Add unstable_check for custom SSP method

[bennibolm]

This PR adds an instability check for the custom SSPRK time integration method.

Currently, when running an unstable simulation (e.g. sedov blast (elixir_euler_sedov_blast_wave_sc_subcell.jl without any limiting) the simulation is already unstable after the second time step

────────────────────────────────────────────────────────────────────────────────────────────────────
 Simulation running 'CompressibleEulerEquations2D' with DGSEM(polydeg=3)
────────────────────────────────────────────────────────────────────────────────────────────────────
 #timesteps:                  2                run time:       3.13681926e-01 s
 Δt:             2.33398471e-03                └── GC time:    0.00000000e+00 s (0.000%)
 sim. time:      5.57222448e-03 (0.557%)       time/DOF/rhs!:  1.80197861e-07 s
                                               PID:            2.37261007e-07 s
 #DOFs per field:         16384                alloc'd memory:        313.587 MiB
 #elements:                1024

 Variable:       rho              rho_v1           rho_v2           rho_e         
 L2 error:             NaN         NaN         NaN         NaN
 Linf error:     1.36283654e-01   4.22437117e-01   4.22437117e-01   6.11449577e-01
 ∑∂S/∂U ⋅ Uₜ :         NaN
────────────────────────────────────────────────────────────────────────────────────────────────────

Nevertheless, the simulation continues to time step 5:

────────────────────────────────────────────────────────────────────────────────────────────────────
 Simulation running 'CompressibleEulerEquations2D' with DGSEM(polydeg=3)
────────────────────────────────────────────────────────────────────────────────────────────────────
 #timesteps:                  5                run time:       8.05128560e-02 s
 Δt:             9.88719304e-01                └── GC time:    0.00000000e+00 s (0.000%)
 sim. time:      1.00000000e+00 (100.000%)     time/DOF/rhs!:  5.24908752e-08 s
                                               PID:            5.65288900e-08 s
 #DOFs per field:         16384                alloc'd memory:        252.749 MiB
 #elements:                1024

 Variable:       rho              rho_v1           rho_v2           rho_e         
 L2 error:             NaN         NaN         NaN         NaN
 Linf error:     2.22044605e-16   3.79470760e-19   3.79470760e-19   6.77626358e-21
 ∑∂S/∂U ⋅ Uₜ :         NaN
────────────────────────────────────────────────────────────────────────────────────────────────────

When using the SSPRK33 method of OrdinaryDiffEqSSPRK, we get the following warning after the second time step and the simulation terminates afterward.

┌ Warning: Instability detected. Aborting
└ @ SciMLBase ~/.julia/packages/SciMLBase/rvXrA/src/integrator_interface.jl:631

from here.

I want to mimic that with this simple instability check.

[github-actions]

Review checklist

This checklist is meant to assist creators of PRs (to let them know what reviewers will typically look for) and reviewers (to guide them in a structured review process). Items do not need to be checked explicitly for a PR to be eligible for merging.

Purpose and scope

• The PR has a single goal that is clear from the PR title and/or description.
• All code changes represent a single set of modifications that logically belong together.
• No more than 500 lines of code are changed or there is no obvious way to split the PR into multiple PRs.

Code quality

• The code can be understood easily.
• Newly introduced names for variables etc. are self-descriptive and consistent with existing naming conventions.
• There are no redundancies that can be removed by simple modularization/refactoring.
• There are no leftover debug statements or commented code sections.
• The code adheres to our conventions and style guide, and to the Julia guidelines.

Documentation

• New functions and types are documented with a docstring or top-level comment.
• Relevant publications are referenced in docstrings (see example for formatting).
• Inline comments are used to document longer or unusual code sections.
• Comments describe intent ("why?") and not just functionality ("what?").
• If the PR introduces a significant change or new feature, it is documented in NEWS.md with its PR number.

Testing

• The PR passes all tests.
• New or modified lines of code are covered by tests.
• New or modified tests run in less then 10 seconds.

Performance

• There are no type instabilities or memory allocations in performance-critical parts.
• If the PR intent is to improve performance, before/after time measurements are posted in the PR.

Verification

• The correctness of the code was verified using appropriate tests.
• If new equations/methods are added, a convergence test has been run and the results
  are posted in the PR.

Created with ❤️ by the Trixi.jl community.

[codecov]

Codecov Report

Attention: Patch coverage is 75.00000% with 3 lines in your changes missing coverage. Please review.

Project coverage is 96.21%. Comparing base (cba534b) to head (6e1fbb7).
Report is 1 commits behind head on main.

Files with missing lines	Patch %	Lines
src/time_integration/methods_SSP.jl	66.67%	2 Missing ⚠️
src/time_integration/time_integration.jl	83.33%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #2445      +/-   ##
==========================================
+ Coverage   94.86%   96.21%   +1.35%     
==========================================
  Files         506      509       +3     
  Lines       41952    42116     +164     
==========================================
+ Hits        39795    40519     +724     
+ Misses       2157     1597     -560     

Flag	Coverage Δ
unittests	96.21% <75.00%> (+1.35%)	⬆️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[DanielDoehring]

Is used correctly in the custom time integrators, i.e., over the ODE.jl implementation?

In other words, no explicit example needed?

[bennibolm]

Is used correctly in the custom time integrators, i.e., over the ODE.jl implementation?

In other words, no explicit example needed?

I'm not really sure what you mean.

Thanks for the reference. I didn't know that ode_unstable_check exists. Unfortunately, in my example dt isn't NaN but only part of the solution itself, so the simulation would still continue.

[bennibolm]

I'm also thinking about adding this to the other custom time integration methods. What do you think?

[ranocha]

It would make sense to unify this. However, the next question is how "correct" it should be.

• all(isfinite, u) will not work with MPI. That's why we have our own ode_unstable_check.
• Shall we have the same flexibility as OrdinaryDiffEq.jl and pass the unstable check function to the options of solve?

[DanielDoehring]

Hm yeah so I also though about this at some point (see the PR) but disregarded this as I did ultimately not care too much if we would take a fistful more timesteps before breaking down.

[ranocha]

If we export it, it becomes part of our public API (although it is not documented...). However, I do not think this is a good idea right now since it can be easily confused with ode_unstable_check. I see two options right now:

• Do something specific for your subcell limiter ecosystem
• Make it general so that it works nicely with the other parts

[bennibolm]

Okay, that's true.
I'm completely fine with the first option. But how would that look like in your opinion? Rename it to something like unstable_check_SSPRK, have an optional parameter unstable_check = nothing for solve() of the custom time integrators and add unstable_check = unstable_check_SSRK to every elixir, where it's needed?
Or really removing the parameter unstable_check from the solve function and hardcode it into the SSPRK method?