Skip to content

Commit 58b847a

Browse files
dlfivefiftydlfivefifty
committed
made it to synthesis test
1 parent 47162aa commit 58b847a

File tree

12 files changed

+58
-34
lines changed

12 files changed

+58
-34
lines changed

‎src/FastTransforms.jl

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -13,9 +13,11 @@ if VERSION < v"0.7-"
1313
else
1414
using FFTW, LinearAlgebra, DSP
1515
import FFTW: r2rFFTWPlan, unsafe_execute!, fftwSingle, fftwDouble, fftwNumber
16-
import FFTW: libfftw, libfftwf, PlanPtr, r2rFFTWPlan
16+
import FFTW: libfftw3, libfftw3f, PlanPtr, r2rFFTWPlan
1717
const LAmul! = LinearAlgebra.mul!
18-
import LinearAlgebra: Factorization
18+
const libfftw = libfftw3
19+
const libfftwf = libfftw3f
20+
import LinearAlgebra: Factorization, transpose, adjoint
1921
flipdim(A,d) = reverse(A; dims=d)
2022
end
2123

‎src/SphericalHarmonics/Butterfly.jl

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -55,7 +55,7 @@ function Butterfly(A::AbstractMatrix{T}, L::Int; isorthogonal::Bool = false, opt
5555
nu = round(Int, ninds[j+1]) - 1
5656
nd = nu-nl+1
5757
if isorthogonal
58-
factors[1][j] = IDPackedV{T}(collect(1:nd),Int[],Array{T}(nd,0))
58+
factors[1][j] = IDPackedV{T}(collect(1:nd),Int[],Array{T}(undef,nd,0))
5959
else
6060
factors[1][j] = idfact!(A[:,nl:nu], LRAOpts)
6161
end
@@ -83,7 +83,7 @@ function Butterfly(A::AbstractMatrix{T}, L::Int; isorthogonal::Bool = false, opt
8383
lc = length(cols)
8484
Av = A[ml:mu,cols]
8585
if maximum(abs, Av) < realmin(real(T))/eps(real(T))
86-
factors[l][j+ctr] = IDPackedV{T}(Int[], collect(1:lc), Array{T}(0,lc))
86+
factors[l][j+ctr] = IDPackedV{T}(Int[], collect(1:lc), Array{T}(undef,0,lc))
8787
else
8888
LRAOpts.rtol = eps(real(T))*max(mu-ml+1, lc)
8989
factors[l][j+ctr] = idfact!(Av, LRAOpts)

‎src/SphericalHarmonics/SphericalHarmonics.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -5,7 +5,7 @@ function *(P::SphericalHarmonicPlan, X::AbstractMatrix)
55
end
66

77
function \(P::SphericalHarmonicPlan, X::AbstractMatrix)
8-
At_mul_B!(zero(X), P, X)
8+
mul!(zero(X), transpose(P), X)
99
end
1010

1111
include("sphfunctions.jl")

‎src/SphericalHarmonics/fastplan.jl

Lines changed: 10 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -16,10 +16,10 @@ function FastSphericalHarmonicPlan(A::Matrix{T}, L::Int; opts...) where T
1616
p2 = plan_normleg12cheb2(A)
1717
p1inv = plan_cheb2normleg(A)
1818
p2inv = plan_cheb22normleg1(A)
19-
B = zeros(A)
20-
Ce = eye(T, M)
21-
Co = eye(T, M)
22-
BF = Vector{Butterfly{T}}(n-2)
19+
B = zero(A)
20+
Ce = Matrix{T}(I, M, M)
21+
Co = Matrix{T}(I, M, M)
22+
BF = Vector{Butterfly{T}}(undef, n-2)
2323
P = Progress(n-2, 0.1, "Pre-computing...", 43)
2424
for j = 1:2:n-2
2525
mul!(Ce, RP.layers[j])
@@ -36,6 +36,12 @@ end
3636

3737
FastSphericalHarmonicPlan(A::Matrix; opts...) = FastSphericalHarmonicPlan(A, floor(Int, log2(size(A, 1)+1)-6); opts...)
3838

39+
if VERSION ≥ v"0.7-"
40+
adjoint(P::FastSphericalHarmonicPlan) = Adjoint(P)
41+
transpose(P::FastSphericalHarmonicPlan) = Transpose(P)
42+
end
43+
44+
3945
function LAmul!(Y::Matrix, FP::FastSphericalHarmonicPlan, X::Matrix)
4046
RP, BF, p1, p2, B = FP.RP, FP.BF, FP.p1, FP.p2, FP.B
4147
fill!(B, zero(eltype(B)))

‎src/SphericalHarmonics/slowplan.jl

Lines changed: 11 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -46,7 +46,7 @@ struct Pnmp2toPlm{T} <: AbstractRotation{T}
4646
end
4747

4848
function Pnmp2toPlm(::Type{T}, n::Int, m::Int) where T
49-
G = Vector{Givens{T}}(n)
49+
G = Vector{Givens{T}}(undef, n)
5050
@inbounds for ℓ = 1:n
5151
c = sqrt(T((2m+2)*(2ℓ+2m+3))/T((ℓ+2m+2)*(ℓ+2m+3)))
5252
s = sqrt(T(ℓ*(ℓ+1))/T((ℓ+2m+2)*(ℓ+2m+3)))
@@ -80,7 +80,7 @@ struct RotationPlan{T} <: AbstractRotation{T}
8080
end
8181

8282
function RotationPlan(::Type{T}, n::Int) where T
83-
layers = Vector{Pnmp2toPlm{T}}(n-1)
83+
layers = Vector{Pnmp2toPlm{T}}(undef, n-1)
8484
@inbounds for m = 0:n-2
8585
layers[m+1] = Pnmp2toPlm(T, n-1-m, m)
8686
end
@@ -99,6 +99,11 @@ function RotationPlan(::Type{T}, n::Int) where T
9999
RotationPlan(layers, snm, cnm)
100100
end
101101

102+
if VERSION ≥ v"0.7-"
103+
adjoint(P::RotationPlan) = Adjoint(P)
104+
transpose(P::RotationPlan) = Transpose(P)
105+
end
106+
102107
function LAmul!(P::RotationPlan, A::AbstractMatrix)
103108
N, M = size(A)
104109
snm = P.snm
@@ -319,7 +324,7 @@ function SlowSphericalHarmonicPlan(A::Matrix{T}) where T
319324
p2 = plan_normleg12cheb2(A)
320325
p1inv = plan_cheb2normleg(A)
321326
p2inv = plan_cheb22normleg1(A)
322-
B = zeros(A)
327+
B = zero(A)
323328
SlowSphericalHarmonicPlan(RP, p1, p2, p1inv, p2inv, B)
324329
end
325330

@@ -359,6 +364,9 @@ if VERSION < v"0.7-"
359364

360365
Base.Ac_mul_B!(Y::Matrix, SP::SlowSphericalHarmonicPlan, X::Matrix) = At_mul_B!(Y, SP, X)
361366
else
367+
transpose(P::SlowSphericalHarmonicPlan) = Transpose(P)
368+
adjoint(P::SlowSphericalHarmonicPlan) = Adjoint(P)
369+
362370
function LinearAlgebra.mul!(Y::Matrix, SPt::Transpose{T,<:SlowSphericalHarmonicPlan}, X::Matrix) where T
363371
SP = parent(SPt)
364372
RP, p1inv, p2inv, B = SP.RP, SP.p1inv, SP.p2inv, SP.B

‎src/SphericalHarmonics/synthesisanalysis.jl

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -224,12 +224,13 @@ function row_analysis!(P, C, vals::Vector{T}) where T
224224
cfs[n÷2+1] *= half(T)
225225
end
226226

227-
negateeven!(reverseeven!(mul!(C, cfs)))
227+
negateeven!(reverseeven!(lmul!(C, cfs)))
228228
end
229229

230+
230231
function row_synthesis!(P, C, cfs::Vector{T}) where T
231232
n = length(cfs)
232-
Ac_mul_B!(C, reverseeven!(negateeven!(cfs)))
233+
lmul!(C, reverseeven!(negateeven!(cfs)))
233234
if iseven(n)
234235
cfs[n÷2+1] *= two(T)
235236
end

‎src/SphericalHarmonics/thinplan.jl

Lines changed: 9 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -20,10 +20,10 @@ function ThinSphericalHarmonicPlan(A::Matrix{T}, L::Int; opts...) where T
2020
p2 = plan_normleg12cheb2(A)
2121
p1inv = plan_cheb2normleg(A)
2222
p2inv = plan_cheb22normleg1(A)
23-
B = zeros(A)
24-
Ce = eye(T, M)
25-
Co = eye(T, M)
26-
BF = Vector{Butterfly{T}}(n-2)
23+
B = zero(A)
24+
Ce = Matrix{T}(I, M, M)
25+
Co = Matrix{T}(I, M, M)
26+
BF = Vector{Butterfly{T}}(undef, n-2)
2727
P = Progress(n-2, 0.1, "Pre-computing...", 43)
2828
for J = 1:2:n-2
2929
mul!(Ce, RP.layers[J])
@@ -40,6 +40,11 @@ end
4040

4141
ThinSphericalHarmonicPlan(A::Matrix; opts...) = ThinSphericalHarmonicPlan(A, floor(Int, log2(size(A, 1)+1)-6); opts...)
4242

43+
if VERSION ≥ v"0.7-"
44+
adjoint(P::ThinSphericalHarmonicPlan) = Adjoint(P)
45+
transpose(P::ThinSphericalHarmonicPlan) = Transpose(P)
46+
end
47+
4348
function LAmul!(Y::Matrix, TP::ThinSphericalHarmonicPlan, X::Matrix)
4449
RP, BF, p1, p2, B = TP.RP, TP.BF, TP.p1, TP.p2, TP.B
4550
copyto!(B, X)

‎test/butterflytests.jl

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -5,6 +5,7 @@ if VERSION < v"0.7-"
55
const mul! = Base.A_mul_B!
66
else
77
Ac_mul_B!(C, A, B) = mul!(C, A', B)
8+
At_mul_B!(C, A, B) = mul!(C, transpose(A), B)
89
end
910

1011
import FastTransforms: Butterfly

‎test/sphericalharmonics/fastplantests.jl

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -16,8 +16,8 @@ import FastTransforms: normalizecolumns!, maxcolnorm
1616
@time FP = FastSphericalHarmonicPlan(A; sketch = :none);
1717
@time SP = SlowSphericalHarmonicPlan(A);
1818

19-
@time A_mul_B!(B, SP, A);
20-
@time A_mul_B!(C, FP, A);
19+
@time mul!(B, SP, A);
20+
@time mul!(C, FP, A);
2121
@time At_mul_B!(D, SP, B);
2222
@time At_mul_B!(E, FP, C);
2323

‎test/sphericalharmonics/slowplantests.jl

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -3,8 +3,9 @@ using Compat.Test
33

44
import FastTransforms: normalizecolumns!, maxcolnorm
55

6+
67
@testset "Slow plan" begin
7-
N = round.([Int],logspace(1,2.5,10))
8+
N = round.([Int],10 .^ range(1,stop=2.5,length=10))
89

910
t = zeros(length(N))
1011
err = zeros(length(N))
@@ -20,7 +21,7 @@ import FastTransforms: normalizecolumns!, maxcolnorm
2021
Ac = copy(A)
2122
B = zero(A)
2223
SP = SlowSphericalHarmonicPlan(A)
23-
A_mul_B!(B, SP, A)
24+
mul!(B, SP, A)
2425
fill!(A, 0.0)
2526
t[j] += @elapsed At_mul_B!(A, SP, B)
2627
nrms[kk] = maxcolnorm(A - Ac)

0 commit comments

Comments
 (0)