GRHD: Make parameters of cons2prim's velocity filter configurable (!140)

This adds the following parameter to the `GRHD` section:
- `c2p_threshold_velocity_limiter`: Threshold `y₀ = ρhW²` below which the velocity is limited by a filtered division following Dumbser+ 2024.

I found emperically that this value strongly influences the reflection and generation of perturbations in the TOV star tests.
In particular, using 1e-4 gives a PSD of the central density evolution in which all peak frequencies are wrong compared to literature (e.g. Three-dimensional numerical general relativistic hydrodynamics. II. Long-term dynamics of single relativistic stars, Font+, 2002, table 1).
Using instead a value of 1e-6 gives better results.
However, this seems to come at the expense of stability, because with lower limiter values we appear to more easily apply surface heating, e.g. epsilon or tau display artifacts at the surface.

The default value of `c2p_threshold_velocity_limiter` is now set to `1e-6` because of the above observation.
This also explains why updated the ref test in this PR.

We did not account for the fact that there a two
data points on an interface as well as two neighbors
(one in each direction).

This adds new parameters to the `Mesh` section:
- `filter_enable`: Enable spectral filtering to be applied after each full time step.
- `filter_parameters`: Parameters for exponential spectral filter, cf. (5.16) in Hesthaven, Warburton 2007.

GRHD: also impose EoS limits when computing final primitives (!137)

This went unnoticed so far, because we have been running only with polytropes and there epsilon is basically unused.

Evolve: Fix t argument for callback_substep for ode_step (!136)

Readd previously removed `load_parameters` (!135)

GRHD: Implement TOV test in 2d axisymmetry using the cartoon formulae (!134, !134)

The `examples/grhd_tov/dg_cartoon.toml` example can run for a good while without immediately crashing. This means that the bdry conditions (reflection symmetry across `x = 0` and `z = 0` axes) aren't so wrong and the surface treatment is doing something.

However, ... I realized that I have been running with smaller time steps than other schemes do. That is, mine scales as `~ N^{-3}`, but it should more be like `~ N^{-2}`.
This explains why I could run with comparatively large CFL numbers some of the 1d tests.
This should be investigated in a separate PR.
Changing to `~N^{-2}` here makes the cartoon test crash shortly after start, because of an issue in the cons2prim recovery. So perhaps we can't get away without any regularization scheme and I have to implement AV too.

Some fixes to dg1d.jl core (!133, !133)

- Adds `setup_env(parfile)` with which you can now directly construct `env, project, bdryconds, prms` from a parfile. These objects are the basic ingredients for an evolution.
- Fixes overspecialization in caching code.
- Fixes a bug where showing `CallbackTrigger` caused an error.

GRHD: Make `cons2prim` work with 2D velocities (!132, !132)

I am glad I have the reftests in place, because this was very delicate to implement correctly.
Part of this struggle was also to not introduce extra allocations from boxing of the closures, due to Julia sharing variables between scopes ...

This PR also applies a general clean up to the `cons2prim` function.

`GRHD`: Implement reflection bdry conditions for `valencia1d, spherical1d, rescaled_spherical1d` formulations (!131, !131)

It works nicely with `spherical1d, rescaled_spherical1d`, but fails miserably with `valencia1d`.
Not sure what the problem there is, because for both `rescaled_spherical1d` and `valencia1d` I need to interpolate onto the origin after unscaling the conservatives,
because the volume factor comes with an extra `r` that makes the conservatives vanish there.
However, for `valencia1d` and `spherical1d` I also have singular source terms, which I don't have in `rescaled_spherical1d`.

`@with_signature`: avoid tuple limit by using `get_variable` unpacking directly in function body (!130, !130)

The implementation of the `doublecartoon` formulation in `GRHD` finally hit the tuple length limit, above which tuples are stack allocated and, thus, incur a drastic performance hit due to extra allocs.

Here we update `@with_signature` so that we no longer have to pass arguments through tuples,
but instead we look them up with the `get_variable` method from `mesh.cache`.
This way we also do the unpacking now directly inside the macro-generated function body.

Its not an ideal solution, because now `@with_signature` can no longer be used without the `get_variable` method from `dg1d`, which increases coupling.
But for the moment I don't bother, and I shall worry about that whenever I decided to move `@with_signature` into a separate pkg, if ever.

---

The `get_variable` method and `Cache` could also use a bit of benchmarking and optimizations.
E.g. atm we perform three lookups
- `isregistered` to check if the variable is available
- `Cache.variables_to_fields` to determine in which variable group the variable can be found
- `getproperty(::Variable, ::Symbol)` to get the variable

Clearly, this could be improved, but we leave that for the future.

---

Obviously, I also updated the tests. But I refused to setup now `Cache`s for each test case.
Instead I just added a dispatch `get_variable(::NamedTuple, ::Symbol)` and a `cache` argument to the expanded function so that the tests didn't so much adjustments.

GRHD: Fix rescaled_spherical1d formulation (!128, !128)

A few issues have silently accumulated because we did not add a reftest when we introduced this formulation.
This PR fixes them and adds a reftest.

GRHD: Add reftest for doublecartoon formulation (!129, !129)

GRHD: Revert changes from !114 regarding spherical1d formulation (!127, !127)

Although we had a reftest added for the TOV test with `spherical1d`, we had not enabled the comparison of data. As such the test was almost useless. This is fixed here for this test. But we should add more reftests for other formulations in separate PRs.

GRHD Fixup atmosphere freezing logic (!126, !126)

We need to account for the fact that global solution is degenerated on cell interfaces.
Because of that we need to cast our data stencil wider when we check for resets on a bdry.

GRHD: Implement reflection bdry conditions for doublecartoon formulation and TOV test. (!125, !125)

I decided to implement the symmetry conditions now differently than I did for the other formulations (`spherical1d, rescaled_spherical1d, valencia1d`).
That is, I now compute the bdry values of the state variable separately before I call `bdryllf_...`.
The values are computed such that
- for variables that are even wrt `r -> -r` trafos, the x derivative of the nodal approximation at the origin vanishes,
- for variables that are odd wrt `r -> -r` trafos, the value of the nodal approximation at the origin vanishes.
Given this data (stored in `bdry_D, bdry_Sx, bdry_tau`), I then impose the bdry conditions weakly by using the LLF flux and the previously computed
bdry values as exterior states.
I also tested computing the exterior states according to `q^+ = -q^- + 2 q_BC`, where `q_BC` are the previously computed values, but did not find that
this has significant effect.

What is interesting is that this implementation (doublecartoon + reflection bdry conds.) can now be evolved without symmetrizing the state vector explicitly
and it appears to be stable, e.g. no blowup at the origin.
I wonder whether adapting this way of imposing bdry conditions for other formulations would also help with stability, although I am not sure anymore
if they suffered from blowups in the reflected case.
Note that all other formulations suffer from the problem of division by zero, which we dynamically fix up by explicit symmetrization. Imposing the bdry
conditions then differently won't fix that issue.

However, for now we keep the bdry conditions in the other formulations as they currently are.

---

More thoughts about bdry values:
- Like I said above, all other formulations use the values from initial data as bdry values, which makes sense because the solution is supposed to be static.
One should check if using instead bdry values such that the symmetry conditions are satisfied actually works better and/or is factually more accurate.
The latter might be plausible, because numerical dissipation (espcially with AV) might just cause a decrease in the central density.
Using instead the initial data as bdry values would try to keep the central density as it is, despite the star solution being altered.
- For the mirrored domain test we can have a grid setup in which the cells faces fall onto the origin. It would make sense to test if using similarly
constructed interface values for the LLF flux at the interfaces that lie at the origin. Maybe that could also cure the stability issue we see in
all the mirrored runs.

* update parameter docstrings

* add minimal description for each formulation

* remove dead code

* split equation.jl into spherical1d.jl, rescaled_spherical1d.jl, valencia1d.jl, doublecartoon.jl

* manually undo the insertion of the refvol factor into spherical1d
formulation

the formulae without the refvol factor should already be correct
in particular, they correspond to the choice where all three
metrics (physical one, conformal one and reference one) all agree
I will document this properly in a follow up PR were I also
disentangle all formulations

* rename llf and bdryllf to refer to spherical1d formulation

* update example parfile

* fix derivatives of induced metric in initial data

* remove unused variables

* fix sign in 4 metric determinant

* fixup minus and fix doublecartoon corrections at x=0

* fix todo conflict

* some todos

* implement symmetrizer around origin, but disable it for now

* fix computation of state variables in llf_doublecartoon

* some debug statements

* implement the doublecartoon correction

* implement the doublecartoon correction

* remove not evolved y,z components and fixup nomenclature

* add doublecartoon implementation

* factor out data reader from initial data setup

* update parfile

* add more tov examples

* fix spherical1d to account for reference volume element

* fixup impose_symmetry! for valencia1d

* reset tau when it falls below 0

this is needed to run with lower atmosphere densities

* impose symmetry on state variables before c2p

* callbacks: fix atm freezing logic

* c2p: todos

* implement valancia1d version that works with tov test

* add grhd_bondi example that uses valencia1d formulation

* implement valencia1d version

* add valencia1d option for formulation parameter

* c2p: more verbose log output

* register more variables for rescaled_spherical1d

* fix dispatch in setup.jl

* verify formulation parameter

* add bdryllf for rescaled_spherical1d

* stop gap for unhandled cases of singularities in impose_symmetry functions

* add simulation time to progress meter output

* implement cartoon derivative expressions for mixed rank (1,1) tensor

* add docstrings

* fix imports

* add cartoon macros

* update 2d reftests to check 2d output

also reduced number of output variables and
times in order to save storage

* IntegrationTests: error when reference data is missing

* fix filename of 2d output

* also test 0d output with fv_muscl_burgers_sine reftest

* fix u_norm2 computation for FVElement

* reftests: update advection_sine_analysis

* implement u_norm2 0d analysis variable

* undo method redefintion

* reftests: update fv_muscl_burgers_sine

note: fv_muscl_burgers_sine gives an asymmetric result when we use
upwinding (muscl_omega=1.0) which is probably wrong

* fixup update_test script

* remove dead code

* utilize limit_slope! method in ScalarEq/rhs.jl

* Implement integration for FVElement meshes with slope limited data

* extract slope limiter from ScalarEq/rhs.jl FV methods

fixup hrsc prms

* remove dead code

* reftest: update fv_muscl_burgers_sine to use aligned_ts for outputs

* bump version

* reftests: replace old output.h5 files with output0d/1d/2d.h5 files

* reftests: update all parfiles to use the new output parameters

* teach IntegrationTests how to deal with new output files

* update unit tests

* SaveCallbacks: handle file format properly

* Let SaveCallback append the file extensions

* for now restrict substep and interpolating output to 1d and h5 format

remove substep_format1d, substep_format2d parameters

* queue in all SaveCallbacks

* only use variables0d/1d/2d/ to determine if that output is enabled

* SaveCallback: adjust terminal value logic for absent format specifier

* Adjust Output parameters

- remove enabled0d/1d/2d
- only output if a non-empty string value of format0d/format1d/format2d
  is provided
- replace variable with variable0d/1d/2d so that output is split into
  three files

* Teach TimeAlignedSaveCallback to adjust time steps
of multiple SaveCallbacks at once

* define Base.append overload to work with AbstractVectors

* type annotate kwargs for SaveCallback

* bump ci julia version to 1.10

* update manifest

* ignore experiments folder

* compat Jacobi

* fix invalidations by upgrading to Jacobi@v0.6.2

requires to dev my fork here:
https://github.com/fatteneder/Jacobi.jl/tree/fa/bump-polynomials

* add advection_sine_analysis as ref test

* hack: also output mesh when 0d output is requested

* add checks for limited SaveCallback functionality

* implement 0d output

* add 0d output parameters

* update advection_sine ref test

* implement analyze_error computation

* add now all parameters to Project type

* add proper register! dispatching

* add new parameters

* wip: stash some tov examples

* update grhd ref tests to use c2p_freeze_atm_reset = false to restore previous behavior

* init c2p_freeze_atm

* fixup logic of freezing atm

* introduce c2p_freeze_atm flag which allows to define when an atmosphere can propagate

* fix time step computation to use maximum of smoothed_mu

* rescale smoothed viscosity by the extrema of D

* add more example parfiles for tov case

* enable c2p_freeze_atm_reset option

* add cons2prim_rescaled_spherical1d_freeze_flags

* add c2p_freeze_atm_reset option

* document perturbation_rho parameter

* update examples

* experiment with partially frozen c2p flags

* don't bernstein in AV computation

* add c2p with partially frozen reset flags

* use buffers for logarithmic AV indicators

* remove isinvalid

* remove dead code

* avoid splat constructor from CallbackSet

* add punctuation to some parameter doc strings

* add more colors to highlight in which files comparisons failed

* IntegrationTests: always use terminal coloring to highlight results

* update update_test script

* include substep outputting in some reftests

* teach IntegrationTests about substep_output.h5

* add unit test

* move substep savecallback logic into convert_to_remote_savecallback!

only add remote savecallback when requested

* enable substep output

* introduce AbstractCallbackTrigger and CallbackTrigger

* add parameters

* remove experiments folder

* grhd_tov: add more slope limiter examples

* add fv_tvb_mirrored.toml TOV example

* add reftest fv_superbee_grhd_bondi_accretion

* implement rhs_slope_limiter!

* need init_ variables also for FV because of bdry conditions

* no need to compute fluxes when setting up initial data

* fixup callbacks

* add fv_grhd_bondi_accretion reftest

* a separate register! for FV

* disambiguate rhs! dispatch

* SRHD: remove unused variables from fv method

* dispatch on mesh for register!

* splice in fv method dispatch

* add new parameters

* SRHD: add reftest srhd_simple_wave_fv_tvb

* SRHD: Implement slope limiter methods

* SRHD: clean up c2p

* SRHD: use Brent's method in c2p because A42 is currently bugged

* only init relevant bdry variables

* register more vars

* wrap existing FV method as rhs_fv_central

* thin out internal rhs! usage

* organize rhs.jl

* a registration for new hrsc method

* rename _register_variables! to register!

* add new parameters

* EulerEq: add more ref tests for FV slope limiters

* EulerEq: rename hrsc option 'muscl' to 'slope_limiter' and implement
more limiters

* TCI: implement maxmod

* EulerEq: compute max speed for each substep

* add ref test fv_muscl_euler_isentropic_flow

* fix doc string of LLF

* TCI: make minmod non-allocating

* log invalid dt

* EulerEq: make muscl scheme run

* EulerEq: PoC muscl scheme

* EulerEq: run callbacks, initialdata only when certain variables are registered

* EulerEq: Add LV dependency

* evolve: init stages before capturing them

* EulerEq: strip down rhs! args

* rename _register_variables! to register! in EulerEq

* add new parameters to EulerEq project

* fix fn name

* remove unused file

* add reftest fv_muscl_burgers_sine

* fixup muscl implementation; also make it work with arbitrary fluxes

* fix parameter check

* rename prm hrsc_omega to muscl_omega

* add ref test for muscl scheme

* implement muscl scheme

* implement new parameters

* fix reftest output coloring when non interactive

* add TODO branches for bdry conds

* fixup check for coordinate signularity on grid

* implement different bdry conditions for bondi and tov problem

also update bondi reftest

* more example parfiles

* c2p: fix computation of W

* make training wheel ID work independent of method

* fix ID interpolation for FV method

* cleanup fv_rhs.jl

* skip symmetry imposition if signularity is not on grid

* mirror TOV initial data for mirrored domains

* add bernstein tov example

* implement bernsteining

* add missing cache wrapper for bdry variables

* reconstruct dynamic variables after viscosity computation

* add perturbation_rho parameter

* add a D dependency to viscosity

* allow to only evolve initial data in star's interior

* add rho h W^2 variable for analysis

* make initial viscosity work

* allow slightly negative energies

* pull out y0 from vr limiter in cons2prim

* add training_wheels options which smooths initial dynamical variables

* update example parfiles

* update tov examples

* wip: tinker with c2p and viscosity computation in atmospheric regions

* filter AV from atmosphere cells

* enable inflow condition on inner boundary

* fix imposition of EoS bounds in cons2prim

also disabling lower bracket adjustment, because it seems to be wrong
with this the FV method can now run

* cleanup cons2prim

* wip: tinkering with cons2prim

* compute AV for initial data

* implement av_drag option

* wip: tinker with c2p routine

* add c2p analysis vars

* stabilize origin of TOV with Neumann bdry conditions

* implement entropy_variables_rescaled_spherical1d

* use log(D) for MDA viscosity because it is more sensitive

* fix maxspeed_rescaled_spherical1d

* example parfile that allows to evolve the tov for quite a while

* apply velocity limiter from dumbser at the end of cons2prim

* fix maxspeed

* add update_test helper script

* add ref test balancedburgers_kink_avmda

* fix time step computation for AV method

* simplify inf clipping

* mda: use local instead of global speed to compute baseline av

* favor likely case in if in mda computation

* implement av dragging

* add av_drag parameter

* spice up IntergrationTests logging with colors

* fix 0th order AV smoother

* implement BalancedBurgers equation

* fix fv stepping to also add dt to source term

* cleanup

* rm old ERK evolutions

* fixup fv ref tests

* rename rk steppers to ode_step!

* for FV methods forward dt from evolution via cache

fixes mismatch in dt values that were used in the time stepper and the
rhs functions
the mismatch was caused by the time stepper always applying the CFL
(and output callbacks could further shorten the proposed time step)
eventually, this lead to wrong propagation speeds, which can be seen
by the test results of fv_advection_sine (we expect same profile as
initial data after two revolutions, but we actually saw 2.5)

only add dt for FV evolutions

* fixup kwarg forwarding in register_variables

* implement a minimal git dirty logging

* restructure variable registration

* implement rhs for av regularization

* update callbacks

* implement ldg and av regularization

* mesh: forward more methods to cache

* dispatch TODO error

* forward Cache methods to Mesh

* nospecialize setup and initialdata

* update examples once more

* remove comments

* cleanup: remove codegen for now

* IntegartionTests: add grhd_tov_spherical1d ref test

* IntegrationTests: remove reltol from grhd_bondi_accretion test

* spherical1d: fix missing shift from maxspeed computation

* cleanup

* make spherical1d formulation work when r=0 is included by removing potential nans in source terms by enforcing symmetry

* make impose_symmetry! dispatch on project

* add tov example that runs with spherical1d formulation

* wip: change maxspeed implementation for tov star initial data

* cons2prim: use non-zero values for rho, eps thresholds

* fix: set atm density in bondi_accretion initial data

* implement tov initial data for spherical1d formulation

* add formulation parameter

* implement rhs for FVElement and rescaled_spherical1d

* implement equations for FVElement and rescaled_spherical1d

* mesh: implement differentiate for FVElement

* add grhd example par files

* update notebooks

* fix missing lapse factor in source terms

* fix typos

* fixup initialdata setup

* derivative computation needs jacobi

* rename cons2prim routines

* add atmosphering to cons2prim

* implement rhs! for rescaled_spherical1d formulation

* specialize timespeed computation based on equations

* add specialized equations for the rescaled GRHD version

* implement initialdata setup for tov case

* spherical accretion:make sure K, Gamma agree with loaded data

* expose spectral differentiation

* wrap SphericalAccretion solver into module; recompute bondi_accretion.h5 with less radial points

* refactor TOV solver to use OrdinaryDiffEq

* initialdata: cleanup dir

* add bondi accretion as ref test

* fix up loading of id file

* fv_rhs: implement fv_update_step! with source term

* add ROOTDIR const global so that parfiles can specify paths relative to project root

* update SphericalAccretion.jl

* update sympy nb

* cleanup initialdata setup

* equations: fix up Tmunu computation

* Revert "fix source term in dg and fv rhs computation"

This reverts commit 161cc034.

* Replace Interpolations.jl with BasicInterpolators.jl

* equation: add missing flux terms to sources; use outflow bdry condition on inner radius

* comments in cons2prim

* implement momentum constraint for debugging

* test SphericalAccretion

- copy Bugner's C version of the solver
- compare my results with Bugner's -- good agreement
- compare against literature results (Papadopoulos 1998,1999) --
  discrepancy in sign of vr, ur, otherwise good agreement

* SphericalAccretion: fix wrong definition of uc^r and Vc

* tweak impl more

* fix source term in dg and fv rhs computation

* fix equations and implement fv method

* equations: combine flux and source computation

* math: update spherical_accretion notebook

* SRHD: export module

* re-implemented GRHD eqs with a non-zero shift

* initialdata: add grhd_source.codegen.jl that was generated with TensorComponents.jl

* deps: need Interpolations.jl for initial data setup

* update math/spherical_accretion.py

* update setup

* start initialdata setup for bondi problem

* rework rhs.jl

* update deps and parameters

* update type defs

* initialdata: add HDF5 dep; write h5 output for SphericalAccretion.jl

* rename cons2prim.jl

* remove unused files

@with_signature: only need one LineNumberNode per @accepts, @returns unpacking block (!100)

* @with_signature: only need one LineNumberNode per @accepts, @returns unpacking block