Skip to content

Commit d46ad69

Browse files
MikaelSlevinskyMikaelSlevinsky
committed
add sph2fourier and fourier2sph API
re-work readme
1 parent 4f39dbf commit d46ad69

File tree

11 files changed

+197
-88
lines changed

11 files changed

+197
-88
lines changed

‎README.md

Lines changed: 92 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -2,24 +2,56 @@
22

33
[![Build Status](https://travis-ci.org/MikaelSlevinsky/FastTransforms.jl.svg?branch=master)](https://travis-ci.org/MikaelSlevinsky/FastTransforms.jl) [![](https://img.shields.io/badge/docs-stable-blue.svg)](https://MikaelSlevinsky.github.io/FastTransforms.jl/stable) [![](https://img.shields.io/badge/docs-latest-blue.svg)](https://MikaelSlevinsky.github.io/FastTransforms.jl/latest)
44

5-
The aim of this package is to provide new classes of fast transforms with low
6-
pre-computation. One approach is based on the use of asymptotic formulae to
7-
relate the transforms to a small number of fast Fourier transforms. Another
8-
approach is based on a Toeplitz-dot-Hankel decomposition of the matrix of
9-
connection coefficients. Both new classes of fast transforms do not
10-
require large pre-computation for fast execution and they are designed
11-
to work on expansions of functions with any degree of regularity.
12-
13-
The Chebyshev—Jacobi transform and its inverse are implemented. This
14-
allows the fast conversion of Chebyshev expansion coefficients to Jacobi expansion
15-
coefficients and back.
5+
The aim of this package is to provide fast orthogonal polynomial transforms that are designed for expansions of functions with any degree of regularity. There are multiple approaches to the classical connection problem, though the user does not need to know the specifics.
6+
7+
One approach is based on the use of asymptotic formulae to relate the transforms to a small number of fast Fourier transforms. Another approach is based on a Toeplitz-dot-Hankel decomposition of the matrix of connection coefficients. Yet another approach uses a hierarchical decomposition of the matrix of connection coefficients.
8+
9+
## The Chebyshev—Legendre Transform
10+
11+
The Chebyshev—Legendre transform allows the fast conversion of Chebyshev expansion coefficients to Legendre expansion coefficients and back.
12+
1613
```julia
1714
julia> Pkg.add("FastTransforms")
1815

1916
julia> using FastTransforms
2017

2118
julia> c = rand(10001);
2219

20+
julia> leg2cheb(c);
21+
22+
julia> cheb2leg(c);
23+
24+
julia> norm(cheb2leg(leg2cheb(c))-c)
25+
5.564168202018823e-13
26+
```
27+
28+
The implementation separates pre-computation into a type of plan. This type is constructed with either `plan_leg2cheb` or `plan_cheb2leg`. Let's see how much faster it is if we pre-compute.
29+
30+
```julia
31+
julia> p1 = plan_leg2cheb(c);
32+
33+
julia> p2 = plan_cheb2leg(c);
34+
35+
julia> @time leg2cheb(c);
36+
0.082615 seconds (11.94 k allocations: 31.214 MiB, 6.75% gc time)
37+
38+
julia> @time p1*c;
39+
0.004297 seconds (6 allocations: 78.422 KiB)
40+
41+
julia> @time cheb2leg(c);
42+
0.110388 seconds (11.94 k allocations: 31.214 MiB, 8.16% gc time)
43+
44+
julia> @time p2*c;
45+
0.004500 seconds (6 allocations: 78.422 KiB)
46+
```
47+
48+
## The Chebyshev—Jacobi Transform
49+
50+
The Chebyshev—Jacobi transform allows the fast conversion of Chebyshev expansion coefficients to Jacobi expansion coefficients and back.
51+
52+
```julia
53+
julia> c = rand(10001);
54+
2355
julia> @time norm(icjt(cjt(c,0.1,-0.2),0.1,-0.2)-c,Inf)
2456
0.258390 seconds (431 allocations: 6.278 MB)
2557
1.4830359162942841e-12
@@ -43,17 +75,58 @@ is valid for the half-open square `(α,β) ∈ (-1/2,1/2]^2`. Therefore, the fas
4375
when the parameters are inside. If the parameters `(α,β)` are not exceptionally beyond the square,
4476
then increment/decrement operators are used with linear complexity (and linear conditioning) in the degree.
4577

46-
The Padua transform and its inverse are also implemented thanks to
47-
[Michael Clarke](https://github.com/MikeAClarke). These are optimized methods
48-
designed for computing the bivariate Chebyshev coefficients by interpolating a
49-
bivariate function at the Padua points on `[-1,1]^2`.
78+
## The Padua Transform
79+
80+
The Padua transform and its inverse are implemented thanks to [Michael Clarke](https://github.com/MikeAClarke). These are optimized methods designed for computing the bivariate Chebyshev coefficients by interpolating a bivariate function at the Padua points on `[-1,1]^2`.
81+
82+
```julia
83+
julia> n=200;
84+
85+
julia> N=div((n+1)*(n+2),2)
86+
87+
julia> v=rand(N); #Length of v is the no. of Padua points
88+
89+
julia> @time norm(ipaduatransform(paduatransform(v))-v)
90+
0.006571 seconds (846 allocations: 1.746 MiB)
91+
3.123637691861415e-14
92+
93+
```
94+
95+
## The Spherical Harmonic Transform
96+
97+
Let `A` be a matrix of spherical harmonic expansion coefficients arranged by increasing alternating order. Then `sph2fourier` converts the representation into a bivariate Fourier series, and `fourier2sph` converts it back.
98+
```julia
99+
julia> A = rand(Float64, 251, 501); FastTransforms.zero_spurious_modes!(A);
100+
101+
julia> B = sph2fourier(A);
102+
103+
julia> C = fourier2sph(B);
104+
105+
julia> norm(A-C)
106+
107+
julia> A = rand(Float64, 1024, 2047); FastTransforms.zero_spurious_modes!(A);
108+
109+
julia> B = sph2fourier(A; sketch = :none);
110+
Pre-computing thin plan...100%|██████████████████████████████████████████████████| Time: 0:00:04
111+
112+
julia> C = fourier2sph(B; sketch = :none);
113+
Pre-computing thin plan...100%|██████████████████████████████████████████████████| Time: 0:00:04
114+
115+
julia> norm(A-C)
116+
1.5062262753260893e-12
117+
118+
```
119+
120+
As with other fast transforms, `plan_sph2fourier` saves effort by caching the pre-computation. Be warned that for dimensions larger than `1000`, this is no small feat!
50121

51122
# References:
52123

53-
1. B. Alpert and V. Rokhlin. <a href="http://dx.doi.org/10.1137/0912009">A fast algorithm for the evaluation of Legendre expansions</a>, *SIAM J. Sci. Stat. Comput.*, **12**:158—179, 1991.
124+
[1] B. Alpert and V. Rokhlin. <a href="http://dx.doi.org/10.1137/0912009">A fast algorithm for the evaluation of Legendre expansions</a>, *SIAM J. Sci. Stat. Comput.*, **12**:158—179, 1991.
125+
126+
[2] N. Hale and A. Townsend. <a href="http://dx.doi.org/10.1137/130932223">A fast, simple, and stable Chebyshev—Legendre transform using and asymptotic formula</a>, *SIAM J. Sci. Comput.*, **36**:A148—A167, 2014.
54127

55-
2. N. Hale and A. Townsend. <a href="http://dx.doi.org/10.1137/130932223">A fast, simple, and stable Chebyshev—Legendre transform using and asymptotic formula</a>, *SIAM J. Sci. Comput.*, **36**:A148—A167, 2014.
128+
[3] R. M. Slevinsky. <a href="https://doi.org/10.1093/imanum/drw070">On the use of Hahn's asymptotic formula and stabilized recurrence for a fast, simple, and stable Chebyshev—Jacobi transform</a>, published online in *IMA J. Numer. Anal.*, 2017.
56129

57-
3. R. M. Slevinsky. <a href="https://doi.org/10.1093/imanum/drw070">On the use of Hahn's asymptotic formula and stabilized recurrence for a fast, simple, and stable Chebyshev—Jacobi transform</a>, published online in *IMA J. Numer. Anal.*, 2017.
130+
[4] R. M. Slevinsky. <a href="https://arxiv.org/abs/">Fast and backward stable transforms between spherical harmonic expansions and bivariate Fourier series</a>, arXiv, 2017.
58131

59-
4. A. Townsend, M. Webb, and S. Olver. <a href="http://arxiv.org/abs/1604.07486">Fast polynomial transforms based on Toeplitz and Hankel matrices</a>, arXiv:1604.07486, 2016.
132+
[5] A. Townsend, M. Webb, and S. Olver. <a href="http://arxiv.org/abs/1604.07486">Fast polynomial transforms based on Toeplitz and Hankel matrices</a>, arXiv:1604.07486, 2016.

‎REQUIRE

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,6 @@
11
julia 0.5
22
ToeplitzMatrices 0.2
3+
HierarchicalMatrices 0.0.1
34
LowRankApprox 0.0.2
45
ProgressMeter 0.3.4
56
Compat 0.18

‎src/FastTransforms.jl

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -3,11 +3,11 @@ module FastTransforms
33

44
using Base, ToeplitzMatrices, HierarchicalMatrices, LowRankApprox, ProgressMeter, Compat
55

6-
import Base: *, size, view, A_mul_B!, At_mul_B!, Ac_mul_B!
6+
import Base: *, \, size, view, A_mul_B!, At_mul_B!, Ac_mul_B!
77
import Base: getindex, setindex!, Factorization, length
88
import Base.LinAlg: BlasFloat, BlasInt
99
import HierarchicalMatrices: HierarchicalMatrix, unsafe_broadcasttimes!
10-
import LowRankApprox: rowperm!, ColPerm
10+
import LowRankApprox: ColPerm
1111

1212
export cjt, icjt, jjt, plan_cjt, plan_icjt
1313
export leg2cheb, cheb2leg, leg2chebu, ultra2ultra, jac2jac
@@ -20,6 +20,7 @@ export paduatransform, ipaduatransform, paduatransform!, ipaduatransform!, padua
2020
export plan_paduatransform!, plan_ipaduatransform!
2121

2222
export SlowSphericalHarmonicPlan, FastSphericalHarmonicPlan, ThinSphericalHarmonicPlan
23+
export sph2fourier, fourier2sph, plan_sph2fourier
2324

2425
# Other module methods and constants:
2526
#export ChebyshevJacobiPlan, jac2cheb, cheb2jac

‎src/SphericalHarmonics/SphericalHarmonics.jl

Lines changed: 22 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,7 +10,29 @@ function zero_spurious_modes!(A::AbstractMatrix)
1010
A
1111
end
1212

13+
@compat abstract type SphericalHarmonicPlan{T} end
14+
15+
function *(P::SphericalHarmonicPlan, X::AbstractMatrix)
16+
A_mul_B!(zero(X), P, X)
17+
end
18+
19+
function \(P::SphericalHarmonicPlan, X::AbstractMatrix)
20+
At_mul_B!(zero(X), P, X)
21+
end
22+
1323
include("slowplan.jl")
1424
include("Butterfly.jl")
1525
include("fastplan.jl")
1626
include("thinplan.jl")
27+
28+
function plan_sph2fourier(A::AbstractMatrix; opts...)
29+
M, N = size(A)
30+
if M ≤ 1022
31+
SlowSphericalHarmonicPlan(A)
32+
else
33+
ThinSphericalHarmonicPlan(A; opts...)
34+
end
35+
end
36+
37+
sph2fourier(A::AbstractMatrix; opts...) = plan_sph2fourier(A; opts...)*A
38+
fourier2sph(A::AbstractMatrix; opts...) = plan_sph2fourier(A; opts...)\A

‎src/SphericalHarmonics/fastplan.jl

Lines changed: 5 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -1,4 +1,4 @@
1-
immutable FastSphericalHarmonicPlan{T}
1+
immutable FastSphericalHarmonicPlan{T} <: SphericalHarmonicPlan{T}
22
RP::RotationPlan{T}
33
BF::Vector{Butterfly{T}}
44
p1::NormalizedLegendreToChebyshevPlan{T}
@@ -14,11 +14,10 @@ function FastSphericalHarmonicPlan{T}(A::Matrix{T}, L::Int; opts...)
1414
@assert N == 2M-1
1515
n = (N+1)÷2
1616
RP = RotationPlan(T, n-1)
17-
a1 = A[:,1]
18-
p1 = plan_normleg2cheb(a1)
19-
p2 = plan_normleg12cheb2(a1)
20-
p1inv = plan_cheb2normleg(a1)
21-
p2inv = plan_cheb22normleg1(a1)
17+
p1 = plan_normleg2cheb(A)
18+
p2 = plan_normleg12cheb2(A)
19+
p1inv = plan_cheb2normleg(A)
20+
p2inv = plan_cheb22normleg1(A)
2221
B = zeros(A)
2322
Ce = eye(T, M)
2423
Co = eye(T, M)

‎src/SphericalHarmonics/slowplan.jl

Lines changed: 5 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -124,7 +124,7 @@ end
124124
Ac_mul_B!(P::RotationPlan, A::AbstractMatrix) = At_mul_B!(P, A)
125125

126126

127-
immutable SlowSphericalHarmonicPlan{T}
127+
immutable SlowSphericalHarmonicPlan{T} <: SphericalHarmonicPlan{T}
128128
RP::RotationPlan{T}
129129
p1::NormalizedLegendreToChebyshevPlan{T}
130130
p2::NormalizedLegendre1ToChebyshev2Plan{T}
@@ -137,11 +137,10 @@ function SlowSphericalHarmonicPlan{T}(A::Matrix{T})
137137
M, N = size(A)
138138
n = (N+1)÷2
139139
RP = RotationPlan(T, n-1)
140-
a1 = A[:,1]
141-
p1 = plan_normleg2cheb(a1)
142-
p2 = plan_normleg12cheb2(a1)
143-
p1inv = plan_cheb2normleg(a1)
144-
p2inv = plan_cheb22normleg1(a1)
140+
p1 = plan_normleg2cheb(A)
141+
p2 = plan_normleg12cheb2(A)
142+
p1inv = plan_cheb2normleg(A)
143+
p2inv = plan_cheb22normleg1(A)
145144
B = zeros(A)
146145
SlowSphericalHarmonicPlan(RP, p1, p2, p1inv, p2inv, B)
147146
end

‎src/SphericalHarmonics/thinplan.jl

Lines changed: 5 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@ const LAYERSKELETON = 64
22

33
checklayer(j::Int) = j÷LAYERSKELETON == j/LAYERSKELETON
44

5-
immutable ThinSphericalHarmonicPlan{T}
5+
immutable ThinSphericalHarmonicPlan{T} <: SphericalHarmonicPlan{T}
66
RP::RotationPlan{T}
77
BF::Vector{Butterfly{T}}
88
p1::NormalizedLegendreToChebyshevPlan{T}
@@ -18,11 +18,10 @@ function ThinSphericalHarmonicPlan{T}(A::Matrix{T}, L::Int; opts...)
1818
@assert N == 2M-1
1919
n = (N+1)÷2
2020
RP = RotationPlan(T, n-1)
21-
a1 = A[:,1]
22-
p1 = plan_normleg2cheb(a1)
23-
p2 = plan_normleg12cheb2(a1)
24-
p1inv = plan_cheb2normleg(a1)
25-
p2inv = plan_cheb22normleg1(a1)
21+
p1 = plan_normleg2cheb(A)
22+
p2 = plan_normleg12cheb2(A)
23+
p1inv = plan_cheb2normleg(A)
24+
p2inv = plan_cheb22normleg1(A)
2625
B = zeros(A)
2726
Ce = eye(T, M)
2827
Co = eye(T, M)

‎src/hierarchical.jl

Lines changed: 30 additions & 47 deletions
Original file line numberDiff line numberDiff line change
@@ -40,8 +40,8 @@ end
4040

4141
size(P::HierarchicalPlanWithParity) = (size(P.even, 1)+size(P.odd, 1), size(P.even, 2)+size(P.odd, 2))
4242

43-
evenlength(v::Vector) = (L = length(v); iseven(L) ? L÷2 : (L+1)÷2)
44-
oddlength(v::Vector) = (L = length(v); iseven(L) ? L÷2 : (L-1)÷2)
43+
evensize(v::AbstractVecOrMat, d) = (L = size(v, d); iseven(L) ? L÷2 : (L+1)÷2)
44+
oddsize(v::AbstractVecOrMat, d) = (L = size(v, d); iseven(L) ? L÷2 : (L-1)÷2)
4545

4646
UpperTriangularHierarchicalMatrix{T}(::Type{T}, f::Function, bd::Int64) = UpperTriangularHierarchicalMatrix(T, f, bd, bd)
4747
UpperTriangularHierarchicalMatrix{T}(::Type{T}, f::Function, b::Int64, d::Int64) = UpperTriangularHierarchicalMatrix(T, f, 1, b, 1, d)
@@ -154,21 +154,6 @@ function getindex(P::ChebyshevToLegendrePlan, i::Int, j::Int)
154154
end
155155
end
156156

157-
LegendreToChebyshevPlan(v::Vector) = LegendreToChebyshevPlan(plan_even_leg2cheb(v), plan_odd_leg2cheb(v))
158-
ChebyshevToLegendrePlan(v::Vector) = ChebyshevToLegendrePlan(plan_even_cheb2leg(v), plan_odd_cheb2leg(v))
159-
160-
plan_leg2cheb(v::Vector) = LegendreToChebyshevPlan(v)
161-
plan_cheb2leg(v::Vector) = ChebyshevToLegendrePlan(v)
162-
163-
plan_even_leg2cheb(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Meven, evenlength(v))
164-
plan_odd_leg2cheb(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Modd, oddlength(v))
165-
166-
plan_even_cheb2leg(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Leven, evenlength(v))
167-
plan_odd_cheb2leg(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Lodd, oddlength(v))
168-
169-
leg2cheb(v::Vector) = plan_leg2cheb(v)*v
170-
cheb2leg(v::Vector) = plan_cheb2leg(v)*v
171-
172157
function A_mul_B!(y::Vector, P::LegendreToChebyshevPlan, x::AbstractVector)
173158
HierarchicalMatrices.A_mul_B!(y, P.even, x, 1, 1, 2, 2)
174159
HierarchicalMatrices.A_mul_B!(y, P.odd, x, 2, 2, 2, 2)
@@ -271,21 +256,6 @@ function getindex(P::ChebyshevToNormalizedLegendrePlan, i::Int, j::Int)
271256
end
272257
end
273258

274-
NormalizedLegendreToChebyshevPlan(v::Vector) = NormalizedLegendreToChebyshevPlan(plan_even_normleg2cheb(v), plan_odd_normleg2cheb(v), eltype(v)[sqrt(j-0.5) for j in 1:length(v)])
275-
ChebyshevToNormalizedLegendrePlan(v::Vector) = ChebyshevToNormalizedLegendrePlan(plan_even_cheb2normleg(v), plan_odd_cheb2normleg(v), eltype(v)[sqrt(i-0.5) for i in 1:length(v)])
276-
277-
plan_normleg2cheb(v::Vector) = NormalizedLegendreToChebyshevPlan(v)
278-
plan_cheb2normleg(v::Vector) = ChebyshevToNormalizedLegendrePlan(v)
279-
280-
plan_even_normleg2cheb(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Mevennorm, evenlength(v))
281-
plan_odd_normleg2cheb(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Moddnorm, oddlength(v))
282-
283-
plan_even_cheb2normleg(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Levennorm, evenlength(v))
284-
plan_odd_cheb2normleg(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Loddnorm, oddlength(v))
285-
286-
normleg2cheb(v::Vector) = plan_normleg2cheb(v)*v
287-
cheb2normleg(v::Vector) = plan_cheb2normleg(v)*v
288-
289259
function A_mul_B!!(y::Vector, P::NormalizedLegendreToChebyshevPlan, x::AbstractVector)
290260
unsafe_broadcasttimes!(x, P.scl)
291261
HierarchicalMatrices.A_mul_B!(y, P.even, x, 1, 1, 2, 2)
@@ -394,21 +364,6 @@ function getindex(P::Chebyshev2ToNormalizedLegendre1Plan, i::Int, j::Int)
394364
end
395365
end
396366

397-
NormalizedLegendre1ToChebyshev2Plan(v::Vector) = NormalizedLegendre1ToChebyshev2Plan(plan_even_normleg12cheb2(v), plan_odd_normleg12cheb2(v), eltype(v)[sqrt((j+0.5)/(j*(j+1))) for j in 1:length(v)])
398-
Chebyshev2ToNormalizedLegendre1Plan(v::Vector) = Chebyshev2ToNormalizedLegendre1Plan(plan_even_cheb22normleg1(v), plan_odd_cheb22normleg1(v), eltype(v)[sqrt(i*(i+1)/(i+0.5)) for i in 1:length(v)])
399-
400-
plan_normleg12cheb2(v::Vector) = NormalizedLegendre1ToChebyshev2Plan(v)
401-
plan_cheb22normleg1(v::Vector) = Chebyshev2ToNormalizedLegendre1Plan(v)
402-
403-
plan_even_normleg12cheb2(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Mnormeven, evenlength(v))
404-
plan_odd_normleg12cheb2(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Mnormodd, oddlength(v))
405-
406-
plan_even_cheb22normleg1(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Lnormeven, evenlength(v))
407-
plan_odd_cheb22normleg1(v::Vector) = UpperTriangularHierarchicalMatrix(eltype(v), Lnormodd, oddlength(v))
408-
409-
normleg12cheb2(v::Vector) = plan_normleg12cheb2(v)*v
410-
cheb22normleg1(v::Vector) = plan_cheb22normleg1(v)*v
411-
412367
function A_mul_B!!(y::Vector, P::NormalizedLegendre1ToChebyshev2Plan, x::AbstractVector)
413368
unsafe_broadcasttimes!(x, P.scl)
414369
HierarchicalMatrices.A_mul_B!(y, P.even, x, 1, 1, 2, 2)
@@ -456,3 +411,31 @@ function A_mul_B_col_J!(Y::Matrix, P::Chebyshev2ToNormalizedLegendre1Plan, X::Ma
456411
HierarchicalMatrices.A_mul_B!(Y, P.odd, X, 2+COLSHIFT, 2+COLSHIFT, 2, 2)
457412
scale_col_J!(P.scl, Y, J)
458413
end
414+
415+
### API
416+
417+
for (f, T, feven, fodd) in ((:leg2cheb, :LegendreToChebyshevPlan, :Meven, :Modd),
418+
(:cheb2leg, :ChebyshevToLegendrePlan, :Leven, :Lodd),
419+
(:normleg2cheb, :NormalizedLegendreToChebyshevPlan, :Mevennorm, :Moddnorm),
420+
(:cheb2normleg, :ChebyshevToNormalizedLegendrePlan, :Levennorm, :Loddnorm),
421+
(:normleg12cheb2, :NormalizedLegendre1ToChebyshev2Plan, :Mnormeven, :Mnormodd),
422+
(:cheb22normleg1, :Chebyshev2ToNormalizedLegendre1Plan, :Lnormeven, :Lnormodd))
423+
plan_f = parse("plan_"*string(f))
424+
plan_even_f = parse("plan_even_"*string(f))
425+
plan_odd_f = parse("plan_odd_"*string(f))
426+
@eval begin
427+
$plan_f(v::VecOrMat) = $T(v)
428+
$f(v::VecOrMat) = $plan_f(v)*v
429+
$plan_even_f(v::VecOrMat) = UpperTriangularHierarchicalMatrix(eltype(v), $feven, evensize(v, 1))
430+
$plan_odd_f(v::VecOrMat) = UpperTriangularHierarchicalMatrix(eltype(v), $fodd, oddsize(v, 1))
431+
end
432+
end
433+
434+
LegendreToChebyshevPlan(v::VecOrMat) = LegendreToChebyshevPlan(plan_even_leg2cheb(v), plan_odd_leg2cheb(v))
435+
ChebyshevToLegendrePlan(v::VecOrMat) = ChebyshevToLegendrePlan(plan_even_cheb2leg(v), plan_odd_cheb2leg(v))
436+
437+
NormalizedLegendreToChebyshevPlan(v::VecOrMat) = NormalizedLegendreToChebyshevPlan(plan_even_normleg2cheb(v), plan_odd_normleg2cheb(v), eltype(v)[sqrt(j-0.5) for j in 1:size(v, 1)])
438+
ChebyshevToNormalizedLegendrePlan(v::VecOrMat) = ChebyshevToNormalizedLegendrePlan(plan_even_cheb2normleg(v), plan_odd_cheb2normleg(v), eltype(v)[sqrt(i-0.5) for i in 1:size(v, 1)])
439+
440+
NormalizedLegendre1ToChebyshev2Plan(v::VecOrMat) = NormalizedLegendre1ToChebyshev2Plan(plan_even_normleg12cheb2(v), plan_odd_normleg12cheb2(v), eltype(v)[sqrt((j+0.5)/(j*(j+1))) for j in 1:size(v, 1)])
441+
Chebyshev2ToNormalizedLegendre1Plan(v::VecOrMat) = Chebyshev2ToNormalizedLegendre1Plan(plan_even_cheb22normleg1(v), plan_odd_cheb22normleg1(v), eltype(v)[sqrt(i*(i+1)/(i+0.5)) for i in 1:size(v, 1)])

0 commit comments

Comments
 (0)