ELBDM Spectral Interpolation  vivi235711/gamer GitHub Wiki
This section provides an overview of the spectral interpolation method in GAMER implemented via the GramFourier continuation algorithm in the Jupyter notebook GramFE.ipynb
located in tool/table_maker/GramFE
. This algorithm is essential for computing extension and interpolation tables.

CompileTime Flag: Ensure GAMER is compiled with
SUPPORT_SPECTRAL_INT
(or configured withspectral_interpolation
). 
Runtime Parameters: Set
OPT__FLU_INT_SCHEME
andOPT__REF_FLU_INT_SCHEME
to8
for enabling spectral interpolation. SetSPEC_INT_TABLE_PATH
to the directory containinginterpolation_tables
andboundary2extension_tables
. EnableSPEC_INT_XY_INSTEAD_DEPHA
to interpolate x = density^0.5*cos( phase/SPEC_INT_WAVELENGTH_MAGNIFIER ), y = density^0.5*sin( phase/SPEC_INT_WAVELENGTH_MAGNIFIER ) instead of density and phase, which has the advantage of being welldefined across vortices.SPEC_INT_WAVELENGTH_MAGNIFIER
is the stretching factor of wavelength; setting it to unity gives x=real part and y=imaginary part.

Script: Use
download_spectral_interpolation_tables.sh
located inexample/test_problem/ELBDM/LSS_Hybrid
.

Script: Use
compute_interpolation_tables.py
located intool/table_maker/GramFE
. 
Execution: Run
mpirun n 16 python3 compute_interpolation_tables.py
for table generation.  Time Consumption: The script may take several hours to execute without additional output.
 Spectral Accuracy: Offers highprecision interpolation results.
 Flexibility: Suitable for various interpolation needs within the GAMER framework.
 Nonconservative and Nonmonotonic: Interpolation results may not always preserve these properties.

Computationally Intensive: Spectral interpolation is slower than polynomial interpolation with a local stencil (i.e. interpolation options
1=MinMod3D, 2=MinMod1D, 3=vanLeer, 4=CQuad, 5=Quad, 6=CQuar, 7=Quar
).
Different algorithms in GAMER introduce different discretisation errors:
 numerical PDE solvers
 interpolation (new patches during refinement and ghost boundaries for neighbouring patches on different AMR levels)
 time interpolation for flexible timesteps (firstorder)
 restriction operation (firstorder)
 fluxfixup operation to ensure mass conservation
To minimise these errors and measure the performance of different interpolation algorithms, you may consider the following settings
 Use vortex pair test problems

Compile time parameters:
ELBDM_SCHEME=ELBDM_WAVE
,WAVE_SCHEME=WAVE_GRAMFE
andGRAMFE_SCHEME == GRAMFE_MAMTUL

Runtime parameters: Turn
OPT_FIXUP_FLUX
,OPT_INIT_RESTRICT
,OPT_FIXUP_RESTRICT
,OPT_PHASE_INT
off, setOPT__DT_LEVEL = 1
 Problem Description: The spectral interpolation method may introduce negative density values. This issue arises because the interpolation is nonmonotonic.

Consequences: Negative density can lead to NaN (Not a Number) values in both density and flux calculations within the fluid scheme. This problem has been observed to potentially cause
MaxPot == 0.0
errors in certain simulations.  Current Investigation: The issue is under active investigation, with more details and ongoing discussions available on GitHub and Slack. Specific instances of this problem have been shared here.

Proposed Solution: A preliminary solution involves applying a density floor when preparing density ghost zones in
Flu_Prepare()
. This adjustment could prevent the generation of negative density values during spectral interpolation. Further details and implementation suggestions can be found in the related GitHub comment thread.
 Users encountering similar issues or observing unexpected results in density and flux calculations should use GAMER's other more reliable interpolation algorithms (i.e.
6=CQuar
(conservative quartic interpolation)).
 The spectral interpolation method in GAMER is implemented through a GramFourier continuation algorithm in the Jupyter notebook
GramFE.ipynb
, located intool/table_maker/GramFE
. The basic idea is to find a periodic extension of a given interpolant. This periodic function can then be interpolated using a DFT. This interpolation is not monotonic and not conservative, but highly accurate.
 Lyon and Bruno's Original Work:

Arbitrary Precision Arithmetic: Utilizes
mpmath
library.  GramSchmidt Orthogonalization: Implements this algorithm for basis computation.
 SVD Algorithm: Used for highprecision continuations of Gram polynomials.
 FFT Computations: Involves forward and backward FFTs in high precision.
The Python code generates two types of tables essential for the spectral interpolation algorithm, stored in interpolation_tables
and boundary2extension_tables
.

Precomputed Operation: The entire interpolation process (
GramFE extension * DFT * shift in k space * IDFT
) is precomputed as a single matrix. 
Storage: The resulting matrices are stored in doubleprecision binary files named
N.bin
, whereN
is the input interpolation size. 
Efficiency: This approach is efficient for small values of
N
, requiringN^2
operations for each matrix multiplication at runtime. 
Equivalent to Polynomial Interpolation: Notably, for
m = N
, this method is equivalent to polynomial interpolation of orderN
, as confirmed by numerical tests.
 Runtime Computation: For large inputs, the extension and DFTs are computed at runtime using the FFT algorithm.

Table Utilization: The tables
%ld_%ld_%d_%d.bin
are used, wherend
,m
,Gamma
,g
represent the table parameters. By applying these tables to the boundary data (size2 * m
), an extension of sizend
is produced. 
Adaptive Size Selection: The size
nd
is adaptively chosen during runtime to optimize the FFT algorithm's performance, favoring sizes with small prime factorizations. 
Matrix and FFT Operations: The interpolation involves
2 * m * m
matrix operations and additionalO(N log N)
operations for the FFT, managed by the FFTW library.
 Determining Factors: The interpolation's accuracy heavily depends on the quality of the periodic extension.

Optimal Parameters: Parameters
Gamma = 150
,g = 63
,nd = 32
, andm = 8
strike a balance between stability and accuracy. 
Artifacts and Behavior: Higher values of
m
can lead to interpolation artifacts, while lower values might reduce accuracy but result in smoother interpolation.
The following plots show a number of accuracy comparisons between quartic interpolation and spectral interpolation:
The last function demonstrates that while the highorder boundary polynomials used in the boundary matching of the GramFE algorithm suffer from the Runge phenomenon, for larger N the numerical properties of the algorithm are determined by the plane wave basis and convergence sets in. Potentially, the large error of the GramFE method in this case can be remedied by switching to lowerorder boundary polynomials for small N at the cost of sacrificing accuracy for more wellbehaved functions.
The following plots show a performance comparison between quartic interpolation and spectral interpolation using matrix multiplication and FFT. It demonstrates that for