Conversion between Fourier spaces

[jishnub]

Fixes JuliaApproximation/ApproxFun.jl#768. After this, the following works:

julia> coefficients([1,2,3], Fourier(0..2pi), Fourier(0..4pi))
5-element Vector{Int64}:
 1
 0
 0
 2
 3

julia> union(Fourier(0..2pi), Fourier(0..4pi))
Fourier(【0.0,12.566370614359172❫)

julia> union(Fourier(0..2pi), Fourier(0..3pi))
Fourier(【0.0,6.283185307179586❫)⨄Fourier(【0.0,9.42477796076938❫)

WIP to allow domains that aren't multiples of each other.

[codecov]

Codecov Report

Patch coverage: 98.33% and project coverage change: +2.67 🎉

Comparison is base (11df9e6) 61.80% compared to head (dce96e7) 64.48%.

Additional details and impacted files

@@            Coverage Diff             @@
##           master      #41      +/-   ##
==========================================
+ Coverage   61.80%   64.48%   +2.67%     
==========================================
  Files          14       14              
  Lines        1063     1115      +52     
==========================================
+ Hits          657      719      +62     
+ Misses        406      396      -10     

Impacted Files	Coverage Δ
src/ApproxFunFourier.jl	71.81% <98.14%> (+5.82%)	⬆️
src/calculus.jl	64.19% <100.00%> (+1.37%)	⬆️

... and 3 files with indirect coverage changes

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report in Codecov by Sentry.
📢 Do you have feedback about the report comment? Let us know in this issue.

[andreasvarga]

My understanding is that the result of

union(Fourier(0..2pi), Fourier(0..4pi))

should be

Fourier(【0.0,12.566370614359172❫)

(i.e.,Fourier(0..4pi)) and not

Fourier(【0.0,6.283185307179586❫)

The later would be the result of the intersection.

Also

union(Fourier(0..2pi), Fourier(0..3pi))

could be

Fourier(【0.0,18.84955592153876❫)

(i.e., (i.e.,Fourier(0..6pi)), which is the smallest common period of both domains.

[jishnub]

My understanding is that the result of

union(Fourier(0..2pi), Fourier(0..4pi))

should be

Fourier(【0.0,12.566370614359172❫)

(i.e.,Fourier(0..4pi)) and not

Fourier(【0.0,6.283185307179586❫)

Apologies, that was a copying error on my part. The result is indeed 0..4pi. I've fixed it now.

Regarding union(Fourier(0..2pi), Fourier(0..3pi)), I see why the LCM makes sense, but am unsure if we should change the domain altogether to one that doesn't correspond to either of the two Funs. The domain may have physical meaning, e.g. it may correspond to the dimensions of a box.

[andreasvarga]

Regarding union(Fourier(0..2pi), Fourier(0..3pi)), I see why the LCM makes sense, but am unsure if we should change the domain altogether to one that doesn't correspond to either of the two Funs. The domain may have physical meaning, e.g. it may correspond to the dimensions of a box.

This was the reason why I wrote in my comment could be.

However, assuming
F1 = Fun(x-> sin(x),Fourier(0..2pi))
F2 = Fun(x-> sin(2*x/3),Fourier(0..3pi))
then to represent F1+F2 as a periodic function, I think the only possible periods are of the form 6kπ, for k = 1, 2, ... . Or it is possible to have a periodic function F1+F2 with the domain
union(Fourier(0..2pi), Fourier(0..3pi)) ?

In any case, even if the domain updating operation cannot be supported or is not possible, I would propose to issue an error, to prevent crashes of stack overflow errors.

[jishnub]

Addition of Fun should work in either case after this PR:

julia> f1 = Fun(sin, Fourier(0..2pi));

julia> f2 = Fun(x->sin((2/3)x), Fourier(0..3pi));

julia> f1 + f2
Fun(Fourier(【0.0,6.283185307179586❫)⨄Fourier(【0.0,9.42477796076938❫), [-2.754298361189885e-18, 1.0, -9.693192265097112e-18, 1.0])

julia> f2 = Fun(sin, Fourier(0..4pi));

julia> f1 + f2
Fun(Fourier(【0.0,12.566370614359172❫), [-4.295736460310579e-17, -2.6593986539776536e-17, -1.0473104888558003e-16, 2.0])

The latter of these will be periodic for a period of 6pi as well, although there might be issues that need to be resolved in this case, as I've not tested this yet.

What doesn't work as yet is the conversion of coefficients between Fourier spaces for domains that aren't multiples of each other, but I'm hoping to fix that as well.

[andreasvarga]

Could you check that for f1+f2 in the first example you get the following values:

julia> F = Fun(x-> sin(x)+sin(2*x/3),Fourier(0..6pi))
Fun(Fourier(【0.0,18.84955592153876❫),[-4.711296368749278e-16, 2.636779683484747e-16, -4.364038231364736e-16, 1.0000000000000002, -2.1109009991810788e-16, 1.0000000000000004])

julia> F.([0,pi,2pi,3pi,4pi,5pi,6pi])
7-element Vector{Float64}:
 -1.1186235599295091e-15
  0.8660254037844379
 -0.8660254037844384
 -1.2335123374182663e-16
  0.8660254037844376
 -0.8660254037844365
 -2.343270359076863e-15

[jishnub]

I think f1 + f2 is broken currently for Fourier sum spaces. I'll have to look into it.

[andreasvarga]

I observed some problems when constructing arrays involving different domains. Here are some examples:

julia> f1 = Fun(sin, Fourier(0..2pi));

julia> f2 = Fun(x->sin((2/3)x), Fourier(0..3pi));

julia> f = f1+f2   # OK
Fun(Fourier(【0.0,18.84955592153876❫), [-1.24475e-17, 0.0, 0.0, 1.0, 0.0, 1.0])

julia> domain(f) # OK
【0.0,18.84955592153876❫

julia> domain([f;]) # OK
【0.0,18.84955592153876❫

julia> domain([f;;]) # OK
1×1 Matrix{PeriodicSegment{Float64}}:
 【0.0,18.84955592153876❫

julia> domain([f]) # error
ERROR: MethodError: no method matching domain(::Vector{Fun{ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}, Float64, Vector{Float64}}})

Closest candidates are:
  domain(::SequenceSpace)
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\4bcPp\src\Space.jl:690
  domain(::Fun)
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\4bcPp\src\Fun.jl:284
  domain(::NormalizedPolynomialSpace)
   @ ApproxFunOrthogonalPolynomials C:\Users\Andreas\.julia\packages\ApproxFunOrthogonalPolynomials\PisRd\src\Spaces\PolynomialSpace.jl:438
  ...

Stacktrace:
 [1] top-level scope
   @ REPL[63]:1

The same if f = f1*f2.

julia> f = [f1;f2;;] # OK
2×1 Matrix{Fun{ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}, Float64, Vector{Float64}}}:
 Fun(Fourier(【0.0,6.283185307179586❫), [-2.7543e-18, 1.0])
 Fun(Fourier(【0.0,9.42477796076938❫), [-9.69319e-18, 1.0])

julia> domain(f)  # OK
2×1 Matrix{PeriodicSegment{Float64}}:
 【0.0,6.283185307179586❫
 【0.0,9.42477796076938❫

julia> f = [f1;f2] # ??
Fun(ApproxFunBase.ArraySpace(ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,6.283185307179586❫), Fourier(【0.0,9.42477796076938❫)]), [-2.7543e-18, 1.0, -9.69319e-18, 1.0])

julia> domain(f) # incomplete
【0.0,6.283185307179586❫

julia> f = [f1 f2] # ??
Fun(ApproxFunBase.ArraySpace(ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,6.283185307179586❫) Fourier(【0.0,9.42477796076938❫)]), [-2.7543e-18, 1.0, -9.69319e-18, 1.0])

julia> domain(f) # incomplete
【0.0,6.283185307179586❫

In my opinion f should be a 1 x 2 matrix, as in the simple case [1 2]. The corresponding domain should be also adjusted.