ladder operators

Author: dlfivefifty

src/rect.jl

@@ -1,3 +1,13 @@
+"""
+    KronPolynomial(A, B, C...)
+
+represents the kronecker product of `A`, `B`, …. In particular, if `K = KronPolynomial(A,B)` and `U` is an infinite
+matrix of coefficients we have
+``
+K[SVector(x,y),:]'DiagTrav(U) == A[x,:]'U*B[y,:]
+``
+"""
+
 struct KronPolynomial{d, T, PP} <: MultivariateOrthogonalPolynomial{d, T}
     args::PP
 end

test/test_rectdisk.jl

@@ -1,6 +1,6 @@
 using MultivariateOrthogonalPolynomials, StaticArrays, BlockArrays, BlockBandedMatrices, ArrayLayouts, Base64,
         QuasiArrays, Test, ClassicalOrthogonalPolynomials, BandedMatrices, FastTransforms, LinearAlgebra
-import MultivariateOrthogonalPolynomials: dunklxu_raising, dunklxu_lowering, AngularMomentum
+import MultivariateOrthogonalPolynomials: dunklxu_raising, dunklxu_lowering, AngularMomentum, coefficients, basis
 using ForwardDiff
 
 @testset "Dunkl-Xu disk" begin
@@ -114,12 +114,53 @@ using ForwardDiff
         x,y = 𝐱 = SVector(0.1,0.2)
         ρ = sqrt(1-x^2)
         ρ′ = -x/ρ
-        n,k = 3,1
-        P = DunklXuDisk()
+        n,k = 4,2
+        β = 0
+        P = DunklXuDisk(β)
         K = Block(n+1)[k+1]
-        @test Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * diff(Legendre())[y/ρ,k+1] ≈ ρ * diff(P,(0,1))[ 𝐱 ,K]
+        @test P[𝐱,K] ≈ Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Legendre()[y/ρ,k+1]
+        @test Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * diff(Legendre())[y/ρ,k+1] ≈ ρ * diff(P,(0,1))[𝐱 ,K]
         @test diff(Jacobi(k+1/2,k+1/2))[x,n-k+1] * ρ^(k+1) * Legendre()[y/ρ,k+1] ≈ -k*ρ′*P[𝐱,K] + y*ρ′*diff(P,(0,1))[𝐱,K] + ρ * diff(P,(1,0))[𝐱,K]
 
-        
+        dunklxudisk(n, k, β, γ, x, y) = jacobip(n-k, k+β+γ+1/2, k+β+γ+1/2, x)* (1-x^2)^(k/2) * jacobip(k, β, β, y/sqrt(1-x^2))
+
+        A = coefficients(diff(Jacobi(k+1/2,k+1/2)))
+        B = Jacobi(1,1)\diff(Legendre())
+        # M₀₁
+        @test diff(P,(0,1))[𝐱 ,K] ≈ Jacobi(k-1+3/2,k-1+3/2)[x,n-k+1] * ρ^(k-1) * Jacobi(1,1)[y/ρ,k]B[k,k+1] ≈ DunklXuDisk(1)[𝐱,Block(n)[k]]B[k,k+1] ≈ dunklxudisk(n-1,k-1,1,0,x,y)B[k,k+1]
+        # M₀₂
+        @test (k+1)*Legendre()[x,k+1] + (1+x)*diff(Legendre())[x,k+1] ≈ (k+1)*Jacobi(1,0)[x,k+1] #L₄
+        @test (k+1)*Legendre()[y/ρ,k+1] + (1+y/ρ)*diff(Legendre())[y/ρ,k+1] ≈ (k+1)*Jacobi(1,0)[y/ρ,k+1]
+        @test (k+1)*P[𝐱,K] + (1+y/ρ)*ρ *diff(P,(0,1))[𝐱 ,K]  ≈ (k+1)*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Jacobi(1,0)[y/ρ,k+1]
+        # M₀₄
+        @test -(1-x)*(k+1)*Legendre()[x,k+1] - (1-x^2)*diff(Legendre())[x,k+1] ≈ 2*(k+1)*Jacobi(-1,0)[x,k+2] #L₄
+        @test -(1-y/ρ)*(k+1)*P[𝐱,K] - (1-y^2/ρ^2)*ρ *diff(P,(0,1))[𝐱 ,K] ≈ 2*(k+1)*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Jacobi(-1,0)[y/ρ,k+2]
+        # M₁₀
+        @test diff(Jacobi(k+1/2,k+1/2))[x,n-k+1] ≈ (n+k+2)/2 * Jacobi(k+3/2,k+3/2)[x,n-k]
+        @test 1/ρ * (-k*ρ′*P[𝐱,K] + y*ρ′*diff(P,(0,1))[𝐱,K] + ρ * diff(P,(1,0))[𝐱,K]) ≈ (n+k+2)/2 * Jacobi(k+3/2,k+3/2)[x,n-k] * ρ^k * Legendre()[y/ρ,k+1] ≈ (n+k+2)/2 * dunklxudisk(n-1,k,0,1,x,y)
+        # M₆₀
+        @test (k+1/2)*Jacobi(k+1/2,k+1/2)[x,n-k+1] + (1+x)*diff(Jacobi(k+1/2,k+1/2))[x,n-k+1] ≈ (n+1/2)*Jacobi(k+3/2,k-1/2)[x,n-k+1] # L₆
+        @test (k+1/2)*Jacobi(k+1/2,k+1/2)[x,n-k+1]* ρ^k * Legendre()[y/ρ,k+1] + (1+x)*diff(Jacobi(k+1/2,k+1/2))[x,n-k+1]* ρ^k * Legendre()[y/ρ,k+1] ≈ (n+1/2)*Jacobi(k+3/2,k-1/2)[x,n-k+1]* ρ^k * Legendre()[y/ρ,k+1]
+        @test (k+1/2)*P[𝐱,K] + (1+x)/ρ*( -k*ρ′*P[𝐱,K] + y*ρ′*diff(P,(0,1))[𝐱,K] + ρ * diff(P,(1,0))[𝐱,K]) ≈ (n+1/2)*Jacobi(k+3/2,k-1/2)[x,n-k+1] * ρ^k * Legendre()[y/ρ,k+1]
+        # M₀₁'
+        @test -(1-x^2)*diff(Legendre())[x,k+1] ≈ 2*(k+1)*Jacobi(-1+eps(),-1+eps())[x,k+2] #L₁'
+        @test -(1-y^2/ρ^2)*diff(Legendre())[y/ρ,k+1] ≈ 2*(k+1)*Jacobi(-1+eps(),-1+eps())[y/ρ,k+2]
+        @test -(1-y^2/ρ^2)* ρ^2 * diff(P,(0,1))[𝐱 ,K] ≈ 2*(k+1)*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^(k+1) *Jacobi(-1+eps(),-1+eps())[y/ρ,k+2] ≈ 2*(k+1)*dunklxudisk(n+1,k+1,-1+eps(),0,x,y)
+        # M₀₂'
+        @test (1-x)*k*Legendre()[x,k+1] - (1-x^2)*diff(Legendre())[x,k+1] ≈ 2*k*Jacobi(-1,0)[x,k+1] #L₂'
+        @test (1-y/ρ)*k*Legendre()[y/ρ,k+1] - (1-y^2/ρ^2)*diff(Legendre())[y/ρ,k+1] ≈ 2*k*Jacobi(-1,0)[y/ρ,k+1]
+        @test (1-y/ρ)*k*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Legendre()[y/ρ,k+1] - (1-y^2/ρ^2)*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * diff(Legendre())[y/ρ,k+1] ≈ 2*k*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Jacobi(-1,0)[y/ρ,k+1] #L₂'
+        @test (1-y/ρ)*k*P[𝐱,K] - (1-y^2/ρ^2)*ρ *diff(P,(0,1))[𝐱 ,K] ≈ 2*k*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Jacobi(-1,0)[y/ρ,k+1]
+        # M₀₄'
+        @test -k*Legendre()[x,k+1] + (1+x)*diff(Legendre())[x,k+1] ≈ k*Jacobi(1,0)[x,k] #L₄'
+        @test -k*P[𝐱,K] + (1+y/ρ)*ρ *diff(P,(0,1))[𝐱 ,K] ≈ k*Jacobi(k+1/2,k+1/2)[x,n-k+1] * ρ^k * Jacobi(1,0)[y/ρ,k]
+        # M₁₀'
+        @test ((1+x)*(n-k+1/2) - (1-x)*(n-k+1/2))*Jacobi(k+1/2,k+1/2)[x,n-k+1] - (1-x^2)*diff(Jacobi(k+1/2,k+1/2))[x,n-k+1] ≈ 2*(n-k+1)*Jacobi(k-1/2,k-1/2)[x,n-k+2] # L₁'
+        @test ((1+x)*(n-k+1/2) - (1-x)*(n-k+1/2))*P[𝐱,K] - (1-x^2)/ρ * (-k*ρ′*P[𝐱,K] + y*ρ′*diff(P,(0,1))[𝐱,K] + ρ * diff(P,(1,0))[𝐱,K]) ≈ 2*(n-k+1)*Jacobi(k-1/2,k-1/2)[x,n-k+2]* ρ^k * Legendre()[y/ρ,k+1]
+        # M₆₀'
+        @test (k+1/2)*Jacobi(k+1/2,k+1/2)[x,n-k+1]-(1-x)*diff(Jacobi(k+1/2,k+1/2))[x,n-k+1] ≈ (n+1/2)*Jacobi(k-1/2,k+3/2)[x,n-k+1] # L₆'
+        @test (k+1/2)*P[𝐱,K]-(1-x)/ρ * (-k*ρ′*P[𝐱,K] + y*ρ′*diff(P,(0,1))[𝐱,K] + ρ * diff(P,(1,0))[𝐱,K]) ≈ (n+1/2)*Jacobi(k-1/2,k+3/2)[x,n-k+1]* ρ^k * Legendre()[y/ρ,k+1]
+
+        # d/dx = M₁₀M₀₂ + M₀₄'M₆₀'
     end
 end