Changing domain for Fourier

andreasvarga · August 5, 2022, 10:11am

I am using ApproxFun.jl as basis to define and manipulate periodic matrices within the PeridicSystems.jl, a collection of tools dedicated for handling periodic systems. I consider ApproxFun.jl a valuable tool for my purposes, being able to address some of computational problems, which I was not able to solve within the Matlab version of this package. So I would like very much to continue using this package further for my developments and I would be very pleased to have some new functionality which would ease using this package for my purposes.

I encountered the following basic issue in handling periodic matrices with different, but still commensurate periods.

If I have two periodic function matrices F1 of period 8π and F2 of period 2π and I would like to compute their sum, this operation is mathematically possible, the sum F1+F2 is simply a periodic function matrix of period 8π. However, this operation is not supported within ApproxFun.jl, leading to either crashing Julia or stack overflow error:

julia> F1
Fun(2×2 ArraySpace:
ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,25.132741228718345❫) Fourier(【0.0,25.132741228718345❫); Fourier(【0.0,25.132741228718345❫) Fourier(【0.0,25.132741228718345❫)],[-5.000560225745643e-18, 2.5363265666181656e-17, 3.497202527569243e-15, -2.720857218197605e-16, 3.497202527569243e-15, -2.720857218197605e-16, 3.552713678800501e-15, 2.3209220587652208e-15, -6.12323399573677e-17, 6.123233995736767e-17  …  -4.73232564246473e-15, -2.3209220587652208e-15, -3.552713678800501e-14, 5.942090786205257e-16, -1.4782794558091705e-16, 1.790234627208065e-15, -1.519318030172918e-16, 1.790234627208065e-15, -1.519318030172918e-16, 1.0000000000000142])

julia> domain(F1)
【0.0,25.132741228718345❫

julia> F2
Fun(2×2 ArraySpace:
ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,6.283185307179586❫) Fourier(【0.0,6.283185307179586❫); Fourier(【0.0,6.283185307179586❫) Fourier(【0.0,6.283185307179586❫)],[0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0289288258043393e-16, 0.0, -0.0, 0.0, -0.0, 0.0, 1.0])

julia> domain(F2)
【0.0,6.283185307179586❫

julia> F1+F2
ERROR: StackOverflowError:
Stacktrace:
 [1] map(f::ApproxFunBase.var"#288#289"{PeriodicSegment{Float64}}, v::Core.SimpleVector)
   @ Base .\essentials.jl:622
 [2] ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}(dom::PeriodicSegment{Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Spaces\SumSpace.jl:45
 [3] SumSpace
   @ C:\Users\Andreas\.julia\packages\ApproxFunFourier\MPNwF\src\ApproxFunFourier.jl:405 [inlined]
 [4] canonicalspace
   @ C:\Users\Andreas\.julia\packages\ApproxFunFourier\MPNwF\src\ApproxFunFourier.jl:97 [inlined]
 [5] map
   @ .\tuple.jl:222 [inlined]
 [6] union_rule(S1sp::ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}, S2sp::ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Spaces\ProductSpaceOperators.jl:186
 [7] union_by_union_rule(a::ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}, b::ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Space.jl:266
 [8] union(a::CosSpace{PeriodicSegment{Float64}, Float64}, b::CosSpace{PeriodicSegment{Float64}, Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Space.jl:279

Note that I already raised basically the same issue in #768 (about four months ago).

I would appreciate very much any idea how to overcome this error (e.g., by converting first F2 to period 8π ?). I was not able to find a simple way to perform this (in my trials I crashed Julia several times). Sorry for my ignorance!

Furthermore, I would be very happy if an automatic conversion feature would be supported in a future release of ApproxFun.jl to perform operations with periodic functions with commensurate periods.

andreasvarga · August 5, 2022, 11:29am

I think I got it. The following works:


julia> F1+Fun(t->F2(t),domain(F1))
Fun(2×2 ArraySpace:
ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,25.132741228718345❫) Fourier(【0.0,25.132741228718345❫); Fourier(【0.0,25.132741228718345❫) Fourier(【0.0,25.132741228718345❫)],[-3.561673020442948e-17, 5.072653133236335e-17, 3.497202527569243e-15, -2.720857218197605e-16, 3.497202527569243e-15, -2.720857218197605e-16, 3.18903060143427e-15, 2.8161797819056647e-15, -1.2246467991473535e-16, 1.224646799147353e-16  …  -4.73232564246473e-15, -2.565902665598115e-15, -3.481523346342104e-14, 1.1257283409025332e-15, -2.9565589116183405e-16, 1.790234627208065e-15, -1.519318030172918e-16, 1.790234627208065e-15, -1.519318030172918e-16, 2.000000000000014])

And

julia> F2+Fun(t->F1(t),domain(F2))
Fun(2×2 ArraySpace:
ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}[Fourier(【0.0,6.283185307179586❫) Fourier(【0.0,6.283185307179586❫); Fourier(【0.0,6.283185307179586❫) Fourier(【0.0,6.283185307179586❫)],[-1.294720671475603e-16, -2.0000000000000004, 2.4384757919829034e-15, -4.181839259808612e-15, 2.4384757919829034e-15, -4.181839259808612e-15, 1.6933257206756878e-14, -3.283993251760897e-14, 6.008022436061937e-16, 5.466932223891189e-17, 2.3058884893441655e-15, 3.493507550589483e-16, 2.3058884893441655e-15, 3.493507550589483e-16, 2.000000000000005])

works too.

DNF · August 5, 2022, 3:11pm

This anonymous function shouldn’t be necessary. You can just write F2, as in

F1+Fun(F2, domain(F1))

This a general observation about functions, nothing to do with ApproxFun.jl

andreasvarga · August 5, 2022, 4:04pm

Unfortunately, this doesn’t work:

julia> F1+Fun(F2, domain(F1))
ERROR: StackOverflowError:
Stacktrace:
 [1] map(f::ApproxFunBase.var"#288#289"{PeriodicSegment{Float64}}, v::Core.SimpleVector)
   @ Base .\essentials.jl:622
 [2] ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64}(dom::PeriodicSegment{Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Spaces\SumSpace.jl:45
 [3] SumSpace
   @ C:\Users\Andreas\.julia\packages\ApproxFunFourier\MPNwF\src\ApproxFunFourier.jl:405 [inlined]
 [4] canonicalspace
   @ C:\Users\Andreas\.julia\packages\ApproxFunFourier\MPNwF\src\ApproxFunFourier.jl:97 [inlined]
 [5] defaultcoefficients(f::Vector{Float64}, a::CosSpace{PeriodicSegment{Float64}, Float64}, b::ApproxFunBase.SumSpace{Tuple{CosSpace{PeriodicSegment{Float64}, Float64}, SinSpace{PeriodicSegment{Float64}, Float64}}, PeriodicSegment{Float64}, Float64})
   @ ApproxFunBase C:\Users\Andreas\.julia\packages\ApproxFunBase\Hasa4\src\Space.jl:332

DNF · August 5, 2022, 4:30pm

Hm. Is it because F2 is a “Fun”, not a function, and t->F2(t) turns it into one? Too bad, x->foo(x) is a bit of an eyesore.

jishnub · August 9, 2022, 5:52am

This requires defining the union of domains properly. I am working on a fix for this. It’s not fully ready yet, but will hopefully be included in ApproxFun soon.