diff --git a/Project.toml b/Project.toml index 554b8d66..410a66fa 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "BandedMatrices" uuid = "aae01518-5342-5314-be14-df237901396f" -version = "1.9.4" +version = "1.9.5" [deps] ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a" diff --git a/src/BandedMatrices.jl b/src/BandedMatrices.jl index 4f2d9f45..f0a20e29 100644 --- a/src/BandedMatrices.jl +++ b/src/BandedMatrices.jl @@ -19,7 +19,7 @@ import LinearAlgebra: _apply_inverse_ipiv_rows!, _apply_ipiv_rows!, _chol!, axpy logabsdet, lu, lu!, mul!, qr, qr!, rmul!, svdvals, svdvals!, tril!, triu! using LinearAlgebra: AbstractTriangular, AdjOrTrans, BlasComplex, BlasFloat, BlasInt, BlasReal, Givens, HermOrSym, - QRPackedQ, RealHermSymComplexHerm, StructuredMatrixStyle, checksquare, chkstride1, ipiv2perm + QRPackedQ, RealHermSymComplexHerm, StructuredMatrixStyle, checksquare, chkstride1, ipiv2perm, AdjointQ using LinearAlgebra.LAPACK @@ -32,7 +32,7 @@ import ArrayLayouts: AbstractTridiagonalLayout, BidiagonalLayout, BlasMatLdivVec hermitianlayout, materialize, materialize!, reflector!, reflectorApply!, rowsupport, sub_materialize, subdiagonaldata, sublayout, supdiagonaldata, symmetricdata, symmetriclayout, symmetricuplo, transposelayout, triangulardata, triangularlayout, zero!, - QRPackedQLayout, AdjQRPackedQLayout + QRPackedQLayout, AdjQRPackedQLayout, LayoutVecOrMats import FillArrays: AbstractFill, getindex_value, _broadcasted_zeros, unique_value, OneElement, RectDiagonal, OneElementMatrix, OneElementVector, ZerosMatrix, ZerosVector diff --git a/src/banded/bandedqr.jl b/src/banded/bandedqr.jl index 4ae5e08f..96751872 100644 --- a/src/banded/bandedqr.jl +++ b/src/banded/bandedqr.jl @@ -1,5 +1,3 @@ -AdjointQType = isdefined(LinearAlgebra, :AdjointQ) ? LinearAlgebra.AdjointQ : Adjoint - _qr(::AbstractBandedLayout, ax, A) = _banded_qr(ax, A) _banded_qr(_, A) = qr!(BandedMatrix{float(eltype(A))}(A, (bandwidth(A,1),bandwidth(A,1)+bandwidth(A,2)))) @@ -65,7 +63,7 @@ function banded_qr_lmul!(A, B) B end -function banded_qr_lmul!(adjA::AdjointQType, B) +function banded_qr_lmul!(adjA::AdjointQ, B) require_one_based_indexing(B) A = parent(adjA) mA, nA = size(A.factors) @@ -121,7 +119,7 @@ function banded_qr_rmul!(A, Q) end A end -function banded_qr_rmul!(A, adjQ::AdjointQType) +function banded_qr_rmul!(A, adjQ::AdjointQ) Q = parent(adjQ) mQ, nQ = size(Q.factors) mA, nA = size(A,1), size(A,2) @@ -150,27 +148,28 @@ function banded_qr_rmul!(A, adjQ::AdjointQType) end banded_lmul!(A::QRPackedQ, B::AbstractVecOrMat) = banded_qr_lmul!(A, B) -banded_lmul!(adjA::AdjointQType{<:Any,<:QRPackedQ}, B::AbstractVecOrMat) = banded_qr_lmul!(adjA, B) +banded_lmul!(adjA::AdjointQ{<:Any,<:QRPackedQ}, B::AbstractVecOrMat) = banded_qr_lmul!(adjA, B) banded_rmul!(A::AbstractMatrix, Q::QRPackedQ) = banded_qr_rmul!(A, Q) -banded_rmul!(A::AbstractMatrix, adjQ::AdjointQType{<:Any,<:QRPackedQ}) = banded_qr_rmul!(A, adjQ) +banded_rmul!(A::AbstractMatrix, adjQ::AdjointQ{<:Any,<:QRPackedQ}) = banded_qr_rmul!(A, adjQ) lmul!(A::QRPackedQ{<:Any,<:AbstractBandedMatrix}, B::AbstractVector) = banded_lmul!(A,B) lmul!(A::QRPackedQ{<:Any,<:AbstractBandedMatrix}, B::AbstractMatrix) = banded_lmul!(A,B) -lmul!(adjA::AdjointQType{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}, B::AbstractVector) = banded_lmul!(adjA,B) -lmul!(adjA::AdjointQType{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}, B::AbstractMatrix) = banded_lmul!(adjA,B) -lmul!(A::QRPackedQ{<:Any,BandedSubBandedMatrix{T,C,R,I1,I2}}, B::AbstractVector) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange} = +lmul!(adjA::AdjointQ{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}, B::AbstractVector) = banded_lmul!(adjA,B) +lmul!(adjA::AdjointQ{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}, B::AbstractMatrix) = banded_lmul!(adjA,B) +lmul!(adjA::AdjointQ{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}, B::LayoutVecOrMats) = banded_lmul!(adjA,B) +lmul!(A::QRPackedQ{<:Any,<:BandedSubBandedMatrix{T,C,R,I1,I2}}, B::AbstractVector) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange} = banded_lmul!(A,B) -lmul!(A::QRPackedQ{<:Any,BandedSubBandedMatrix{T,C,R,I1,I2}}, B::AbstractMatrix) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange} = +lmul!(A::QRPackedQ{<:Any,<:BandedSubBandedMatrix{T,C,R,I1,I2}}, B::AbstractMatrix) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange} = banded_lmul!(A,B) -lmul!(adjA::AdjointQType{T,<:QRPackedQ{T,<:BandedSubBandedMatrix{T,C,R,I1,I2,t}}}, B::AbstractVector) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange,t} = +lmul!(adjA::AdjointQ{T,<:QRPackedQ{T,<:BandedSubBandedMatrix{T,C,R,I1,I2,t}}}, B::AbstractVector) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange,t} = banded_lmul!(adjA,B) -lmul!(adjA::AdjointQType{T,<:QRPackedQ{T,<:BandedSubBandedMatrix{T,C,R,I1,I2,t}}}, B::AbstractMatrix) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange,t} = +lmul!(adjA::AdjointQ{T,<:QRPackedQ{T,<:BandedSubBandedMatrix{T,C,R,I1,I2,t}}}, B::AbstractMatrix) where {T,C,R,I1<:AbstractUnitRange,I2<:AbstractUnitRange,t} = banded_lmul!(adjA,B) -# rmul!(A::AbstractMatrix, adjQ::AdjointQType{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}) = banded_rmul!(A, adjA) -# rmul!(A::StridedMatrix, adjQ::AdjointQType{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}) = banded_rmul!(A, adjA) +# rmul!(A::AbstractMatrix, adjQ::AdjointQ{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}) = banded_rmul!(A, adjA) +# rmul!(A::StridedMatrix, adjQ::AdjointQ{<:Any,<:QRPackedQ{<:Any,<:AbstractBandedMatrix}}) = banded_rmul!(A, adjA) rmul!(A::StridedVecOrMat{T}, Q::QRPackedQ{T,B}) where {T<:BlasFloat,B<:AbstractBandedMatrix{T}} = banded_rmul!(A, Q) -rmul!(A::StridedVecOrMat{T}, adjQ::AdjointQType{<:Any,<:QRPackedQ{T,B}}) where {T<:BlasComplex,B<:AbstractBandedMatrix{T}} = banded_rmul!(A, adjQ) -rmul!(A::StridedVecOrMat{T}, adjQ::AdjointQType{<:Any,<:QRPackedQ{T,B}}) where {T<:BlasReal,B<:AbstractBandedMatrix{T}} = banded_rmul!(A, adjQ) +rmul!(A::StridedVecOrMat{T}, adjQ::AdjointQ{<:Any,<:QRPackedQ{T,B}}) where {T<:BlasComplex,B<:AbstractBandedMatrix{T}} = banded_rmul!(A, adjQ) +rmul!(A::StridedVecOrMat{T}, adjQ::AdjointQ{<:Any,<:QRPackedQ{T,B}}) where {T<:BlasReal,B<:AbstractBandedMatrix{T}} = banded_rmul!(A, adjQ) function _banded_widerect_ldiv!(A::QR{T}, B) where T diff --git a/test/test_bandedqr.jl b/test/test_bandedqr.jl index a2dbf355..c8c7e15a 100644 --- a/test/test_bandedqr.jl +++ b/test/test_bandedqr.jl @@ -150,6 +150,17 @@ Random.seed!(0) Q = BandedMatrices.QRPackedQ(BandedMatrix(-2 => -(s ./ (c .+ 1))), c .+ 1) @test Q'Q ≈ Matrix(Q)'Matrix(Q) ≈ I end + + @testset "bandedsubbanded qr " begin + A = brand(100,100,1,1) + B = BandedMatrix(A,(1,2)) # pad + b = randn(5) + c = randn(5,3) + @test lmul!(qr!(view(copy(B),1:5,1:6)).Q, copy(b)) ≈ qr(A[1:5,1:6]).Q * b + @test lmul!(qr!(view(copy(B),1:5,1:6)).Q, copy(c)) ≈ qr(A[1:5,1:6]).Q * c + @test lmul!(qr!(view(copy(B),1:5,1:6)).Q', copy(b)) ≈ qr(A[1:5,1:6]).Q' * b + @test lmul!(qr!(view(copy(B),1:5,1:6)).Q', copy(c)) ≈ qr(A[1:5,1:6]).Q' * c + end end end # module