Wrong tensor product between Fourier and Chebyshev

[mzaffalon]

using ApproxFun
f₁ = Fun(Fourier(), [zeros(2); 1]);
f₂ = Fun(Chebyshev(), [zeros(3); 1]);
f = f₁ ⊗ f₂;
coefficients(f) # all coefficients are zero

with ApproxFun v0.13.25 and Julia v1.10.0. The problem seems to be associated with Chebyshev: increasing its order past n=2 makes all coefficients zero.

The evaluation of the individual terms work:

ϕs = range(0, 2π, 51);
zs = range(-1, 1, 31);
[f₁(ϕ) * f₂(z) for ϕ in ϕs, z in zs]

[mzaffalon]

The problem is with the tensor product with two Fourier spaces:

using ApproxFun
f = Fun(Fourier(), [0; 0; 1]);
kron(f,f).coefficients

This if-statement breaks too early. According to the documentation,, the ordering in the function is assumed to be (1,1), (1,2), (2,1), (1,3), (2,2)... but for Fourier spaces it appears that the ordering is instead (1,1), (1,2), (2,1), (2,2), (1,3).

What does TODO: generalize refer too?