# A Julia Linear Operator Package

Operators behave like matrices (with exceptions) but are defined by their effect when applied to a vector. They can be transposed, conjugated, or combined with other operators cheaply. The costly operation is deferred until multiplied with a vector.

## Compatibility

Julia 1.3 and up.

## How to Install

```
pkg> add LinearOperators
pkg> test LinearOperators
```

## Operators Available

Operator | Description |
---|---|

`LinearOperator` | Base class. Useful to define operators from functions |

`TimedLinearOperator` | Linear operator instrumented with timers from `TimerOutputs` |

`BlockDiagonalOperator` | Block-diagonal linear operator |

`opEye` | Identity operator |

`opOnes` | All ones operator |

`opZeros` | All zeros operator |

`opDiagonal` | Square (equivalent to `diagm()` ) or rectangular diagonal operator |

`opInverse` | Equivalent to `\` |

`opCholesky` | More efficient than `opInverse` for symmetric positive definite matrices |

`opLDL` | Similar to `opCholesky` , for general sparse symmetric matrices |

`opHouseholder` | Apply a Householder transformation `I-2hh'` |

`opHermitian` | Represent a symmetric/hermitian operator based on the diagonal and strict lower triangle |

`opRestriction` | Represent a selection of "rows" when composed on the left with an existing operator |

`opExtension` | Represent a selection of "columns" when composed on the right with an existing operator |

`LBFGSOperator` | Limited-memory BFGS approximation in operator form (damped or not) |

`InverseLBFGSOperator` | Inverse of a limited-memory BFGS approximation in operator form (damped or not) |

`LSR1Operator` | Limited-memory SR1 approximation in operator form |

`kron` | Kronecker tensor product in linear operator form |

## Utility Functions

Function | Description |
---|---|

`check_ctranspose` | Cheap check that `A'` is correctly implemented |

`check_hermitian` | Cheap check that `A = A'` |

`check_positive_definite` | Cheap check that an operator is positive (semi-)definite |

`diag` | Extract the diagonal of an operator |

`Matrix` | Convert an abstract operator to a dense array |

`hermitian` | Determine whether the operator is Hermitian |

`push!` | For L-BFGS or L-SR1 operators, push a new pair {s,y} |

`reset!` | For L-BFGS or L-SR1 operators, reset the data |

`shape` | Return the size of a linear operator |

`show` | Display basic information about an operator |

`size` | Return the size of a linear operator |

`symmetric` | Determine whether the operator is symmetric |

`normest` | Estimate the 2-norm |

## Other Operations on Operators

Operators can be transposed (`A.'`

), conjugated (`conj(A)`

) and conjugate-transposed (`A'`

). Operators can be sliced (`A[:,3]`

, `A[2:4,1:5]`

, `A[1,1]`

), but unlike matrices, slices always return operators (see differences).

## Differences

Unlike matrices, an operator never reduces to a vector or a number.

```
A = rand(5,5)
opA = LinearOperator(A)
A[:,1] * 3 # Vector
```

5-element Vector{Float64}: 1.4471654905502307 2.0720315298108285 0.09493729757887093 0.6397410183094988 2.502579300790608

`opA[:,1] * 3 # LinearOperator`

Linear operator nrow: 5 ncol: 1 eltype: Float64 symmetric: false hermitian: false nprod: 0 ntprod: 0 nctprod: 0

`# A[:,1] * [3] # ERROR`

`opA[:,1] * [3] # Vector`

5-element Vector{Float64}: 1.4471654905502307 2.0720315298108285 0.09493729757887093 0.6397410183094988 2.502579300790608

This is also true for `A[i,:]`

, which returns vectors on Julia 0.6, and for the scalar `A[i,j]`

. Similarly, `opA[1,1]`

is an operator of size (1,1):"

`(opA[1,1] * [3])[1] - A[1,1] * 3`

0.0

In the same spirit, the operator `Matrix`

always returns a matrix.

`Matrix(opA[:,1])`

5×1 Matrix{Float64}: 0.4823884968500769 0.6906771766036095 0.03164576585962364 0.21324700610316627 0.834193100263536

## Other Operators

- LLDL features a limited-memory LDLᵀ factorization operator that may be used as preconditioner in iterative methods
- MUMPS.jl features a full distributed-memory factorization operator that may be used to represent the preconditioner in, e.g., constraint-preconditioned Krylov methods.

## Testing

`julia> Pkg.test("LinearOperators")`