W'W and W'T mass matrices (#99)

Author: dlfivefifty

src/classical/chebyshev.jl

@@ -162,6 +162,21 @@ end
 
 
 
+@simplify function *(A::QuasiAdjoint{<:Any,<:Weighted{<:Any,<:ChebyshevU}}, B::Weighted{<:Any,<:ChebyshevU})
+    T = promote_type(eltype(A), eltype(B))
+    f = (k,j) ->  isodd(j-k) ? zero(T) : -((2*(one(T) + (-1)^(j + k))*(1 + j)*(1 + k))/((-1 + j - k)*(1 + j - k)*(1 + j + k)*(3 + j + k)))
+    BroadcastMatrix{T}(f, 0:∞, (0:∞)')
+end
+
+
+@simplify function *(A::QuasiAdjoint{<:Any,<:Weighted{<:Any,<:ChebyshevU}}, B::ChebyshevT)
+    T = promote_type(eltype(A), eltype(B))
+    W = parent(A)
+    U = ChebyshevU{T}()
+    (W'U) * (U\B)
+end
+
+
 ##########
 # Derivatives
 ##########

test/test_chebyshev.jl

@@ -460,13 +460,18 @@ import ContinuumArrays: MappedWeightedBasisLayout, Map, WeightedBasisLayout
     @testset "innerproduct" begin
         T = ChebyshevT()
         U = ChebyshevU()
+        W = Weighted(U)
         x = axes(T,1)
 
         @test (T'U)[1:10,1:10] ≈ (U'T)[1:10,1:10]'
         @test (T'T)[1:10,1:10] ≈ ((T'U)*(U\T))[1:10,1:10]
         @test_broken (U'U)[1:10,1:10] ≈ ((U'T)*(T\U))[1:10,1:10]
         f,g = (T/T\exp.(x)),(U/U\exp.(x))
         @test f'g ≈ g'f ≈ f'f ≈ dot(f,g) ≈ dot(f,f) ≈ exp(2)/2 - exp(-2)/2
+
+        h = W * (U\exp.(x))
+        @test h'h ≈ cosh(2) - sinh(2)/2
+        @test f'h ≈ h'f ≈ 2.498566528528904 # 1/2 π BesselI[1, 2]
     end
 end