Skip to content

Commit 8dce3b8

Browse files
dlfivefiftydlfivefifty
committed
work on inf qr
1 parent 51914c6 commit 8dce3b8

File tree

4 files changed

+164
-158
lines changed

4 files changed

+164
-158
lines changed

‎test/runtests.jl

Lines changed: 1 addition & 73 deletions
Original file line numberDiff line numberDiff line change
@@ -41,78 +41,6 @@ import BandedMatrices: bandeddata, _BandedMatrix
4141
@test Ac[1:10,1:10] ≈ (A')[1:10,1:10] ≈ A[1:10,1:10]'
4242
@test At[1:10,1:10] ≈ transpose(A)[1:10,1:10] ≈ transpose(A[1:10,1:10])
4343
end
44-
4544
include("test_hessenbergq.jl")
45+
include("test_infql.jl")
4646

47-
@testset "PertTriToeplitz QL" begin
48-
A = Tridiagonal(Vcat(Float64[], Fill(2.0,∞)),
49-
Vcat(Float64[2.0], Fill(0.0,∞)),
50-
Vcat(Float64[], Fill(0.5,∞)))
51-
for λ in (-2.1-0.01im,-2.1+0.01im,-2.1+eps()im,-2.1-eps()im,-2.1+0.0im,-2.1-0.0im,-1.0+im,-3.0+im,-3.0-im)
52-
Q, L = ql(A - λ*I)
53-
@test Q[1:10,1:12]*L[1:12,1:10] ≈ A[1:10,1:10] - λ*I
54-
end
55-
end
56-
57-
@testset "Pert Hessenberg Toeplitz" begin
58-
a = [1,2,5,0.5]
59-
Random.seed!(0)
60-
A = _BandedMatrix(Hcat(randn(4,2), reverse(a) * Ones(1,∞)), ∞, 2, 1)
61-
@test A isa PertToeplitz
62-
@test BandedMatrix(A, (3,1))[1:10,1:10] == A[1:10,1:10]
63-
Q,L = ql(A)
64-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
65-
66-
a = [0.1,1,2,3,0.5]
67-
A = _BandedMatrix(Hcat([0.5 0.5; -1 3; 2 2; 1 1; 0.1 0.1], reverse(a) * Ones(1,∞)), ∞, 3, 1)
68-
@test A isa PertToeplitz
69-
@test BandedMatrix(A, (3,1))[1:10,1:10] == A[1:10,1:10]
70-
71-
B = BandedMatrix(A, (3,1))
72-
@test B[1:10,1:10] == A[1:10,1:10]
73-
Q,L = ql(A)
74-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
75-
76-
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)),
77-
0 => Vcat([0.0], Fill(1/2,∞)),
78-
1 => Vcat([1/4], Fill(1/4,∞)))
79-
Q,L = ql(A)
80-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
81-
a = [-0.1,0.2,0.3]
82-
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)), 0 => Vcat(a, Fill(0,∞)), 1 => Vcat([1/4], Fill(1/4,∞)))
83-
λ = 0.5+0.1im
84-
85-
B = BandedMatrix(A-λ*I, (3,1))
86-
T = toeptail(B)
87-
Q,L = ql(T)
88-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
89-
90-
Q,L = ql(A-λ*I)
91-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ (A-λ*I)[1:10,1:10]
92-
93-
a = [-0.1,0.2,0.3]
94-
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)), 0 => Vcat(a, Fill(0,∞)), 1 => Vcat([1/4], Fill(1/4,∞)))
95-
for λ in (0.1+0im,0.1-0im, 3.0, 3.0+1im, 3.0-im, -0.1+0im, -0.1-0im)
96-
Q, L = ql(A-λ*I)
97-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ (A-λ*I)[1:10,1:10]
98-
end
99-
end
100-
101-
@testset "Pert faux-periodic QL" begin
102-
a = [0.5794879759059747 + 0.0im,0.538107104952824 - 0.951620830938543im,-0.19352887774167749 - 0.3738926065520737im,0.4314153362874331,0.0]
103-
T = _BandedMatrix(a*Ones{ComplexF64}(1,∞), ∞, 3,1)
104-
Q,L = ql(T)
105-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
106-
Qn,Ln = ql(T[1:1001,1:1001])
107-
@test Qn[1:10,1:10] * diagm(0=>[-1; (-1).^(1:9)]) ≈ Q[1:10,1:10]
108-
@test diagm(0=>[-1; (-1).^(1:9)]) * Ln[1:10,1:10] ≈ L[1:10,1:10]
109-
110-
f = [0.0+0.0im 0.522787+0.0im ; 0.59647-1.05483im 0.538107-0.951621im; -0.193529-0.373893im -0.193529-0.373893im;
111-
0.431415+0.0im 0.431415+0.0im; 0.0+0.0im 0.0+0.0im]
112-
A = _BandedMatrix(Hcat(f, a*Ones{ComplexF64}(1,∞)), ∞, 3,1)
113-
Q,L = ql(A)
114-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
115-
Qn,Ln = ql(A[1:1000,1:1000])
116-
@test Qn[1:10,1:10] * diagm(0 => [Ones(5); -(-1).^(1:5)]) ≈ Q[1:10,1:10]
117-
@test diagm(0 => [Ones(5); -(-1).^(1:5)]) * Ln[1:10,1:10] ≈ L[1:10,1:10]
118-
end

‎test/test_hessenbergq.jl

Lines changed: 0 additions & 84 deletions
Original file line numberDiff line numberDiff line change
@@ -62,87 +62,3 @@ import InfiniteLinearAlgebra: UpperHessenbergQ, LowerHessenbergQ, tail_de, _Band
6262
end
6363
end
6464
end
65-
66-
67-
@testset "Toeplitz QLHessenberg" begin
68-
@testset "Derivation" begin
69-
A = BandedMatrix(-1 => Fill(2,∞), 0 => Fill(5,∞), 1 => Fill(0.5,∞))
70-
a = reverse(A.data.args[1])
71-
d,e = tail_de(a)
72-
X = [transpose(a); 0 d e]
73-
Q = LowerHessenbergQ(Fill(ql!(X).Q,∞))
74-
L = _BandedMatrix(Hcat([e; X[2,2]; X[2,1]], X[2,end:-1:1] * Ones{Float64}(1,∞)), ∞, 2, 0)
75-
Qn,Ln = ql(A[1:1000,1:1000])
76-
@test Q[1:10,1:10] ≈ Qn[1:10,1:10]
77-
@test (Q'A)[1:10,1:10] ≈ Ln[1:10,1:10] ≈ L[1:10,1:10]
78-
79-
A = BandedMatrix(-1 => Fill(2,∞), 0 => Fill(5+im,∞), 1 => Fill(0.5,∞))
80-
a = reverse(A.data.args[1])
81-
d,e = tail_de(a)
82-
X = [transpose(a); 0 d e]
83-
q = ql!(X).Q
84-
Q = LowerHessenbergQ(Fill(q,∞))
85-
L = _BandedMatrix(Hcat([e; X[2,2]; X[2,1]], X[2,end:-1:1] * Ones{Float64}(1,∞)), ∞, 2, 0)
86-
Qn,Ln = ql(A[1:1000,1:1000])
87-
@test Q[1:10,1:10] ≈ Qn[1:10,1:10] * diagm(0 => [1; -Ones(9)] )
88-
@test (Q'A)[1:10,1:10] ≈ diagm(0 => [1; -Ones(9)] ) * Ln[1:10,1:10] ≈ L[1:10,1:10]
89-
end
90-
91-
@testset "Tridiagonal Toeplitz" begin
92-
for (Z,A,B) in ((2.0,5.1,0.5), (2.0,2.2,0.5), (2.0,-2.2,0.5), (2.0,-5.1,0.5),
93-
(0.5,5.1,2.0), (0.5,-5.1,2.0))
94-
n = 100_000; T = Tridiagonal(Fill(Z,∞), Fill(A,∞), Fill(B,∞)); Qn,Ln = ql(BandedMatrix(T)[1:n,1:n]);
95-
Q,L = ql(T)
96-
@test L[1:10,1:10] ≈ Ln[1:10,1:10]
97-
@test Q[1:10,1:10] ≈ Qn[1:10,1:10]
98-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
99-
end
100-
101-
for (Z,A,B) in ((2,5.0im,0.5),
102-
(2,2.1+0.1im,0.5),
103-
(2,2.1-0.1im,0.5),
104-
(2,1.5+0.1im,0.5),
105-
(2,1.5+0.0im,0.5),
106-
(2,1.5-0.0im,0.5),
107-
(2,0.0+0.0im,0.5),
108-
(2,-0.1+0.1im,0.5),
109-
(2,-1.1+0.1im,0.5),
110-
(2,-0.1-0.1im,0.5))
111-
T = Tridiagonal(Fill(ComplexF64(Z),∞), Fill(ComplexF64(A),∞), Fill(ComplexF64(B),∞))
112-
Q,L = ql(T)
113-
@test Q[1:10,1:12] * L[1:12,1:10] ≈ T[1:10,1:10]
114-
end
115-
116-
for (Z,A,B) in ((0.5,2.1,2.0),(0.5,-1.1,2.0))
117-
@test_throws DomainError ql(Tridiagonal(Fill(Z,∞), Fill(A,∞), Fill(B,∞)))
118-
end
119-
end
120-
121-
@testset "Hessenberg Toeplitz" begin
122-
a = [1,2,3,0.5]
123-
T = _BandedMatrix(reverse(a) * Ones(1,∞), ∞, 2, 1)
124-
F = ql(T)
125-
@test F.Q[1:10,1:11]*F.L[1:11,1:10] ≈ T[1:10,1:10]
126-
127-
a = [1,2,3+im,0.5]
128-
T = _BandedMatrix(reverse(a) * Ones{eltype(a)}(1,∞), ∞, 2, 1)
129-
Q,L = ql(T)
130-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
131-
@test T isa InfToeplitz
132-
133-
T = BandedMatrix(-2 => Fill(1,∞), 0 => Fill(0.5+eps()im,∞), 1 => Fill(0.25,∞))
134-
Q,L = ql(T)
135-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
136-
137-
for T in (BandedMatrix(-2 => Fill(1,∞), 0 => Fill(0.5,∞), 1 => Fill(0.25,∞)),
138-
BandedMatrix(-2 => Fill(1/4,∞), 1 => Fill(1,∞))-im*I)
139-
Q,L = ql(T)
140-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
141-
end
142-
143-
a = [ -2.531640004434771-0.0im , 0.36995310821558014+2.5612894011525276im, -0.22944284364953327+0.39386202384951985im, -0.2700241133710857 + 0.8984628598798804im, 4.930380657631324e-32 + 0.553001215633963im ]
144-
T = _BandedMatrix(a * Ones{ComplexF64}(1,∞), ∞, 3, 1)
145-
Q,L = ql(T)
146-
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
147-
end
148-
end

‎test/test_infql.jl

Lines changed: 158 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,158 @@
1+
2+
@testset "Inf QL" begin
3+
@testset "Toeplitz QLHessenberg" begin
4+
@testset "Derivation" begin
5+
A = BandedMatrix(-1 => Fill(2,∞), 0 => Fill(5,∞), 1 => Fill(0.5,∞))
6+
a = reverse(A.data.args[1])
7+
d,e = tail_de(a)
8+
X = [transpose(a); 0 d e]
9+
Q = LowerHessenbergQ(Fill(ql!(X).Q,∞))
10+
L = _BandedMatrix(Hcat([e; X[2,2]; X[2,1]], X[2,end:-1:1] * Ones{Float64}(1,∞)), ∞, 2, 0)
11+
Qn,Ln = ql(A[1:1000,1:1000])
12+
@test Q[1:10,1:10] ≈ Qn[1:10,1:10]
13+
@test (Q'A)[1:10,1:10] ≈ Ln[1:10,1:10] ≈ L[1:10,1:10]
14+
15+
A = BandedMatrix(-1 => Fill(2,∞), 0 => Fill(5+im,∞), 1 => Fill(0.5,∞))
16+
a = reverse(A.data.args[1])
17+
d,e = tail_de(a)
18+
X = [transpose(a); 0 d e]
19+
q = ql!(X).Q
20+
Q = LowerHessenbergQ(Fill(q,∞))
21+
L = _BandedMatrix(Hcat([e; X[2,2]; X[2,1]], X[2,end:-1:1] * Ones{Float64}(1,∞)), ∞, 2, 0)
22+
Qn,Ln = ql(A[1:1000,1:1000])
23+
@test Q[1:10,1:10] ≈ Qn[1:10,1:10] * diagm(0 => [1; -Ones(9)] )
24+
@test (Q'A)[1:10,1:10] ≈ diagm(0 => [1; -Ones(9)] ) * Ln[1:10,1:10] ≈ L[1:10,1:10]
25+
end
26+
27+
@testset "Tridiagonal Toeplitz" begin
28+
for (Z,A,B) in ((2.0,5.1,0.5), (2.0,2.2,0.5), (2.0,-2.2,0.5), (2.0,-5.1,0.5),
29+
(0.5,5.1,2.0), (0.5,-5.1,2.0))
30+
n = 100_000; T = Tridiagonal(Fill(Z,∞), Fill(A,∞), Fill(B,∞)); Qn,Ln = ql(BandedMatrix(T)[1:n,1:n]);
31+
Q,L = ql(T)
32+
@test L[1:10,1:10] ≈ Ln[1:10,1:10]
33+
@test Q[1:10,1:10] ≈ Qn[1:10,1:10]
34+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
35+
end
36+
37+
for (Z,A,B) in ((2,5.0im,0.5),
38+
(2,2.1+0.1im,0.5),
39+
(2,2.1-0.1im,0.5),
40+
(2,1.5+0.1im,0.5),
41+
(2,1.5+0.0im,0.5),
42+
(2,1.5-0.0im,0.5),
43+
(2,0.0+0.0im,0.5),
44+
(2,-0.1+0.1im,0.5),
45+
(2,-1.1+0.1im,0.5),
46+
(2,-0.1-0.1im,0.5))
47+
T = Tridiagonal(Fill(ComplexF64(Z),∞), Fill(ComplexF64(A),∞), Fill(ComplexF64(B),∞))
48+
Q,L = ql(T)
49+
@test Q[1:10,1:12] * L[1:12,1:10] ≈ T[1:10,1:10]
50+
end
51+
52+
for (Z,A,B) in ((0.5,2.1,2.0),(0.5,-1.1,2.0))
53+
@test_throws DomainError ql(Tridiagonal(Fill(Z,∞), Fill(A,∞), Fill(B,∞)))
54+
end
55+
end
56+
57+
@testset "Hessenberg Toeplitz" begin
58+
a = [1,2,3,0.5]
59+
T = _BandedMatrix(reverse(a) * Ones(1,∞), ∞, 2, 1)
60+
F = ql(T)
61+
@test F.Q[1:10,1:11]*F.L[1:11,1:10] ≈ T[1:10,1:10]
62+
63+
a = [1,2,3+im,0.5]
64+
T = _BandedMatrix(reverse(a) * Ones{eltype(a)}(1,∞), ∞, 2, 1)
65+
Q,L = ql(T)
66+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
67+
@test T isa InfToeplitz
68+
69+
T = BandedMatrix(-2 => Fill(1,∞), 0 => Fill(0.5+eps()im,∞), 1 => Fill(0.25,∞))
70+
Q,L = ql(T)
71+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
72+
73+
for T in (BandedMatrix(-2 => Fill(1,∞), 0 => Fill(0.5,∞), 1 => Fill(0.25,∞)),
74+
BandedMatrix(-2 => Fill(1/4,∞), 1 => Fill(1,∞))-im*I)
75+
Q,L = ql(T)
76+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
77+
end
78+
79+
a = [ -2.531640004434771-0.0im , 0.36995310821558014+2.5612894011525276im, -0.22944284364953327+0.39386202384951985im, -0.2700241133710857 + 0.8984628598798804im, 4.930380657631324e-32 + 0.553001215633963im ]
80+
T = _BandedMatrix(a * Ones{ComplexF64}(1,∞), ∞, 3, 1)
81+
Q,L = ql(T)
82+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
83+
end
84+
end
85+
86+
@testset "PertTriToeplitz QL" begin
87+
A = Tridiagonal(Vcat(Float64[], Fill(2.0,∞)),
88+
Vcat(Float64[2.0], Fill(0.0,∞)),
89+
Vcat(Float64[], Fill(0.5,∞)))
90+
for λ in (-2.1-0.01im,-2.1+0.01im,-2.1+eps()im,-2.1-eps()im,-2.1+0.0im,-2.1-0.0im,-1.0+im,-3.0+im,-3.0-im)
91+
Q, L = ql(A - λ*I)
92+
@test Q[1:10,1:12]*L[1:12,1:10] ≈ A[1:10,1:10] - λ*I
93+
end
94+
end
95+
96+
@testset "Pert Hessenberg Toeplitz" begin
97+
a = [1,2,5,0.5]
98+
Random.seed!(0)
99+
A = _BandedMatrix(Hcat(randn(4,2), reverse(a) * Ones(1,∞)), ∞, 2, 1)
100+
@test A isa PertToeplitz
101+
@test BandedMatrix(A, (3,1))[1:10,1:10] == A[1:10,1:10]
102+
Q,L = ql(A)
103+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
104+
105+
a = [0.1,1,2,3,0.5]
106+
A = _BandedMatrix(Hcat([0.5 0.5; -1 3; 2 2; 1 1; 0.1 0.1], reverse(a) * Ones(1,∞)), ∞, 3, 1)
107+
@test A isa PertToeplitz
108+
@test BandedMatrix(A, (3,1))[1:10,1:10] == A[1:10,1:10]
109+
110+
B = BandedMatrix(A, (3,1))
111+
@test B[1:10,1:10] == A[1:10,1:10]
112+
Q,L = ql(A)
113+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
114+
115+
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)),
116+
0 => Vcat([0.0], Fill(1/2,∞)),
117+
1 => Vcat([1/4], Fill(1/4,∞)))
118+
Q,L = ql(A)
119+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
120+
a = [-0.1,0.2,0.3]
121+
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)), 0 => Vcat(a, Fill(0,∞)), 1 => Vcat([1/4], Fill(1/4,∞)))
122+
λ = 0.5+0.1im
123+
124+
B = BandedMatrix(A-λ*I, (3,1))
125+
T = toeptail(B)
126+
Q,L = ql(T)
127+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
128+
129+
Q,L = ql(A-λ*I)
130+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ (A-λ*I)[1:10,1:10]
131+
132+
a = [-0.1,0.2,0.3]
133+
A = BandedMatrix(-2 => Vcat([1], Fill(1,∞)), 0 => Vcat(a, Fill(0,∞)), 1 => Vcat([1/4], Fill(1/4,∞)))
134+
for λ in (0.1+0im,0.1-0im, 3.0, 3.0+1im, 3.0-im, -0.1+0im, -0.1-0im)
135+
Q, L = ql(A-λ*I)
136+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ (A-λ*I)[1:10,1:10]
137+
end
138+
end
139+
140+
@testset "Pert faux-periodic QL" begin
141+
a = [0.5794879759059747 + 0.0im,0.538107104952824 - 0.951620830938543im,-0.19352887774167749 - 0.3738926065520737im,0.4314153362874331,0.0]
142+
T = _BandedMatrix(a*Ones{ComplexF64}(1,∞), ∞, 3,1)
143+
Q,L = ql(T)
144+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ T[1:10,1:10]
145+
Qn,Ln = ql(T[1:1001,1:1001])
146+
@test Qn[1:10,1:10] * diagm(0=>[-1; (-1).^(1:9)]) ≈ Q[1:10,1:10]
147+
@test diagm(0=>[-1; (-1).^(1:9)]) * Ln[1:10,1:10] ≈ L[1:10,1:10]
148+
149+
f = [0.0+0.0im 0.522787+0.0im ; 0.59647-1.05483im 0.538107-0.951621im; -0.193529-0.373893im -0.193529-0.373893im;
150+
0.431415+0.0im 0.431415+0.0im; 0.0+0.0im 0.0+0.0im]
151+
A = _BandedMatrix(Hcat(f, a*Ones{ComplexF64}(1,∞)), ∞, 3,1)
152+
Q,L = ql(A)
153+
@test Q[1:10,1:11]*L[1:11,1:10] ≈ A[1:10,1:10]
154+
Qn,Ln = ql(A[1:1000,1:1000])
155+
@test Qn[1:10,1:10] * diagm(0 => [Ones(5); -(-1).^(1:5)]) ≈ Q[1:10,1:10]
156+
@test diagm(0 => [Ones(5); -(-1).^(1:5)]) * Ln[1:10,1:10] ≈ L[1:10,1:10]
157+
end
158+
end

‎test/test_infqr.jl

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,9 +1,13 @@
11
using InfiniteLinearAlgebra, LinearAlgebra, BandedMatrices, InfiniteArrays, FillArrays
2-
import BandedMatrices: _BandedMatrix
2+
import BandedMatrices: _BandedMatrix, colsupport
33

44

55
A = _BandedMatrix(Vcat(Ones(1,∞), (1:∞)', Ones(1,∞)), ∞, 1, 1)
66
C = cache(A)
7+
V = view(C.data,1:10,1:11)
8+
V isa BandedMatrices.BandedSubBandedMatrix
9+
qr!(V)
10+
qr!(C.data)
711

812
F = QR(
913

0 commit comments

Comments
 (0)