Skip to content

Commit 437a207

Browse files
authored
Support tolerance in adaptive QR (#109)
* Support tolerance in adaptive QR * Update for LazyBandedMatrices v0.8 * Update Project.toml * Update test_infqr.jl
1 parent 63d8586 commit 437a207

File tree

6 files changed

+29
-28
lines changed

6 files changed

+29
-28
lines changed

Project.toml

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
name = "InfiniteLinearAlgebra"
22
uuid = "cde9dba0-b1de-11e9-2c62-0bab9446c55c"
3-
version = "0.6.8"
3+
version = "0.6.9"
44

55
[deps]
66
ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"
@@ -17,15 +17,15 @@ MatrixFactorizations = "a3b82374-2e81-5b9e-98ce-41277c0e4c87"
1717
SemiseparableMatrices = "f8ebbe35-cbfb-4060-bf7f-b10e4670cf57"
1818

1919
[compat]
20-
ArrayLayouts = "0.8.6"
20+
ArrayLayouts = "0.8.9"
2121
BandedMatrices = "0.17"
2222
BlockArrays = "0.16.14"
2323
BlockBandedMatrices = "0.11.5"
2424
DSP = "0.7"
25-
FillArrays = "0.12, 0.13"
25+
FillArrays = "0.13"
2626
InfiniteArrays = "0.12"
2727
LazyArrays = "0.22"
28-
LazyBandedMatrices = "0.7.15, 0.8"
28+
LazyBandedMatrices = "0.8"
2929
MatrixFactorizations = "0.9"
3030
SemiseparableMatrices = "0.3"
3131
julia = "1.6"

src/InfiniteLinearAlgebra.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -28,7 +28,7 @@ import BlockArrays: AbstractBlockVecOrMat, sizes_from_blocks, _length, BlockedUn
2828

2929
import BandedMatrices: BandedMatrix, bandwidths, AbstractBandedLayout, _banded_qr!, _banded_qr, _BandedMatrix, banded_chol!
3030

31-
import LazyBandedMatrices: ApplyBandedLayout, BroadcastBandedLayout, _krontrav_axes, _block_interlace_axes, LazyBandedLayout,AbstractLazyBandedLayout
31+
import LazyBandedMatrices: ApplyBandedLayout, BroadcastBandedLayout, _krontrav_axes, _block_interlace_axes, LazyBandedLayout,AbstractLazyBandedLayout, OneToCumsum
3232

3333
import BlockBandedMatrices: _BlockSkylineMatrix, _BandedMatrix, _BlockSkylineMatrix, blockstart, blockstride,
3434
BlockSkylineSizes, BlockSkylineMatrix, BlockBandedMatrix, _BlockBandedMatrix, BlockTridiagonal,

src/blockbanded/blockbanded.jl

Lines changed: 2 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,4 @@
11
const OneToInfCumsum = InfiniteArrays.RangeCumsum{Int,OneToInf{Int}}
2-
const OneToCumsum = InfiniteArrays.RangeCumsum{Int,OneTo{Int}}
32

43
BlockArrays.sortedunion(::AbstractVector{<:PosInfinity}, ::AbstractVector{<:PosInfinity}) = [∞]
54
function BlockArrays.sortedunion(::AbstractVector{<:PosInfinity}, b)
@@ -12,11 +11,7 @@ function BlockArrays.sortedunion(b, ::AbstractVector{<:PosInfinity})
1211
b
1312
end
1413
BlockArrays.sortedunion(a::OneToInfCumsum, ::OneToInfCumsum) = a
15-
BlockArrays.sortedunion(a::OneToCumsum, ::OneToCumsum) = a
16-
function BlockArrays.sortedunion(a::RangeCumsum{<:Any,<:AbstractRange}, b::RangeCumsum{<:Any,<:AbstractRange})
17-
@assert a == b
18-
a
19-
end
14+
2015

2116

2217
function BlockArrays.sortedunion(a::Vcat{Int,1,<:Tuple{Union{Int,AbstractVector{Int}},<:AbstractRange}},
@@ -78,5 +73,4 @@ BroadcastStyle(::Type{<:PseudoBlockArray{T,N,<:AbstractArray{T,N},<:NTuple{N,Blo
7873
# KronTrav
7974
###
8075

81-
_krontrav_axes(A::NTuple{N,OneToInf{Int}}, B::NTuple{N,OneToInf{Int}}) where N =
82-
@. blockedrange(oneto(length(A)))
76+
_krontrav_axes(A::OneToInf{Int}, B::OneToInf{Int}) where N = blockedrange(oneto(length(A)))

src/infqr.jl

Lines changed: 13 additions & 14 deletions
Original file line numberDiff line numberDiff line change
@@ -273,12 +273,11 @@ function materialize!(M::MatLmulVec{<:QRPackedQLayout{<:AdaptiveLayout{<:Abstrac
273273
B
274274
end
275275

276-
function materialize!(M::MatLmulVec{<:AdjQRPackedQLayout{<:AdaptiveLayout{<:AbstractBlockBandedLayout}},<:PaddedLayout})
276+
function materialize!(M::MatLmulVec{<:AdjQRPackedQLayout{<:AdaptiveLayout{<:AbstractBlockBandedLayout}},<:PaddedLayout}; tolerance=1E-30)
277277
adjA,B_in = M.A,M.B
278278
A = adjA.parent
279279
T = eltype(M)
280280
COLGROWTH = 300 # rate to grow columns
281-
tol = 1E-30
282281
ax1 = axes(A.factors.data.data,1)
283282
B = PseudoBlockVector(B_in, (ax1,))
284283

@@ -296,7 +295,7 @@ function materialize!(M::MatLmulVec{<:AdjQRPackedQLayout{<:AdaptiveLayout{<:Abst
296295
resizedata!(B, CS_max)
297296
mx = maximum(abs,view(B,J:last(blockcolsupport(A.factors.data.data.array,J))))
298297
isnan(mx) && error("Not-a-number encounted")
299-
if J > SB && mx tol
298+
if J > SB && mx tolerance
300299
break
301300
end
302301
partialqr!(A.factors.data, CS_max)
@@ -307,19 +306,19 @@ function materialize!(M::MatLmulVec{<:AdjQRPackedQLayout{<:AdaptiveLayout{<:Abst
307306
JR = last(JR)+1:findblock(ax1,last(jr)+COLGROWTH)
308307
end
309308
end
310-
resizedata_chop!(B, tol)
309+
resizedata_chop!(B, tolerance)
311310
end
312311

313312

314-
function _lmul_copymutable(A::AbstractMatrix{T}, x::AbstractVector{S}) where {T,S}
313+
function _lmul_copymutable(A::AbstractMatrix{T}, x::AbstractVector{S}; kwds...) where {T,S}
315314
TS = promote_op(matprod, T, S)
316-
lmul!(A, Base.copymutable(convert(AbstractVector{TS},x)))
315+
lmul!(A, Base.copymutable(convert(AbstractVector{TS},x)); kwds...)
317316
end
318317

319-
(*)(A::QRPackedQ{T,<:AdaptiveQRFactors}, x::AbstractVector) where {T} = _lmul_copymutable(A, x)
320-
(*)(A::Adjoint{T,<:QRPackedQ{T,<:AdaptiveQRFactors}}, x::AbstractVector) where {T} = _lmul_copymutable(A, x)
321-
(*)(A::QRPackedQ{T,<:AdaptiveQRFactors}, x::LayoutVector) where {T} = _lmul_copymutable(A, x)
322-
(*)(A::Adjoint{T,<:QRPackedQ{T,<:AdaptiveQRFactors}}, x::LayoutVector) where {T} = _lmul_copymutable(A, x)
318+
(*)(A::QRPackedQ{T,<:AdaptiveQRFactors}, x::AbstractVector; kwds...) where {T} = _lmul_copymutable(A, x; kwds...)
319+
(*)(A::Adjoint{T,<:QRPackedQ{T,<:AdaptiveQRFactors}}, x::AbstractVector; kwds...) where {T} = _lmul_copymutable(A, x; kwds...)
320+
(*)(A::QRPackedQ{T,<:AdaptiveQRFactors}, x::LayoutVector; kwds...) where {T} = _lmul_copymutable(A, x; kwds...)
321+
(*)(A::Adjoint{T,<:QRPackedQ{T,<:AdaptiveQRFactors}}, x::LayoutVector; kwds...) where {T} = _lmul_copymutable(A, x; kwds...)
323322

324323
function ldiv!(R::UpperTriangular{<:Any,<:AdaptiveQRFactors}, B::CachedVector{<:Any,<:Any,<:Zeros{<:Any,1}})
325324
n = B.datasize[1]
@@ -343,10 +342,10 @@ end
343342

344343
ldiv!(dest::AbstractVector, F::QR{<:Any,<:AdaptiveQRFactors}, b::AbstractVector) =
345344
ldiv!(F, copyto!(dest, b))
346-
ldiv!(F::QR{<:Any,<:AdaptiveQRFactors}, b::AbstractVector) = ldiv!(F.R, lmul!(F.Q',b))
347-
ldiv!(F::QR{<:Any,<:AdaptiveQRFactors}, b::LayoutVector) = ldiv!(F.R, lmul!(F.Q',b))
348-
\(F::QR{<:Any,<:AdaptiveQRFactors}, B::AbstractVector) = ldiv!(F.R, F.Q'B)
349-
\(F::QR{<:Any,<:AdaptiveQRFactors}, B::LayoutVector) = ldiv!(F.R, F.Q'B)
345+
ldiv!(F::QR{<:Any,<:AdaptiveQRFactors}, b::AbstractVector; kwds...) = ldiv!(F.R, lmul!(F.Q',b; kwds...))
346+
ldiv!(F::QR{<:Any,<:AdaptiveQRFactors}, b::LayoutVector; kwds...) = ldiv!(F.R, lmul!(F.Q',b; kwds...))
347+
\(F::QR{<:Any,<:AdaptiveQRFactors}, B::AbstractVector; kwds...) = ldiv!(F.R, *(F.Q', B; kwds...))
348+
\(F::QR{<:Any,<:AdaptiveQRFactors}, B::LayoutVector; kwds...) = ldiv!(F.R, *(F.Q', B; kwds...))
350349

351350

352351
factorize(A::BandedMatrix{<:Any,<:Any,<:OneToInf}) = qr(A)

test/runtests.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -100,7 +100,7 @@ include("test_infbanded.jl")
100100
A = KronTrav- 2I, Eye(∞))
101101
@test axes(A, 1) isa InfiniteLinearAlgebra.OneToInfBlocks
102102
V = view(A, Block.(Base.OneTo(3)), Block.(Base.OneTo(3)))
103-
@test MemoryLayout(V) isa BlockBandedMatrices.BandedBlockBandedLayout
103+
@test MemoryLayout(V) isa LazyBandedMatrices.KronTravBandedBlockBandedLayout
104104

105105
u = A * [1; zeros(∞)]
106106
@test u[1:3] == A[1:3, 1]

test/test_infqr.jl

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -157,6 +157,14 @@ import SemiseparableMatrices: AlmostBandedLayout, VcatAlmostBandedLayout
157157
= BandedMatrix(0 => 1:∞, 1=> Ones(∞), -1=> Ones(∞))
158158
@test qr(A).R[1:10,1:10] qr(Ã).R[1:10,1:10]
159159
end
160+
161+
@testset "AbstractVector b" begin
162+
A = _BandedMatrix(Vcat(Ones(1,∞), (1:∞)', Ones(1,∞)), ℵ₀, 1, 1)
163+
F = qr(A);
164+
b = [1.; 2; 3; zeros(∞)]
165+
@test F\b F\view(b,:)
166+
@test_broken F\b ldiv!(F, view(copy(b),:))
167+
end
160168
end
161169

162170
@testset "almost-banded" begin

0 commit comments

Comments
 (0)