Skip to content

Commit cdaaf23

Browse files
committed
remove assert_applicable as FFTW does that for us
1 parent b25c099 commit cdaaf23

File tree

1 file changed

+22
-44
lines changed

1 file changed

+22
-44
lines changed

src/chebyshevtransform.jl

Lines changed: 22 additions & 44 deletions
Original file line numberDiff line numberDiff line change
@@ -6,14 +6,6 @@ size(P::ChebyshevPlan) = isdefined(P, :plan) ? size(P.plan) : (0,)
66
length(P::ChebyshevPlan) = isdefined(P, :plan) ? length(P.plan) : 0
77

88

9-
# Check whether a ChebyshevPlan is applicable to a given input array, and
10-
# throw an informative error if not:
11-
function assert_applicable(p::ChebyshevPlan{T}, X::StridedArray{T}) where T
12-
if size(X) != size(p)
13-
throw(ArgumentError("Chebyshev plan applied to wrong-size array"))
14-
end
15-
end
16-
179
const FIRSTKIND = 5
1810
const SECONDKIND = 3
1911

@@ -36,7 +28,7 @@ function plan_chebyshevtransform!(x::AbstractArray{T,N}, ::Val{1}, dims...; kws.
3628
ChebyshevTransformPlan{T,1,kindtuple(FIRSTKIND,N,dims...)}(FFTW.plan_r2r!(x, FFTW.REDFT10, dims...; kws...))
3729
end
3830
end
39-
function plan_chebyshevtransform!(x::AbstractArray{T}, ::Val{2}, dims...; kws...) where T<:fftwNumber
31+
function plan_chebyshevtransform!(x::AbstractArray{T,N}, ::Val{2}, dims...; kws...) where {T<:fftwNumber,N}
4032
length(x) 1 && throw(ArgumentError("Array must contain at least 2 entries"))
4133
ChebyshevTransformPlan{T,2,kindtuple(SECONDKIND,N,dims...)}(FFTW.plan_r2r!(x, FFTW.REDFT00, dims...; kws...))
4234
end
@@ -49,7 +41,7 @@ function plan_chebyshevtransform(x::AbstractArray{T,N}, ::Val{1}, dims...; kws..
4941
ChebyshevTransformPlan{T,1,kindtuple(FIRSTKIND,N,dims...)}(FFTW.plan_r2r(x, FFTW.REDFT10, dims...; kws...))
5042
end
5143
end
52-
function plan_chebyshevtransform(x::AbstractArray{T}, ::Val{2}, dims...; kws...) where T<:fftwNumber
44+
function plan_chebyshevtransform(x::AbstractArray{T,N}, ::Val{2}, dims...; kws...) where {T<:fftwNumber,N}
5345
length(x) 1 && throw(ArgumentError("Vector must contain at least 2 entries"))
5446
ChebyshevTransformPlan{T,2,kindtuple(SECONDKIND,N,dims...)}(FFTW.plan_r2r(x, FFTW.REDFT00, dims...; kws...))
5547
end
@@ -83,7 +75,6 @@ end
8375

8476
function *(P::ChebyshevTransformPlan{T,1,K,true}, x::AbstractArray{T}) where {T,K}
8577
n = length(x)
86-
assert_applicable(P, x)
8778
n == 0 && return x
8879

8980
y = P.plan*x # will be === x if in-place
@@ -93,7 +84,6 @@ end
9384
function mul!(y::AbstractArray{T}, P::ChebyshevTransformPlan{T,1,K,false}, x::AbstractArray) where {T,K}
9485
n = length(x)
9586
length(y) == n || throw(DimensionMismatch("output must match dimension"))
96-
assert_applicable(P, x)
9787
n == 0 && return y
9888
_plan_mul!(y, P.plan, x)
9989
_cheb1_rescale!(P.plan.region, y)
@@ -138,25 +128,25 @@ chebyshevtransform(x, dims...; kws...) = plan_chebyshevtransform(x, dims...; kws
138128
## Inverse transforms take Chebyshev coefficients and produce values at Chebyshev points of the first and second kinds
139129

140130

141-
const IFIRSTKIND = (4,)
131+
const IFIRSTKIND = 4
142132

143-
struct IChebyshevTransformPlan{T,kind,inplace,N,R} <: ChebyshevPlan{T}
144-
plan::FFTW.r2rFFTWPlan{T,kind,inplace,N,R}
145-
IChebyshevTransformPlan{T,kind,inplace,N,R}(plan) where {T,kind,inplace,N,R} = new{T,kind,inplace,N,R}(plan)
146-
IChebyshevTransformPlan{T,kind,inplace,N,R}() where {T,kind,inplace,N,R} = new{T,kind,inplace,N,R}()
133+
struct IChebyshevTransformPlan{T,kind,K,inplace,N,R} <: ChebyshevPlan{T}
134+
plan::FFTW.r2rFFTWPlan{T,K,inplace,N,R}
135+
IChebyshevTransformPlan{T,kind,K,inplace,N,R}(plan) where {T,kind,K,inplace,N,R} = new{T,kind,K,inplace,N,R}(plan)
136+
IChebyshevTransformPlan{T,kind,K,inplace,N,R}() where {T,kind,K,inplace,N,R} = new{T,kind,K,inplace,N,R}()
147137
end
148138

149-
IChebyshevTransformPlan{T,kind}(F::FFTW.r2rFFTWPlan{T,kind,inplace,N,R}) where {T,kind,inplace,N,R} =
150-
IChebyshevTransformPlan{T,kind,inplace,N,R}(F)
139+
IChebyshevTransformPlan{T,kind,K}(F::FFTW.r2rFFTWPlan{T,K,inplace,N,R}) where {T,kind,K,inplace,N,R} =
140+
IChebyshevTransformPlan{T,kind,K,inplace,N,R}(F)
151141

152142
size(P::IChebyshevTransformPlan) = isdefined(P, :plan) ? size(P.plan) : (0,)
153143
length(P::IChebyshevTransformPlan) = isdefined(P, :plan) ? length(P.plan) : 0
154144

155145

156146
# second kind Chebyshev transforms share a plan with their inverse
157147
# so we support this via inv
158-
inv(P::ChebyshevTransformPlan{T,SECONDKIND}) where {T} = IChebyshevTransformPlan{T,SECONDKIND}(P.plan)
159-
inv(P::IChebyshevTransformPlan{T,SECONDKIND}) where {T} = ChebyshevTransformPlan{T,SECONDKIND}(P.plan)
148+
inv(P::ChebyshevTransformPlan{T,2,K}) where {T,K} = IChebyshevTransformPlan{T,2,K}(P.plan)
149+
inv(P::IChebyshevTransformPlan{T,2,K}) where {T,K} = ChebyshevTransformPlan{T,2,K}(P.plan)
160150

161151

162152
\(P::ChebyshevTransformPlan, x::AbstractArray) = inv(P) * x
@@ -165,9 +155,9 @@ inv(P::IChebyshevTransformPlan{T,SECONDKIND}) where {T} = ChebyshevTransformPlan
165155

166156
function plan_ichebyshevtransform!(x::AbstractArray{T,N}, ::Val{1}, dims...; kws...) where {T<:fftwNumber,N}
167157
if isempty(x)
168-
IChebyshevTransformPlan{T,IFIRSTKIND,true,N,isempty(dims) ? UnitRange{Int} : typeof(dims)}()
158+
IChebyshevTransformPlan{T,1,kindtuple(IFIRSTKIND,N,dims...),true,N,isempty(dims) ? UnitRange{Int} : typeof(dims)}()
169159
else
170-
IChebyshevTransformPlan{T,IFIRSTKIND}(FFTW.plan_r2r!(x, FFTW.REDFT01, dims...; kws...))
160+
IChebyshevTransformPlan{T,1,kindtuple(IFIRSTKIND,N,dims...)}(FFTW.plan_r2r!(x, FFTW.REDFT01, dims...; kws...))
171161
end
172162
end
173163

@@ -177,9 +167,9 @@ end
177167

178168
function plan_ichebyshevtransform(x::AbstractArray{T,N}, ::Val{1}, dims...; kws...) where {T<:fftwNumber,N}
179169
if isempty(x)
180-
IChebyshevTransformPlan{T,IFIRSTKIND,false,N,isempty(dims) ? UnitRange{Int} : typeof(dims)}()
170+
IChebyshevTransformPlan{T,1,kindtuple(IFIRSTKIND,N,dims...),false,N,isempty(dims) ? UnitRange{Int} : typeof(dims)}()
181171
else
182-
IChebyshevTransformPlan{T,IFIRSTKIND}(FFTW.plan_r2r(x, FFTW.REDFT01, dims...; kws...))
172+
IChebyshevTransformPlan{T,1,kindtuple(IFIRSTKIND,N,dims...)}(FFTW.plan_r2r(x, FFTW.REDFT01, dims...; kws...))
183173
end
184174
end
185175

@@ -191,20 +181,18 @@ plan_ichebyshevtransform!(x::AbstractArray, dims...; kws...) = plan_ichebyshevtr
191181
plan_ichebyshevtransform(x::AbstractArray, dims...; kws...) = plan_ichebyshevtransform(x, Val(1), dims...; kws...)
192182

193183

194-
function *(P::IChebyshevTransformPlan{T,IFIRSTKIND,true}, x::AbstractVector{T}) where T<:fftwNumber
184+
function *(P::IChebyshevTransformPlan{T,1,K,true}, x::AbstractVector{T}) where {T<:fftwNumber,K}
195185
n = length(x)
196-
assert_applicable(P, x)
197186
n == 0 && return x
198187

199188
x[1] *= 2
200189
x = lmul!(one(T)/2, P.plan*x)
201190
x
202191
end
203192

204-
function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,IFIRSTKIND,false}, x::AbstractVector{T}) where T<:fftwNumber
193+
function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,1,K,false}, x::AbstractVector{T}) where {T<:fftwNumber,K}
205194
n = length(x)
206195
length(y) == n || throw(DimensionMismatch("output must match dimension"))
207-
assert_applicable(P, x)
208196
n == 0 && return y
209197

210198
x[1] *= 2 # Todo: don't mutate x
@@ -213,20 +201,18 @@ function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,IFIRSTKIND,fals
213201
lmul!(one(T)/2, y)
214202
end
215203

216-
function *(P::IChebyshevTransformPlan{T,SECONDKIND, true}, x::AbstractVector{T}) where T<:fftwNumber
204+
function *(P::IChebyshevTransformPlan{T,2,K, true}, x::AbstractVector{T}) where {T<:fftwNumber,K}
217205
n = length(x)
218-
assert_applicable(P, x)
219206

220207
x[1] *= 2; x[end] *= 2
221208
x = inv(P)*x
222209
x[1] *= 2; x[end] *= 2
223210
lmul!(convert(T,n-1)/2,x)
224211
end
225212

226-
function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,SECONDKIND,false}, x::AbstractVector{T}) where T<:fftwNumber
213+
function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,2,K,false}, x::AbstractVector{T}) where {T<:fftwNumber,K}
227214
n = length(x)
228215
length(y) == n || throw(DimensionMismatch("output must match dimension"))
229-
assert_applicable(P, x)
230216

231217
x[1] *= 2; x[end] *= 2
232218
_plan_mul!(y, inv(P), x)
@@ -235,13 +221,9 @@ function mul!(y::AbstractVector{T}, P::IChebyshevTransformPlan{T,SECONDKIND,fals
235221
lmul!(convert(T,(n-1))/2,y)
236222
end
237223

238-
*(P::IChebyshevTransformPlan{T,k,false},x::AbstractVector{T}) where {T,k} =
239-
mul!(similar(x), P, convert(Array,x))
240-
241-
ichebyshevtransform!(x::AbstractVector{T}, kind=Val(1)) where T =
242-
plan_ichebyshevtransform!(x, kind)*x
243-
244-
ichebyshevtransform(x, kind=Val(1)) = ichebyshevtransform!(Array(x), kind)
224+
*(P::IChebyshevTransformPlan{T,kind,K,false},x::AbstractVector{T}) where {T,kind,K} = mul!(similar(x), P, convert(Array,x))
225+
ichebyshevtransform!(x::AbstractArray, dims...; kwds...) = plan_ichebyshevtransform!(x, dims...; kwds...)*x
226+
ichebyshevtransform(x, dims...; kwds...) = plan_ichebyshevtransform(x, dims...; kwds...)*x
245227

246228

247229
## Chebyshev U
@@ -290,7 +272,6 @@ plan_chebyshevutransform(x::AbstractVector) = plan_chebyshevutransform(x, Val(1)
290272

291273
function *(P::ChebyshevUTransformPlan{T,1,true},x::AbstractVector{T}) where T
292274
n = length(x)
293-
assert_applicable(P, x)
294275
n  1 && return x
295276

296277
for k=1:n # sqrt(1-x_j^2) weight
@@ -301,7 +282,6 @@ end
301282

302283
function *(P::ChebyshevUTransformPlan{T,2,true}, x::AbstractVector{T}) where T
303284
n = length(x)
304-
assert_applicable(P, x)
305285
n  1 && return x
306286

307287
c = one(T)/ (n+1)
@@ -370,7 +350,6 @@ plan_ichebyshevutransform(x::AbstractVector) = plan_ichebyshevutransform(x, Val(
370350

371351
function *(P::IChebyshevUTransformPlan{T,1,true}, x::AbstractVector{T}) where T<:fftwNumber
372352
n = length(x)
373-
assert_applicable(P, x)
374353
n  1 && return x
375354

376355
x = P.plan * x
@@ -384,7 +363,6 @@ end
384363

385364
function *(P::IChebyshevUTransformPlan{T,2,true}, x::AbstractVector{T}) where T<:fftwNumber
386365
n = length(x)
387-
assert_applicable(P, x)
388366
n  1 && return x
389367

390368
c = one(T)/ (n+1)

0 commit comments

Comments
 (0)