Version history
What's new in v3.11
- The
trfunction fromLinearAlgebra.jlis now overloaded both for genericLinearMaptypes and specialized for most providedLinearMaptypes. In the generic case, this is computationally as expensive as computing the whole matrix representation, though the latter is, of course, not stored.
What's new in v3.10
- A new
MulStyletrait calledTwoArghas been added. It should be used forLinearMaps that do not admit a mutating multiplication à la (3-arg or 5-arg)mul!, but only out-of-place multiplication à laA * x. Products (akaCompositeMaps) and sums (akaLinearCombinations) ofTwoArg-LinearMaps now have memory-optimized multiplication kernels. For instance,A*B*C*xfor threeTwoArg-LinearMapsA,BandCnow allocates onlyy = C*x,z = B*yand the result ofA*z. - The construction of function-based
LinearMaps, typedFunctionMap, has been rearranged. Additionally to the convenience constructorLinearMap{T=Float64}(f, [fc,] M, N=M; kwargs...), the newly exported constructorFunctionMap{T,iip}(f, [fc], M, N; kwargs...)is readily available. Here,iipis eithertrueorfalse, and encodes whetherf(andfcif present) are mutating functions. In the convenience constructor, this is determined via theBoolkeyword argumentismutatingand may not be fully inferred.
What's new in v3.9
- The application of
LinearMaps to vectors operation, i.e.,(A,x) -> A*x = A(x), is now differentiable w.r.t. to the inputxfor integration with machine learning frameworks such asFlux.jl. The reverse differentiation rule makesA::LinearMapusable as a static, i.e., non-trainable, layer in a network, and requires the adjointA'ofAto be defined. - New map types called
KhatriRaoMapandFaceSplittingMapare introduced. These correspond to lazy representations of the column-wise Kronecker product and the row-wise Kronecker product (or "transposed Khatri-Rao product"), respectively. They can be constructed from two matricesAandBviakhatrirao(A, B)andfacesplitting(A, B), respectively. The first is particularly efficient as it makes use of the vec-trick for Kronecker products and computesy = khatrirao(A, B) * xfor a vectorxasy = vec(B * Diagonal(x) * transpose(A)). As such, the Khatri-Rao product can actually be built for generalLinearMaps, including function-based types. Even for moderate sizes of 5 or more columns, this map-vector product is faster than first creating the explicit Khatri-Rao product in memory and then multiplying with the vector; not to mention the memory savings. Unfortunately, similar efficiency cannot be achieved for the face-splitting product.
What's new in v3.8
- A new map called
InverseMapis introduced. Letting anInverseMapact on a vector is equivalent to solving the linear system, i.e.InverseMap(A) * bis the same asA \ b. The default solver isldiv!, but can be specified with thesolverkeyword argument to the constructor (see the docstring for details). Note thatAmust be compatible with the solver:Acan, for example, be a factorization, or anotherLinearMapin combination with an iterative solver. - New constructors for lazy representations of Kronecker products (
squarekron) and sums (sumkronsum) for square factors and summands, respectively, are introduced. They target cases with 3 or more factors/summands, and benchmarking intended use cases for comparison withKroneckerMap(constructed viaBase.kron) andKroneckerSumMap(constructed viakronsum) is recommended.
What's new in v3.7
mul!(M::AbstractMatrix, A::LinearMap, s::Number, a, b)methods are provided, mimicking similar methods inBase.LinearAlgebra. This version allows for the memory efficient implementation of in-place addition and conversion of aLinearMaptoMatrix. Efficient specialisations forWrappedMap,ScaledMap, andLinearCombinationare provided. If users supply the corresponding_unsafe_mul!method for their custom maps, conversion, construction, and inplace addition will benefit from this supplied efficient implementation. If no specialisation is supplied, a generic fallback is used that is based on feeding the canonical basis of unit vectors to the linear map.A new map type called
EmbeddedMapis introduced. It is a wrapper of a "small"LinearMap(or a suitably convertedAbstractVecOrMat) embedded into a "larger" zero map. Hence, the "small" map acts only on a subset of the coordinates and maps to another subset of the coordinates of the "large" map. The large mapLcan therefore be imagined as the composition of a sampling/projection mapP, of the small mapA, and of an embedding mapE:L = E ⋅ A ⋅ P. It is implemented, however, by acting on a view of the vectorxand storing the result into a view of the result vectory. Such maps can be constructed by the new methods:LinearMap(A::MapOrVecOrMat, dims::Dims{2}, index::NTuple{2, AbstractVector{Int}}), wheredimsis the dimensions of the "large" map and index is a tuple of thex- andy-indices that interact withA, respectively;LinearMap(A::MapOrVecOrMat, dims::Dims{2}; offset::Dims{2}), where the keyword argumentoffsetdetermines the dimension of a virtual upper-left zero block, to whichAgets (virtually) diagonally appended.
An often requested new feature has been added: slicing (i.e., non-scalar indexing) any
LinearMapobject viaBase.getindexoverloads. Note, however, that only rather efficient complete slicing operations are implemented:A[:,j],A[:,J], andA[:,:], wherej::IntegerandJis either of typeAbstractVector{<:Integer>}or anAbstractVector{Bool}of appropriate length ("logical slicing"). Partial slicing operations such asA[I,j]andA[I,J]whereIis asJabove are disallowed.Scalar indexing
A[i::Integer,j::Integer]as well as other indexing operations that fall back on scalar indexing such as logical indexing by someAbstractMatrix{Bool}, or indexing by vectors of (linear or Cartesian) indices are not supported; as an exception,getindexcalls on wrappedAbstractVecOrMats is forwarded to correspondinggetindexmethods fromBaseand therefore allow any type of usual indexing/slicing. If scalar indexing is really required, consider usingA[:,j][i]which is as efficient as a reasonable generic implementation forLinearMaps can be.Furthermore, (predominantly) horizontal slicing operations require the adjoint operation of the
LinearMaptype to be defined, or will fail otherwise. Important note:LinearMapobjects are meant to model objects that act on vectors efficiently, and are in general not backed up by storage-like types likeArrays. Therefore, slicing ofLinearMaps is potentially slow, and it may require the (repeated) allocation of standard unit vectors. As a consequence, generic algorithms relying heavily on indexing and/or slicing are likely to run much slower than expected forAbstractArrays. To avoid repeated indexing operations which may involve redundant computations, it is strongly recommended to considerconvertingLinearMap-typed objects toMatrixorSparseMatrixCSCfirst, if memory permits.
What's new in v3.6
- Support for Julia versions below v1.6 has been dropped.
Block[Diagonal]Map,CompositeMap,KroneckerMapandLinearCombinationtype objects can now be backed by aVectorofLinearMap-type elements. This can be beneficial in cases where these higher-orderLinearMaps are constructed from many maps where a tuple backend may get inefficient or impose hard work for the compiler at construction. The default behavior, however, does not change, and construction of vector-basedLinearMaps requires usage of the unexported constructors ("expert usage"), except for constructions likesum([A, B, C])orprod([A, B, C])(== C*B*A), whereA,BandCare of someLinearMaptype.
What's new in v3.5
WrappedMap,ScaledMap,LinearCombination,AdjointMap,TransposeMapandCompositeMap, instead of using the defaultaxes(A) = map(oneto, size(A)), now forward calls toaxesto the underlying wrapped linear map. This allows allocating operations such as*to determine the appropriate storage and axes type of their outputs. For example, linear maps that wrapBlockArrayswill, upon multiplicative action, produce aBlockArrays.PseudoBlockVectorwith block structure inherited from the operator's output axesaxes(A,1).
What's new in v3.4
- In
WrappedMapconstructors, as implicitly called in addition and multiplication ofLinearMaps andAbstractMatrixobjects, (conjugate) symmetry and positive definiteness are only determined for matrix types for which these checks are expected to be very cheap or even known at compile time based on the concrete type. The default forLinearMapsubtypes is to call, for instance,issymmetric, because symmetry properties are either stored or easily obtained from constituting maps. For custom matrix types, define corresponding methodsLinearMaps._issymmetric,LinearMaps._ishermitianandLinearMaps._isposdefto hook into the property checking mechanism.
What's new in v3.3
AbstractVectors can now be wrapped by aLinearMapjust likeAbstractMatrixtyped objects. Upon wrapping, they are not implicitly reshaped to matrices. This feature might be helpful, for instance, in the lazy representation of rank-1 operatorskron(LinearMap(u), v') == ⊗(u, v') == u ⊗ v'for vectorsuandv. The action on vectors,(u⊗v')x, is implemented optimally viau(v'x)`.
What's new in v3.2
- In-place left-multiplication
mul!(Y, X, A::LinearMap)is now allowed forX::AbstractMatrixand implemented via the adjoint equationY' = A'X'.
What's new in v3.1
- In Julia v1.3 and above,
LinearMap-typed objects are callable onAbstractVectors: ForL::LinearMapandx::AbstractVector,L(x) = L*x.
What's new in v3.0
- BREAKING change: Internally, any dependence on former
A*_mul_B!methods is abandoned. For customLinearMapsubtypes, there are now two options:- In case your type is invariant under adjoint/transposition (i.e.,
adjoint(L::MyLinearMap)::MyLinearMapsimilar to, for instance,LinearCombinations orCompositeMaps),At_mul_B!andAc_mul_B!do not require any replacement! Rather, multiplication byL'is, in this case, handled bymul!(y, L::MyLinearMap, x[, α, β]). - Otherwise, you will need to define
mul!methods with the signaturemul!(y, L::TransposeMap{<:Any,MyLinearMap}, x[, α, β])andmul!(y, L::AdjointMap{<:Any,MyLinearMap}, x[, α, β]).
- In case your type is invariant under adjoint/transposition (i.e.,
- Left multiplying by a transpose or adjoint vector (e.g.,
y'*A) produces a transpose or adjoint vector output, rather than a compositeLinearMap. - Block concatenation now handles matrices and vectors directly by internal promotion to
LinearMaps. For[h/v/hc]catit suffices to have aLinearMapobject anywhere in the list of arguments. For the block-diagonal concatenation viaSparseArrays.blockdiag, aLinearMapobject has to appear among the first 8 arguments. This restriction, however, does not apply to block-diagonal concatenation viaBase.cat(As...; dims=(1,2)). - Introduction of more expressive and visually appealing
showmethods, replacing the fallback to the genericshow.
What's new in v2.7
- Potential reduction of memory allocations in multiplication of
LinearCombinations,BlockMaps, and real- or complex-scaledLinearMaps. For the latter, a new internal typeScaledMaphas been introduced. - Multiplication code for
CompositeMaps has been refactored to facilitate to provide memory for storage of intermediate results by directly calling helper functions.
What's new in v2.6
- New feature: "lazy" Kronecker product, Kronecker sums, and powers thereof for
LinearMaps.AbstractMatrixobjects are promoted toLinearMaps if one of the first 8 Kronecker factors is aLinearMapobject. - Compatibility with the generic multiply-and-add interface (a.k.a. 5-arg
mul!) introduced in julia v1.3
What's new in v2.5
- New feature: concatenation of
LinearMaps objects withUniformScalings, consistent with (h-, v-, and hc-)concatenation of matrices. Note, matricesAmust be wrapped asLinearMap(A),UniformScalings are promoted toLinearMaps automatically.
What's new in v2.4
- Support restricted to Julia v1.0+.
What's new in v2.3
- Fully Julia v0.7/v1.0/v1.1 compatible.
- Full support of noncommutative number types such as quaternions.
What's new in v2.2
- Fully Julia v0.7/v1.0 compatible.
- A
convert(SparseMatrixCSC, A::LinearMap)function, that calls thesparsematrix generating function.
What's new in v2.1
- Fully Julia v0.7 compatible; dropped compatibility for previous versions of Julia from LinearMaps.jl v2.0.0 on.
- A 5-argument version for
mul!(y, A::LinearMap, x, α=1, β=0), which computesy := α * A * x + β * yand implements the usual 3-argumentmul!(y, A, x)for the defaultαandβ. - Synonymous
convert(Matrix, A::LinearMap)andconvert(Array, A::LinearMap)functions, that call theMatrixconstructor and return the matrix representation ofA. - Multiplication with matrices, interpreted as a block row vector of vectors:
mul!(Y::AbstractArray, A::LinearMap, X::AbstractArray, α=1, β=0): appliesAto each column ofXand stores the result in-place in the corresponding column ofY;- for the out-of-place multiplication, the approach is to compute
convert(Matrix, A * X); this is equivalent to applyingAto each column ofX. In generic code which handles bothA::AbstractMatrixandA::LinearMap, the additional call toconvertis a noop whenAis a matrix.
- Full compatibility with Arpack.jl's
eigsandsvds; previously onlyeigswas working. For more, nicely collaborating packages see the Example section.