Fix calcuation of n_elements in resizing limiter-specific containers

[bennibolm]

This fixes an error in the resizing of limiter-specific containers. Before, length(u_ode) was used for the new number of elements, although length(u_ode) = n_variables * n_nodes^n_dims * n_elements.
So, now it is possible to use AMR and initially nonconforming meshes for subcell limiting.

Nevertheless, respecting the bounds is not guaranteed for subcell IDP limiting at non-conforming subcell interfaces. Therefore, hanging nodes and mortars are ignored in the limiting process right now. However, often the simulation is still running nearly as good and with less computing effort.

[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.

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: 2 lines in your changes are missing coverage. Please review.

Comparison is base (0cf3e67) 95.46% compared to head (fafb9ef) 96.11%.
Report is 1 commits behind head on main.

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #1655      +/-   ##
==========================================
+ Coverage   95.46%   96.11%   +0.65%     
==========================================
  Files         418      418              
  Lines       34247    34248       +1     
==========================================
+ Hits        32692    32915     +223     
+ Misses       1555     1333     -222     

Flag	Coverage Δ
unittests	96.11% <0.00%> (+0.65%)	⬆️

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

Files	Coverage Δ
src/time_integration/methods_SSP.jl	83.87% <0.00%> (-0.91%)	⬇️

... and 9 files with indirect coverage changes

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

[sloede]

So, now it is possible to use AMR and initially nonconforming meshes for subcell limiting.

🎉

Nevertheless, respecting the bounds is not guaranteed for subcell IDP limiting at non-conforming subcell interfaces. Therefore, hanging nodes and mortars are ignored in the limiting process right now. However, often the simulation is still running nearly as good and with less computing effort.

I don't understand this, especially the "respecting the bounds is not guaranteed" part. If this may or may not produce correct results for refined meshes, I think it is too dangerous to activate. In that case it might produce unexpected errors or exhibit otherwise undefined behavior, which is extremely annoying to debug - then, I'd rather not enable it at all.

If it only means that subcell limiting is not as effective at refinement boundaries as in uniform boundaries but otherwise work correctly (including no out of bounds memory access), I'd leave it up to you to decide whether it is "good enough" to use. However, I still think a proper mortar support would be better - is that possible?

[andrewwinters5000]

Nevertheless, respecting the bounds is not guaranteed for subcell IDP limiting at non-conforming subcell interfaces.

So this means that the IDP on non-conforming meshes might violate something like positivity of density, correct? IIRC, this is also the case if one does not respect the (quite restrictive) time step criterion of the IDP schemes. Am I interpreting this correctly IDP crew @bennibolm @amrueda @gregorgassner ?

Is there any existing literature of IDP + non-conforming meshes? I did a cursory check on Nazarov's webpage and some quick searches but did not see anything.

[amrueda]

To some of your questions:

So this means that the IDP on non-conforming meshes might violate something like positivity of density, correct?

Yes. Unfortunately, the methods that we have developed and @bennibolm has implemented are not fully functional on non-conforming meshes. The main reason for this is that the implementation assumes that the low- and high-order schemes have the same surface terms. This is perfectly fine for conforming meshes, but not for non-conforming meshes. The low-order scheme with a high-order mortar is not invariant domain preserving.

As @bennibolm points out, the subcell IDP positivity limiter still works in many situations, where the low-order scheme + DG mortar produces "positive" solutions although the method is not really provably positivity preserving. Similar to when a time-step size higher than the IDP time-step size is used.

Other subcell limiters (e.g., the one that enforces local minimum and maximum principles) are less functional, as we currently do not have a strategy to define the bounds for nodes that sit at a non-confoming interface.

I still think a proper mortar support would be better - is that possible?

That would be great, but it is unfortunately not trivial. It would be nice to explore that in the future though!

In my opinion, the possibility to use subcell positivity limiting in non-conforming meshes is a nice feature to have, even though the method can still fail and the simulation crash. So I would add this feature, but I would suggest adding a warning message when the mesh is / can become non-conforming.

Is there any existing literature of IDP + non-conforming meshes?

I know of some literature in the context of algebraic flux correction for FEM. However, we have not experimented with any of these methods. Here some references:

• A. Jha and P. Knobloch, Adaptive Grids in the Context of Algebraic Stabilizations for Convection-Diffusion-Reaction Equations, arXiv:2007.08405v4 [math.NA] https://arxiv.org/abs/2007.08405
• J. Bonilla and S. Badia, Monotonicity-preserving finite element schemes with adaptive mesh refinement for hyperbolic problems. Journal of Computational Physics 416 (2020) 109522.
• J. Bonilla, Monotonicity-preserving finite element methods for hyperbolic problems. PhD thesis, UPC Barcelona, 2019. http://hdl.handle.net/10803/668242
• M. Bittl and D. Kuzmin, An hp-adaptive flux-corrected transport algorithm for continuous finite elements, Computing, 95 (2013), pp. 27-48. https://doi.org/10.1007/s00607-012-0223-y

[sloede]

Thanks for the clarifications, @amrueda! If this essentially means that IDP does work from an implementation point of view, I'm fine to merge it as is. I was just concerned about out of bounds memory access. Not being able to prove positivity on refined meshes (at least for now) is absolutely acceptable if it means we can continue to use this stabilization technique for a larger set of problems.

It would be good, however, if this were documented in one of the docstrings for the limiter (preferably the most high-level one), such that users are aware of this formal loss of positivity for refined meshes.

[andrewwinters5000]

Thanks for the detailed answer and references @amrueda !

[sloede]

It would be good, however, if this were documented in one of the docstrings for the limiter (preferably the most high-level one), such that users are aware of this formal loss of positivity for refined meshes

@bennibolm it would be great if you could put this into one of your next PRs 🙏

[bennibolm]

It would be good, however, if this were documented in one of the docstrings for the limiter (preferably the most high-level one), such that users are aware of this formal loss of positivity for refined meshes

@bennibolm it would be great if you could put this into one of your next PRs 🙏

Yes, sure. I already added a TODO to a PR.