ConcreteConversion in Chebyshev-Jacobi (#89)

* ConcreteConversion in Chebyshev-Jacobi

* add tests for deriv

* tests for Chebyshev->Jacobi conversions

* tests for different domains

Author: jishnub

src/Spaces/Jacobi/JacobiOperators.jl

@@ -300,14 +300,14 @@ function Conversion(A::Jacobi,B::Chebyshev)
     elseif A.a == A.b == 0
         ConversionWrapper(
             SpaceOperator(
-                Conversion(Ultraspherical(1//2),B),
+                ConcreteConversion(Ultraspherical(1//2),B),
                 A,B))
     elseif A.a == A.b
         US = Ultraspherical(A)
-        ConversionWrapper(Conversion(US,B)*Conversion(A,US))
+        ConversionWrapper(Conversion(US,B)*ConcreteConversion(A,US))
     else
         J = Jacobi(B)
-        Conversion(J,B)*Conversion(A,J)
+        ConcreteConversion(J,B)*Conversion(A,J)
     end
 end
 
@@ -317,21 +317,22 @@ function Conversion(A::Chebyshev,B::Jacobi)
     elseif B.a == B.b == 0
         ConversionWrapper(
             SpaceOperator(
-                Conversion(A,Ultraspherical(1//2,domain(B))),
+                ConcreteConversion(A,Ultraspherical(1//2,domain(B))),
                 A,B))
     elseif B.a == B.b
         US = Ultraspherical(B)
-        ConversionWrapper(Conversion(US,B)*Conversion(A,US))
+        ConcreteConversion(US,B) * Conversion(A,US)
     else
         J = Jacobi(A)
-        Conversion(J,B)*Conversion(A,J)
+        Conversion(J,B)*ConcreteConversion(A,J)
     end
 end
 
 
 function Conversion(A::Jacobi,B::Ultraspherical)
     if A.a == A.b == -0.5
-        ConversionWrapper(Conversion(Chebyshev(domain(A)),B)*Conversion(A,Chebyshev(domain(A))))
+        ConversionWrapper(Conversion(Chebyshev(domain(A)),B)*
+            ConcreteConversion(A,Chebyshev(domain(A))))
     elseif A.a == A.b == order(B)-0.5
         ConcreteConversion(A,B)
     elseif A.a == A.b == 0
@@ -341,16 +342,17 @@ function Conversion(A::Jacobi,B::Ultraspherical)
                 A,B))
     elseif A.a == A.b
         US = Ultraspherical(A)
-        ConversionWrapper(Conversion(US,B)*Conversion(A,US))
+        ConversionWrapper(Conversion(US,B)*ConcreteConversion(A,US))
     else
         J = Jacobi(B)
-        Conversion(J,B)*Conversion(A,J)
+        ConcreteConversion(J,B)*Conversion(A,J)
     end
 end
 
 function Conversion(A::Ultraspherical,B::Jacobi)
     if B.a == B.b == -0.5
-        ConversionWrapper(Conversion(Chebyshev(domain(A)),B)*Conversion(A,Chebyshev(domain(A))))
+        ConversionWrapper(ConcreteConversion(Chebyshev(domain(A)),B)*
+            Conversion(A,Chebyshev(domain(A))))
     elseif B.a == B.b == order(A)-0.5
         ConcreteConversion(A,B)
     elseif B.a == B.b == 0
@@ -360,10 +362,10 @@ function Conversion(A::Ultraspherical,B::Jacobi)
                 A,B))
     elseif B.a == B.b
         US = Ultraspherical(B)
-        ConversionWrapper(Conversion(US,B)*Conversion(A,US))
+        ConversionWrapper(ConcreteConversion(US,B)*Conversion(A,US))
     else
         J = Jacobi(A)
-        Conversion(J,B)*Conversion(A,J)
+        Conversion(J,B)*ConcreteConversion(A,J)
     end
 end
 

src/Spaces/Ultraspherical/UltrasphericalOperators.jl

@@ -123,8 +123,8 @@ function Conversion(A::Chebyshev,B::Ultraspherical)
     else
         d=domain(A)
         US=Ultraspherical(order(B)-1,d)
-        ConversionWrapper(TimesOperator(ConcreteConversion(US,B),
-                                        Conversion(Chebyshev(d),US)))
+        ConversionWrapper(TimesOperator(
+            ConcreteConversion(US,B), Conversion(A,US)))
     end
 end
 
@@ -148,8 +148,9 @@ function Conversion(A::Ultraspherical,B::Ultraspherical)
         ConcreteConversion(A,B)
     elseif b-a > 1
         d=domain(A)
-        ConversionWrapper(TimesOperator(ConcreteConversion(Ultraspherical(b-1,d),B),
-                                        Conversion(A,Ultraspherical(b-1,d))))
+        US=Ultraspherical(b-1,d)
+        ConversionWrapper(TimesOperator(
+            ConcreteConversion(US,B), Conversion(A,US)))
     else
         throw(ArgumentError("Cannot convert from $A to $B"))
     end

test/JacobiTest.jl

@@ -55,12 +55,40 @@ import ApproxFunOrthogonalPolynomials: jacobip
                 end
             end
         end
+
+        @testset for d in [-1..1, 0..1]
+            f = Fun(x->x^2, Chebyshev(d))
+            C = space(f)
+            for J1 = Any[Jacobi(-0.5, -0.5, d), Legendre(d),
+                            Jacobi(0.5, 0.5, d), Jacobi(2.5, 1.5, d)]
+                for J in [J1, NormalizedPolynomialSpace(J1)]
+                    g = Fun(f, J)
+                    if !any(isnan, coefficients(g))
+                        @test Conversion(C, J) * f ≈ g
+                    end
+                end
+            end
+        end
     end
 
     @testset "Derivative" begin
         D=Derivative(Jacobi(0.,1.,Segment(1.,0.)))
         @time testbandedoperator(D)
 
+
+        @testset for d in [-1..1, 0..1]
+            f = Fun(x->x^2, Chebyshev(d))
+            C = space(f)
+            for J = Any[Jacobi(-0.5, -0.5, d), Legendre(d)]
+                g = (Derivative(J) * Conversion(C, J)) * f
+                h = Derivative(C) * f
+                @test g ≈ h
+
+                g = (Derivative(C) * Conversion(J, C)) * f
+                h = Derivative(J) * f
+                @test g ≈ h
+            end
+        end
         @testset for S1 in Any[Jacobi(0,0),
             Jacobi(0,0,1..2), Jacobi(2,2,1..2), Jacobi(0.5,2.5,1..2)],
                 S in Any[S1, NormalizedPolynomialSpace(S1)]