example with `Chebyshev`

[rveltz]

Hi,

Is it possible to make this example work with ChebyshevT ? Right now, it seems you cannot square differentiation_operator(B)?

Thank you

[daanhb]

Thanks for your interest. I'm afraid the example is not easily modified due to hardcoded limitations. The PDE part of the FrameFun package is not well developed. If you are interested in Chebyshev based methods, perhaps you can consider ApproxFun instead?

I'll elaborate a little. The first line that fails, squaring the differentiation operator, does so because by default the (matrix representation of the) operator has size 40x41. That is because differentiating a polynomial lowers its degree. This can be seen, and fixed, as follows:

julia> B = ChebyshevT(41);

julia> D1 = differentiation_operator(B)
Chebyshev differentiation matrix of order 1 and size 41

julia> size(D1)
(40, 41)

julia> src(D1)
Chebyshev polynomials (first kind)
    ↳ length = 41
    ↳ Float64 -> Float64
    ↳ support = -1.0..1.0 (Chebyshev)

julia> dest(D1)
Chebyshev polynomials (first kind)
    ↳ length = 40
    ↳ Float64 -> Float64
    ↳ support = -1.0..1.0 (Chebyshev)

julia> D2 = differentiation_operator(B, B, 1);

julia> size(D2)
(41, 41)

julia> D2*D2
Operator C * C

C   :   Chebyshev differentiation matrix of order 1 and size 41

The second operator D2 above asks for a differentiation operator that maps an expansion in B to another expansion in B, rather than a smaller version of B. (The third argument 1 is the order which seems to be a necessary argument here. A bit clumsy.)

A bigger problem is that the solver (in the invocation of solve) currently assumes that the discretization of the differentiation operator has a known pseudo inverse. That is the case for Fourier, because the operator for any constant coefficient differential operator is diagonal, but it is not the case for Chebyshev. This assumption is hardcoded, I'm afraid. The code was used to illustrate some parts of the PhD thesis of Roel Matthysen and we haven't further developed it since.

One way around this would be to manually construct the discretization matrix. All the necessary ingredients are there, but it would require a fair amount of coding, and it would not result in an efficient solver. It would also yield a fairly unconventional PDE solver, for which there is no theory yet. We will eventually do this, in the not too distant future even if all goes well, but I'm assuming that that is not what you're after right now.

[rveltz]

Thank you for your detailed answer! I thought it was a nice example to try ChebyshevT.

If you are interested in Chebyshev based methods, perhaps you can consider ApproxFun instead?

I tried but it is too slow. I was hoping for a tone down version of ApproxFun where the basis size is fixed, this would be super useful to discretize PDE. This seems to be the realm of BasisFunctions.jl although not entirely in practice. So I'll wait :D