better support for rdiv with KronExpansionLayout

Project.toml

@@ -1,6 +1,6 @@
 name = "ContinuumArrays"
 uuid = "7ae1f121-cc2c-504b-ac30-9b923412ae5c"
-version = "0.20.6"
+version = "0.20.7"
 
 [deps]
 AbstractFFTs = "621f4979-c628-5d54-868e-fcf4e3e8185c"

src/ContinuumArrays.jl

@@ -13,7 +13,7 @@ import LinearAlgebra: pinv, inv, dot, norm2, ldiv!, mul!
 import BandedMatrices: AbstractBandedLayout, _BandedMatrix
 import BlockArrays: block, blockindex, unblock, blockedrange, _BlockedUnitRange, _BlockArray, BlockIndexRange, _maybetail
 import FillArrays: AbstractFill, getindex_value, SquareEye
-import ArrayLayouts: mul, ldiv, ZerosLayout, ScalarLayout, AbstractStridedLayout, check_mul_axes, check_ldiv_axes
+import ArrayLayouts: mul, ldiv, ZerosLayout, ScalarLayout, AbstractStridedLayout, check_mul_axes, check_ldiv_axes, Rdiv
 import QuasiArrays: cardinality, checkindex, QuasiAdjoint, QuasiTranspose, Inclusion, SubQuasiArray,
                     QuasiDiagonal, MulQuasiArray, MulQuasiMatrix, MulQuasiVector, QuasiMatMulMat, QuasiArrayLayout,
                     ApplyQuasiArray, ApplyQuasiMatrix, LazyQuasiArrayApplyStyle, AbstractQuasiArrayApplyStyle, AbstractQuasiLazyLayout,
@@ -25,6 +25,9 @@ import QuasiArrays: cardinality, checkindex, QuasiAdjoint, QuasiTranspose, Inclu
 import InfiniteArrays: Infinity, InfAxes
 import AbstractFFTs: Plan
 
+const LazyArraysBandedMatricesExt = Base.get_extension(LazyArrays, :LazyArraysBandedMatricesExt)
+const BandedLazyLayouts = LazyArraysBandedMatricesExt.BandedLazyLayouts
+
 export Spline, LinearSpline, HeavisideSpline, DiracDelta, Derivative, ℵ₁, Inclusion, Basis, grid, plotgrid, affine, .., transform, expand, plan_transform, basis, coefficients,
         weaklaplacian, laplacian, Laplacian, AbsLaplacian, abslaplacian
 

src/bases/bases.jl

@@ -118,6 +118,18 @@ end
 # default to transform for expanding weights
 copy(L::Ldiv{<:AbstractBasisLayout,<:AbstractWeightLayout}) = transform_ldiv(L.A, L.B)
 
+simplifiable(L::Rdiv{<:AbstractLazyLayout,<:AdjointBasisLayout}) = simplifiable(\, L.B', L.A')
+@inline function copy(P::Rdiv{<:AbstractLazyLayout,<:AdjointBasisLayout})
+    A, B = P.A, P.B
+    (B' \ A')'
+end
+
+@inline function copy(P::Rdiv{ApplyLayout{typeof(*)},<:AdjointBasisLayout})
+    A, B = P.A, P.B
+    (B' \ A')'
+end
+
+
 # multiplication operators, reexpand in basis A
 @inline function _broadcast_mul_ldiv(::Tuple{Any,AbstractBasisLayout}, A, B)
     a,b = arguments(B)
@@ -455,6 +467,7 @@ LazyArrays._mul_arguments(::ExpansionLayout, A) = LazyArrays._mul_arguments(Appl
 copy(L::Ldiv{Lay,<:ExpansionLayout}) where Lay<:AbstractBasisLayout = copy(Ldiv{Lay,ApplyLayout{typeof(*)}}(L.A, L.B))
 copy(L::Ldiv{Lay,<:ExpansionLayout}) where Lay<:MappedBasisLayouts = copy(Ldiv{Lay,ApplyLayout{typeof(*)}}(L.A, L.B))
 copy(L::Mul{<:ExpansionLayout,Lay}) where Lay = copy(Mul{ApplyLayout{typeof(*)},Lay}(L.A, L.B))
+copy(L::Mul{<:ExpansionLayout,Lay}) where Lay<:BandedLazyLayouts = copy(Mul{ApplyLayout{typeof(*)},Lay}(L.A, L.B))
 copy(L::Mul{<:ExpansionLayout,Lay}) where Lay<:AbstractLazyLayout = copy(Mul{ApplyLayout{typeof(*)},Lay}(L.A, L.B))
 
 function broadcastbasis_layout(::typeof(+), _, _, a, b)

src/bases/basiskron.jl

@@ -22,4 +22,8 @@ function sub_coefficients_layout(::KronExpansionLayout, P, j)
     X * Bt[:,j]
 end
 
-sum_layout(::KronExpansionLayout, F, dims...) = sum_layout(ApplyLayout{typeof(*)}(), F, dims...)
+sum_layout(::KronExpansionLayout, F, dims...) = sum_layout(ApplyLayout{typeof(*)}(), F, dims...)
+
+diff_layout(::KronExpansionLayout, F, order...; dims...) = diff_layout(ApplyLayout{typeof(*)}(), F, order...; dims...)
+
+copy(L::Ldiv{Lay,<:KronExpansionLayout}) where Lay<:AbstractBasisLayout = copy(Ldiv{Lay,ApplyLayout{typeof(*)}}(L.A, L.B))

test/test_basiskron.jl

@@ -22,4 +22,9 @@ end
     @test F[:,[0.1,0.2]][0.3,:] ≈ F[0.3,[0.1,0.2]] ≈ F[0.3,:][[0.1,0.2]]
     @test F[[0.1,0.2],:][:,0.3] ≈ F[[0.1,0.2],0.3] ≈ F[:,0.3][[0.1,0.2]]
     @test F[[0.1,0.2],[0.3,0.4]] ≈ F[[0.1,0.2],:][:,[0.3,0.4]] ≈ F[:,[0.3,0.4]][[0.1,0.2],:]
+
+    @test diff(F)[0.1,0.2] ≈ diff(F;dims=1)[0.1,0.2]  ≈ diff(L)[0.1,:]'*C*L[0.2,:]
+    @test diff(F;dims=2)[0.1,0.2]  ≈ L[0.1,:]'*C*diff(L)[0.2,:]
+
+    @test L\F/L' == (L\F)/L' == L\(F/L') == C
 end

test/test_splines.jl

@@ -1,6 +1,6 @@
 using ContinuumArrays, LinearAlgebra, Base64, FillArrays, QuasiArrays, BandedMatrices, BlockArrays, StatsBase, Random, Test
 using QuasiArrays: ApplyQuasiArray, ApplyStyle, MemoryLayout, mul, MulQuasiMatrix, Vec
-import LazyArrays: MulStyle, LdivStyle, arguments, applied, apply, simplifiable
+import LazyArrays: MulStyle, LdivStyle, arguments, applied, apply, simplifiable, ApplyArray, Rdiv
 import ContinuumArrays: basis, AdjointBasisLayout, ExpansionLayout, BasisLayout, SubBasisLayout, AdjointMappedBasisLayouts, MappedBasisLayout, plan_grid_transform, weaklaplacian
 
 Random.seed!(24543)
@@ -701,4 +701,17 @@ Random.seed!(24543)
         g = expand(L[affine(0..1,0..1), :], cos)
         @test (f + g)[1/9] ≈ exp(1/9) + cos(1/9)
     end
+
+    @testset "expansion * lazy banded" begin
+        L = LinearSpline(range(0,1,10))
+        F = L * rand(10, 11)
+        B = brand(11,11,1,1)
+        @test F * ApplyArray(*, B, B) == F * B*B
+    end
+
+    @testset "adj and rdiv" begin
+        L = LinearSpline(range(0,1,10))
+        @test L' / L' ≡ Eye(10)
+        @test simplifiable(Rdiv(L', L')) ≡ Val(true)
+    end
 end