Skip to content

Commit aba1cc1

Browse files
dlfivefiftydlfivefifty
authored
Use Mul to determine MulAdd, Lmul, or Rmul (#30)
* Fix BandedMatrices #188 * mul(Eye(n), A) == copy(A) * _mul for Eye, return type for ldiv!(::QRPackedQ, _) * v0.3.8 * _mul -> Mul * Update mul.jl * copy for Diagonal, LayoutVector * overloads, DualLayout * generic triangular Lmul * Remove unused copy for Diagonal fills * Simplify Diagonal mul copy * Improve triangular copy * Add mulreduce * getindex for Mul * Update mul.jl * Generic copyto! * support ContinuumArrays.jl in Mul getindex * Add docs * Update mul.jl * SparseArray copyto! * bi/tri/symtridiagonal layouts * Fill *DiagLayout * Update memorylayout.jl * Increase coverage * Update ArrayLayouts.jl * Increase coverage * Dot tests * increase cov
1 parent 96c427f commit aba1cc1

20 files changed

+654
-158
lines changed

‎.travis.yml

Lines changed: 9 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -9,9 +9,17 @@ julia:
99
- 1.4
1010
- 1.5
1111
- nightly
12-
matrix:
12+
jobs:
1313
allow_failures:
1414
- julia: nightly
15+
- stage: "Documentation"
16+
julia: 1.4
17+
os: linux
18+
script:
19+
- julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd()));
20+
Pkg.instantiate()'
21+
- julia --project=docs/ docs/make.jl
22+
after_success: skip
1523
notifications:
1624
email: false
1725
after_success:

‎Project.toml

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1,14 +1,15 @@
11
name = "ArrayLayouts"
22
uuid = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
33
authors = ["Sheehan Olver <[email protected]>"]
4-
version = "0.3.8"
4+
version = "0.4.0"
55

66
[deps]
77
FillArrays = "1a297f60-69ca-5386-bcde-b61e274b549b"
88
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
9+
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
910

1011
[compat]
11-
FillArrays = "0.8.10, 0.9"
12+
FillArrays = "0.9.1"
1213
julia = "1"
1314

1415
[extras]

‎README.md

Lines changed: 22 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -18,4 +18,25 @@ julia> @time mul!(y, A, x); # Julia does not recognize that V is symmetric
1818

1919
julia> @time muladd!(1.0, V, x, 0.0, y); # ArrayLayouts does and is 3x faster as it calls BLAS
2020
0.017677 seconds (4 allocations: 160 bytes)
21-
```
21+
```
22+
23+
## Internal design
24+
25+
The package is based on assigning a `MemoryLayout` to every array, which is used for dispatch. For example,
26+
```julia
27+
julia> MemoryLayout(A) # Each column of A is column major, and columns stored in order
28+
DenseColumnMajor()
29+
30+
julia> MemoryLayout(view(A, 1:3,:)) # Each column of A is column major
31+
ColumnMajor()
32+
33+
julia> MemoryLayout(V) # A symmetric version, whose storage is DenseColumnMajor
34+
SymmetricLayout{DenseColumnMajor}()
35+
```
36+
37+
This is then used by `muladd!(α, A, B, β, C)`, `ArrayLayouts.lmul!(A, B)`, and `ArrayLayouts.rmul!(A, B)` to
38+
lower to the correct BLAS calls via lazy objects `MulAdd(α, A, B, β, C)`, `Lmul(A, B)`, `Rmul(A, B)` which are materialized,
39+
in analogy to `Base.Broadcasted`.
40+
41+
Note there is also a higher level function `mul(A, B)` that materializes via `Mul(A, B)`, which uses the layout of `A` and `B`
42+
to further reduce to either `MulAdd`, `Lmul`, and `Rmul`.

‎docs/.gitignore

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,2 @@
1+
build/
2+
site/

‎docs/Project.toml

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,2 @@
1+
[deps]
2+
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"

‎docs/make.jl

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
using Documenter
2+
using ArrayLayouts
3+
4+
makedocs(
5+
sitename = "ArrayLayouts",
6+
format = Documenter.HTML(),
7+
modules = [ArrayLayouts]
8+
)
9+
10+
# Documenter can also automatically deploy documentation to gh-pages.
11+
# See "Hosting Documentation" and deploydocs() in the Documenter manual
12+
# for more information.
13+
#=deploydocs(
14+
repo = "<repository url>"
15+
)=#

‎docs/src/index.md

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,3 @@
1+
# ArrayLayouts.jl
2+
3+
Documentation for ArrayLayouts.jl

‎src/ArrayLayouts.jl

Lines changed: 19 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
11
module ArrayLayouts
2-
using Base, Base.Broadcast, LinearAlgebra, FillArrays
2+
using Base, Base.Broadcast, LinearAlgebra, FillArrays, SparseArrays
33
import LinearAlgebra.BLAS
44

55
import Base: AbstractArray, AbstractMatrix, AbstractVector,
@@ -33,8 +33,9 @@ import Base.Broadcast: BroadcastStyle, AbstractArrayStyle, Broadcasted, broadcas
3333
combine_eltypes, DefaultArrayStyle, instantiate, materialize,
3434
materialize!, eltypes
3535

36-
import LinearAlgebra: AbstractTriangular, AbstractQ, checksquare, pinv, fill!, tilebufsize, Abuf, Bbuf, Cbuf, dot, factorize, qr, lu, cholesky,
37-
norm2, norm1, normInf, normMinusInf, qr, lu, qr!, lu!, AdjOrTrans, HermOrSym
36+
import LinearAlgebra: AbstractTriangular, AbstractQ, checksquare, pinv, fill!, tilebufsize, Abuf, Bbuf, Cbuf, factorize, qr, lu, cholesky,
37+
norm2, norm1, normInf, normMinusInf, qr, lu, qr!, lu!, AdjOrTrans, HermOrSym, copy_oftype,
38+
AdjointAbsVec, TransposeAbsVec
3839

3940
import LinearAlgebra.BLAS: BlasFloat, BlasReal, BlasComplex
4041

@@ -47,10 +48,12 @@ else
4748
import Base: require_one_based_indexing
4849
end
4950

50-
export materialize, materialize!, MulAdd, muladd!, Ldiv, Rdiv, Lmul, Rmul, lmul, rmul, ldiv, rdiv, mul, MemoryLayout, AbstractStridedLayout,
51+
export materialize, materialize!, MulAdd, muladd!, Ldiv, Rdiv, Lmul, Rmul, Dot,
52+
lmul, rmul, mul, ldiv, rdiv, mul, MemoryLayout, AbstractStridedLayout,
5153
DenseColumnMajor, ColumnMajor, ZerosLayout, FillLayout, AbstractColumnMajor, RowMajor, AbstractRowMajor, UnitStride,
52-
DiagonalLayout, ScalarLayout, SymTridiagonalLayout, HermitianLayout, SymmetricLayout, TriangularLayout,
53-
UnknownLayout, AbstractBandedLayout, ApplyBroadcastStyle, ConjLayout, AbstractFillLayout,
54+
DiagonalLayout, ScalarLayout, SymTridiagonalLayout, TridiagonalLayout, BidiagonalLayout,
55+
HermitianLayout, SymmetricLayout, TriangularLayout,
56+
UnknownLayout, AbstractBandedLayout, ApplyBroadcastStyle, ConjLayout, AbstractFillLayout, DualLayout,
5457
colsupport, rowsupport, layout_getindex, QLayout, LayoutArray, LayoutMatrix, LayoutVector
5558

5659
struct ApplyBroadcastStyle <: BroadcastStyle end
@@ -94,6 +97,7 @@ function unsafe_convert(::Type{ConjPtr{T}}, V::SubArray{T,2}) where {T,N,P}
9497
end
9598

9699
include("memorylayout.jl")
100+
include("mul.jl")
97101
include("muladd.jl")
98102
include("lmul.jl")
99103
include("ldiv.jl")
@@ -115,6 +119,7 @@ macro _layoutgetindex(Typ)
115119
@inline Base.getindex(A::$Typ, kr::AbstractUnitRange, jr::AbstractUnitRange) = ArrayLayouts.layout_getindex(A, kr, jr)
116120
@inline Base.getindex(A::$Typ, kr::AbstractVector, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr)
117121
@inline Base.getindex(A::$Typ, kr::Colon, jr::AbstractVector) = ArrayLayouts.layout_getindex(A, kr, jr)
122+
@inline Base.getindex(A::$Typ, kr::Colon, jr::Integer) = ArrayLayouts.layout_getindex(A, kr, jr)
118123
@inline Base.getindex(A::$Typ, kr::AbstractVector, jr::Colon) = ArrayLayouts.layout_getindex(A, kr, jr)
119124
end)
120125
end
@@ -162,8 +167,13 @@ copyto!(dest::SubArray{<:Any,N,<:LayoutArray}, src::AbstractArray{<:Any,N}) wher
162167
_copyto!(MemoryLayout(dest), MemoryLayout(src), dest, src)
163168
copyto!(dest::AbstractArray{<:Any,N}, src::SubArray{<:Any,N,<:LayoutArray}) where N =
164169
_copyto!(MemoryLayout(dest), MemoryLayout(src), dest, src)
165-
166-
170+
# ambiguity from sparsematrix.jl
171+
if VERSION ≥ v"1.5"
172+
copyto!(dest::LayoutMatrix, src::SparseArrays.AbstractSparseMatrixCSC) =
173+
_copyto!(MemoryLayout(dest), MemoryLayout(src), dest, src)
174+
copyto!(dest::SubArray{<:Any,2,<:LayoutMatrix}, src::SparseArrays.AbstractSparseMatrixCSC) =
175+
_copyto!(MemoryLayout(dest), MemoryLayout(src), dest, src)
176+
end
167177

168178
zero!(A::AbstractArray{T}) where T = fill!(A,zero(T))
169179
function zero!(A::AbstractArray{<:AbstractArray})
@@ -245,7 +255,7 @@ Base.replace_in_print_matrix(A::Union{LayoutVector,
245255
AdjOrTrans{<:Any,<:LayoutMatrix},
246256
HermOrSym{<:Any,<:LayoutMatrix},
247257
SubArray{<:Any,2,<:LayoutMatrix}}, i::Integer, j::Integer, s::AbstractString) =
248-
layout_replace_in_print_matrix(MemoryLayout(A), A, i, j, s)
258+
layout_replace_in_print_matrix(MemoryLayout(A), A, i, j, s)
249259

250260
Base.print_matrix_row(io::IO,
251261
X::Union{LayoutMatrix,

‎src/diagonal.jl

Lines changed: 11 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,10 @@
33
# Lmul
44
####
55

6+
mulreduce(M::Mul{<:DiagonalLayout,<:DiagonalLayout}) = Lmul(M)
7+
mulreduce(M::Mul{<:DiagonalLayout}) = Lmul(M)
8+
mulreduce(M::Mul{<:Any,<:DiagonalLayout}) = Rmul(M)
9+
610
# Diagonal multiplication never changes structure
711
similar(M::Lmul{<:DiagonalLayout}, ::Type{T}, axes) where T = similar(M.B, T, axes)
812
# equivalent to rescaling
@@ -11,10 +15,13 @@ function materialize!(M::Lmul{<:DiagonalLayout{<:AbstractFillLayout}})
1115
M.B
1216
end
1317

14-
copy(M::Lmul{<:DiagonalLayout,<:DiagonalLayout}) = Diagonal(M.A.diag .* M.B.diag)
15-
copy(M::Lmul{<:DiagonalLayout}) = M.A.diag .* M.B
16-
copy(M::Lmul{<:DiagonalLayout{<:AbstractFillLayout}}) = getindex_value(M.A.diag) * M.B
17-
copy(M::Lmul{<:DiagonalLayout{<:AbstractFillLayout},<:DiagonalLayout}) = Diagonal(getindex_value(M.A.diag) * M.B.diag)
18+
19+
copy(M::Lmul{<:DiagonalLayout,<:DiagonalLayout}) = Diagonal(diagonaldata(M.A) .* diagonaldata(M.B))
20+
copy(M::Lmul{<:DiagonalLayout}) = diagonaldata(M.A) .* M.B
21+
copy(M::Rmul{<:Any,<:DiagonalLayout}) = M.A .* permutedims(diagonaldata(M.B))
22+
23+
24+
1825

1926
# Diagonal multiplication never changes structure
2027
similar(M::Rmul{<:Any,<:DiagonalLayout}, ::Type{T}, axes) where T = similar(M.A, T, axes)
@@ -24,8 +31,6 @@ function materialize!(M::Rmul{<:Any,<:DiagonalLayout{<:AbstractFillLayout}})
2431
M.A
2532
end
2633

27-
copy(M::Rmul{<:Any,<:DiagonalLayout}) = M.A .* permutedims(M.B.diag)
28-
copy(M::Rmul{<:Any,<:DiagonalLayout{<:AbstractFillLayout}}) = M.A .* getindex_value(M.B.diag)
2934

3035
function materialize!(M::Ldiv{<:DiagonalLayout})
3136
M.B .= M.A.diag .\ M.B

‎src/factorizations.jl

Lines changed: 7 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -81,8 +81,10 @@ function materialize!(Ldv::Ldiv{<:QRPackedLayout,<:Any,<:Any,<:AbstractMatrix{T}
8181
end
8282
return B
8383
end
84-
materialize!(Ldv::Ldiv{<:QRPackedLayout,<:Any,<:Any,<:AbstractVector{T}}) where T =
85-
ldiv!(Ldv.A, reshape(Ldv.B, length(Ldv.B), 1))[:]
84+
function materialize!(Ldv::Ldiv{<:QRPackedLayout,<:Any,<:Any,<:AbstractVector{T}}) where T
85+
ldiv!(Ldv.A, reshape(Ldv.B, length(Ldv.B), 1))
86+
Ldv.B
87+
end
8688

8789

8890

@@ -101,7 +103,9 @@ adjointlayout(::Type, ::QRPackedQLayout{SLAY,TLAY}) where {SLAY,TLAY} = AdjQRPac
101103
adjointlayout(::Type, ::QRCompactWYQLayout{SLAY,TLAY}) where {SLAY,TLAY} = AdjQRCompactWYQLayout{SLAY,TLAY}()
102104

103105
copy(M::Lmul{<:AbstractQLayout}) = copyto!(similar(M), M)
104-
106+
mulreduce(M::Mul{<:AbstractQLayout,<:AbstractQLayout}) = MulAdd(M)
107+
mulreduce(M::Mul{<:AbstractQLayout}) = Lmul(M)
108+
mulreduce(M::Mul{<:Any,<:AbstractQLayout}) = Rmul(M)
105109

106110
function copyto!(dest::AbstractArray{T}, M::Lmul{<:AbstractQLayout}) where T
107111
A,B = M.A,M.B

0 commit comments

Comments
 (0)