Skip to content

Commit 453e9c6

Browse files
committed
Added leg2cheb, leg2chebu
1 parent 5a71929 commit 453e9c6

File tree

4 files changed

+50
-25
lines changed

4 files changed

+50
-25
lines changed

REQUIRE

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,2 +1,2 @@
11
julia 0.4
2-
ToeplitzMatrices
2+
ToeplitzMatrices 0.1

src/FastTransforms.jl

Lines changed: 10 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -6,6 +6,7 @@ using Base, ToeplitzMatrices
66
import Base: *
77

88
export cjt, icjt, jjt, plan_cjt, plan_icjt
9+
export leg2cheb,leg2chebu
910
export gaunt
1011

1112
# Other module methods and constants:
@@ -38,4 +39,13 @@ include("cjt.jl")
3839

3940
include("gaunt.jl")
4041

42+
43+
include("toeplitzhankel.jl")
44+
45+
leg2cheb(x...)=th_leg2cheb(x...)
46+
leg2chebu(x...)=th_leg2chebu(x...)
47+
48+
49+
50+
4151
end # module

src/toeplitzhankel.jl

Lines changed: 34 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -41,17 +41,26 @@ end
4141

4242

4343
# Plan a multiply by DL*(T.*H)*DR
44-
immutable ToepltizHankelPlan{T}
45-
T::TriangularToeplitz{T}
46-
C::Matrix{T} # A cholesky factorization of H: H=CC'
47-
DL::Vector{T}
48-
DR::Vector{T}
44+
immutable ToeplitzHankelPlan{TT}
45+
T::TriangularToeplitz{TT}
46+
C::Matrix{TT} # A cholesky factorization of H: H=CC'
47+
DL::Vector{TT}
48+
DR::Vector{TT}
49+
50+
ToeplitzHankelPlan(T,C,DL,DR)=new(T,C,DL,DR)
4951
end
5052

51-
ToepltizHankelPlan(T::TriangularToeplitz,C::Matrix)=ToepltizHankelPlan(T,C,ones(size(T,1)),ones(size(T,2)))
52-
ToepltizHankelPlan(T::TriangularToeplitz,H::Hankel,D...)=ToepltizHankelPlan(T,partialchol(H),D...)
5353

54-
*(P::ToepltizHankelPlan,v::Vector)=P.DL.*toeplitzcholmult(P.T,P.C,P.DR.*v)
54+
function ToeplitzHankelPlan(T::TriangularToeplitz,C::Matrix,DL::AbstractVector,DR::AbstractVector)
55+
TT=promote_type(eltype(T),eltype(C),eltype(DL),eltype(DR))
56+
ToeplitzHankelPlan{TT}(T,C,collect(TT,DL),collect(TT,DR))
57+
end
58+
ToeplitzHankelPlan(T::TriangularToeplitz,C::Matrix) =
59+
ToeplitzHankelPlan(T,C,ones(size(T,1)),ones(size(T,2)))
60+
ToeplitzHankelPlan(T::TriangularToeplitz,H::Hankel,D...) =
61+
ToeplitzHankelPlan(T,partialchol(H),D...)
62+
63+
*(P::ToeplitzHankelPlan,v::Vector)=P.DL.*toeplitzcholmult(P.T,P.C,P.DR.*v)
5564

5665

5766

@@ -87,16 +96,16 @@ end
8796
Λ(z::AbstractArray)=Λ(eltype(z),z)
8897

8998

90-
function leg2chebuTH(n)
91-
λ=Λ(0:0.5:n-1)
92-
t=zeros(n)
99+
function leg2chebuTH{TT}(::Type{TT},n)
100+
λ=Λ(TT,0:0.5:n-1)
101+
t=zeros(TT,n)
93102
t[1:2:end]=λ[1:2:n]./(((1:2:n)-2))
94103
T=TriangularToeplitz(-2/π*t,:U)
95104
H=Hankel(λ[1:n]./((1:n)+1),λ[n:end]./((n:2n-1)+1))
96105
T,H
97106
end
98107

99-
leg2chebuplan(n)=ToepltizHankelPlan(leg2chebuTH(n)...,1:n,ones(n))
108+
th_leg2chebuplan{TT}(::Type{TT},n)=ToeplitzHankelPlan(leg2chebuTH(TT,n)...,1:n,ones(TT,n))
100109

101110
function leg2chebTH{TT}(::Type{TT},n)
102111
λ=Λ(TT,0:0.5:n-1)
@@ -109,18 +118,19 @@ function leg2chebTH{TT}(::Type{TT},n)
109118
T,H,DL
110119
end
111120

112-
leg2chebplan{TT}(::Type{TT},n)=ToepltizHankelPlan(leg2chebTH(TT,n)...,ones(TT,n))
121+
th_leg2chebplan{TT}(::Type{TT},n)=ToeplitzHankelPlan(leg2chebTH(TT,n)...,ones(TT,n))
113122

114-
function leg2chebTHslow(n)
115-
λ=map(Λ,0:0.5:n-1)
116-
t=zeros(n)
117-
t[1:2:end]=λ[1:2:n]
118-
T=TriangularToeplitz(2/π*t,:U)
119-
H=Hankel(λ[1:n],λ[n:end])
120-
DL=ones(n)
121-
DL[1]*=0.5
122-
T,H,DL
123-
end
123+
# function leg2chebTHslow(n)
124+
# λ=map(Λ,0:0.5:n-1)
125+
# t=zeros(n)
126+
# t[1:2:end]=λ[1:2:n]
127+
# T=TriangularToeplitz(2/π*t,:U)
128+
# H=Hankel(λ[1:n],λ[n:end])
129+
# DL=ones(n)
130+
# DL[1]*=0.5
131+
# T,H,DL
132+
# end
124133

125134

126-
leg2cheb(v)=leg2chebplan(eltype(v),length(v))*v
135+
th_leg2cheb(v)=th_leg2chebplan(eltype(v),length(v))*v
136+
th_leg2chebu(v)=th_leg2chebuplan(eltype(v),length(v))*v

test/runtests.jl

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -175,3 +175,8 @@ include("gaunttest.jl")
175175
println("Testing BigFloat support of FFT and DCT methods")
176176

177177
include("fftBigFloattest.jl")
178+
179+
r=rand(1000)
180+
@test_approx_eq leg2cheb(r) cjt(r,0.,0.)
181+
182+
@test_approx_eq leg2chebu([1.0,2,3,4,5]) [0.546875,0.5,0.5390625,1.25,1.3671875]

0 commit comments

Comments
 (0)