Discrepancy in number of points and number of values on tensor spaces

[michakraus]

When trying to use ApproxFun to interpolate a function on Chebyshev-points and compute derivatives, I came across the following behaviour:

Define two spaces of Chebyshev and ultraspherical polynomials:

using ApproxFun
C2 = Chebyshev(0..1)^2
U2 = Ultraspherical(1, 0..1)^2

Project a function onto the Chebyshev space:

> fc = Fun((x,y) -> cos(π*x) * cos(π*y), C2)
> (length(points(fc)), length(values(fc)), ncoefficients(fc))
(190, 190, 185)

Transform to ultraspherical space:

> fu = Conversion(C2, U2) * fc
> (length(points(fu)), length(values(fu)), ncoefficients(fu))
(196, 361, 190)

According to the docs, values(f::Fun) returns the f::Fun evaluated on points(f), but that does not seem to be the case as the number of points and values are very different. I am wondering what exactly is happening here.

The result is similar when projecting onto the ultraspherical space right away

> fu2 = Fun((x,y) -> cos(π*x) * cos(π*y), U2)
> (length(points(fu2)), length(values(fu2)), ncoefficients(fu2))
(196, 361, 185)

I am also not quite sure why ncoefficients(f) is smaller than length(points(f)).

[dlfivefifty]

This is caused by inconsistency in the ordering by polynomial degree in the coefficients vs the tensor transforms. It will require some thought on a "smart" way to do this.

[MikaelSlevinsky]

We have to contend with the duality of the array (from tensor product) and vector (from diagonal traversal for degree enumeration) representations.

What if we had a coefficient type that stored both a matrix, a vector, a diagonal traversal type (to indicate the ragged stride configuration like the differences between disk/triangle/sphere OPs), and a hash value (or integer or Boolean) to tell us which backend is active. Then, one can work with the vector type for differential operators, and use the matrix type for computing other things like transforms & evaluation. The "transpose" only takes place if necessary to perform an operation.

[dlfivefifty]

That exists in LazyBandedMatrices as DiagTrav

[MikaelSlevinsky]

Looks like that's only valid for triangles? or does it work for traversals other than anti-diagonal.

[dlfivefifty]

Isn’t it the exact same on the square with Padua points?

Note triangle points double-sample so that’s morally the same as using the tensor grid for the square