Merge pull request #9 from JuliaApproximation/jp/cache

Jp/cache

.github/workflows/ci.yml

@@ -12,7 +12,7 @@ jobs:
           - '1'
         os:
           - ubuntu-latest
-          - macOS-latest
+          # - macOS-latest
         arch:
           - x64
     steps:

Project.toml

@@ -1,7 +1,7 @@
 name = "AnnuliOrthogonalPolynomials"
 uuid = "de1797fd-24c3-4035-91a2-b52aecdcfb01"
 authors = ["john.papad "]
-version = "0.0.4"
+version = "0.0.5"
 
 [deps]
 BlockArrays = "8e7c35d0-a365-5155-bbbb-fb81a777f24e"

src/annulus.jl

@@ -131,8 +131,14 @@ function \(A::ZernikeAnnulus{T}, B::Weighted{V,ZernikeAnnulus{V}}) where {T,V}
 
     ρ = convert(TV, A.ρ); t=inv(one(TV)-ρ^2)
 
-    L₁ = Weighted.(SemiclassicalJacobi{real(TV)}.(t,zero(TV),zero(TV),zero(TV):∞)) .\ Weighted.(SemiclassicalJacobi{real(TV)}.(t,one(TV),one(TV),zero(TV):∞))
-    L₂ = SemiclassicalJacobi{real(TV)}.(t,one(TV),one(TV),zero(TV):∞) .\ SemiclassicalJacobi{real(TV)}.(t,zero(TV),zero(TV),zero(TV):∞)
+    # L₁ = Weighted.(SemiclassicalJacobi{real(TV)}.(t,zero(TV),zero(TV),zero(TV):∞)) .\ Weighted.(SemiclassicalJacobi{real(TV)}.(t,one(TV),one(TV),zero(TV):∞))
+    # L₂ = SemiclassicalJacobi{real(TV)}.(t,one(TV),one(TV),zero(TV):∞) .\ SemiclassicalJacobi{real(TV)}.(t,zero(TV),zero(TV),zero(TV):∞)
+
+    Q₀₀ = SemiclassicalJacobi{real(TV)}.(t,0,0,0:∞)
+    Q₁₁ = SemiclassicalJacobi{real(TV)}.(t,1,1,0:∞)
+
+    L₁ = Weighted.(Q₀₀) .\ Weighted.(Q₁₁)
+    L₂ = Q₁₁ .\ Q₀₀
 
     # L = (one(TV)-ρ^2)^2 .* (L₂ .* L₁)
     # Workaround for broken lazy multiplication
@@ -210,6 +216,7 @@ function plotvalues(u::ApplyQuasiVector{T,typeof(*),<:Tuple{Weighted{<:Any,<:Zer
     w[x] .* U
 end
 
+
 # FastTransforms uses orthonormalized annulus OPs so we need to correct the normalization
 # Gives a vector with the weights of each m-mode.
 function normalize_mmodes(w::AnnulusWeight{T}) where T
@@ -238,6 +245,35 @@ function denormalize_annulus(A::AbstractVector, a, b, c, ρ, analysis=true)
     A ./ d # divide vector by denormalization if synthesis
 end
 
+# # FastTransforms uses orthonormalized annulus OPs so we need to correct the normalization
+# # Gives a vector with the weights of each m-mode.
+# function normalize_mmodes(w::AnnulusWeight{T}, l::Int) where T
+#     (ρ,a,b) = convert(T,w.ρ), convert(T,w.a), convert(T,w.b)
+#     t = inv(one(T)-ρ^2)
+#     # contribution from non-normalized semiclassical Jacobi
+#     jw = sum.(SemiclassicalJacobiWeight.(t,b,a,0:(l-1)))
+#     # π comes from contribution of Harmonic polynomials and t^(a+b+one(T)) from the change of variables
+#     m₀ =  sqrt( convert(T,π) / t^(a+b+one(T)) ) * sqrt(jw[1])
+#     vcat([m₀], sqrt.( convert(T,π) ./ (2 * t.^( (a+b+2*one(T)):(convert(T,l)+a+b)))) .* sqrt.(jw[2:end]))
+# end
+
+# # Creates vector correctly interlacing the denormalization
+# # constants for each m-mode.
+# function denormalize_annulus(A::AbstractVector, a, b, c, ρ, analysis=true)
+#     supp = Int(last(blockcolsupport(A)))
+#     l = length(A[Block.(1:supp)])
+#     bl = [findblockindex(blockedrange(oneto(∞)), j) for j in 1:l]
+#     ℓ = [bl[j].I[1]-1 for j in 1:l] # degree
+#     k = [bl[j].α[1] for j in 1:l]   # index of degree
+#     m = [iseven(ℓ[j]) ? k[j]-isodd(k[j]) : k[j]-iseven(k[j]) for j in 1:l] # m-mode
+#     s = (-1).^Int.((ℓ .- m) ./ 2) # denormalization includes negatives (-1)^(degree - m)/2
+#     w = AnnulusWeight(ρ, a, b)
+#     constants = normalize_mmodes(w, l) # m-mode constants
+#     d = [inv(constants[mm+1]*ss) for (mm, ss) in zip(m, s)] # multiply by relevant (-1)
+#     analysis && return [d.*A[1:l];A[l+1:end]] # multiply vector by denormalization if analysis
+#     [A[1:l]./d;A[l+1:end]]  # divide vector by denormalization if synthesis
+# end
+
 struct ZernikeAnnulusTransform{T} <: Plan{T}
     N::Int
     ann2cxf::FastTransforms.FTPlan{T,2,FastTransforms.ANNULUS}