Skip to content

WIP: Add non-conservative 2D Euler equations with gravity potential#79

Draft
amrueda wants to merge 29 commits intomainfrom
arr/euler_gravity_geopotential
Draft

WIP: Add non-conservative 2D Euler equations with gravity potential#79
amrueda wants to merge 29 commits intomainfrom
arr/euler_gravity_geopotential

Conversation

@amrueda
Copy link
Copy Markdown
Collaborator

@amrueda amrueda commented Mar 18, 2025

This set of equations extends the CompressibleEulerEquations2D by incorporating geopotential energy into the total energy rho_e. This ensures that the energy equation remains a conservation law, even in the presence of gravitational effects. Additionally, the source terms in the momentum equations are treated as non-conservative terms rather than point-wise source terms, following:

  • Souza, A. N., He, J., Bischoff, T., Waruszewski, M., Novak, L., Barra, V., ... & Schneider, T. (2023). The flux‐differencing discontinuous galerkin method applied to an idealized fully compressible nonhydrostatic dry atmosphere. Journal of Advances in Modeling Earth Systems, 15(4), e2022MS003527. https://doi.org/10.1029/2022MS003527.
  • Waruszewski, M., Kozdon, J. E., Wilcox, L. C., Gibson, T. H., & Giraldo, F. X. (2022). Entropy stable discontinuous Galerkin methods for balance laws in non-conservative form: Applications to the Euler equations with gravity. Journal of Computational Physics, 468, 111507. https://doi.org/10.1016/j.jcp.2022.111507.

We implement the non-conservative term flux_nonconservative_waruszewski, which provides well-balancedness (to machine precision accuracy) for isothermal hydrostatic equilibria:
image

Other standard benchmark tests are provided as well:

RisingBubble_afterBugFix_lmars.mp4

image

@amrueda amrueda marked this pull request as draft March 18, 2025 13:28
@codecov
Copy link
Copy Markdown

codecov bot commented Mar 18, 2025

Codecov Report

❌ Patch coverage is 84.77612% with 51 lines in your changes missing coverage. Please review.
✅ Project coverage is 96.31%. Comparing base (aa057e7) to head (2825d7e).
⚠️ Report is 12 commits behind head on main.

Files with missing lines Patch % Lines
src/equations/compressible_euler_gravity_2d.jl 84.77% 51 Missing ⚠️
Additional details and impacted files
@@            Coverage Diff             @@
##             main      #79      +/-   ##
==========================================
- Coverage   97.08%   96.31%   -0.78%     
==========================================
  Files          37       38       +1     
  Lines        5010     5345     +335     
==========================================
+ Hits         4864     5148     +284     
- Misses        146      197      +51     

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

Co-authored-by: Benedict Geihe <bgeihe@uni-koeln.de>
@amrueda amrueda force-pushed the arr/euler_gravity_geopotential branch from e5b5534 to 57ebc3b Compare March 18, 2025 13:42
@tristanmontoya
Copy link
Copy Markdown
Member

@amrueda What still needs to be done for this PR to be ready?

Copy link
Copy Markdown
Collaborator

@MarcoArtiano MarcoArtiano left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks, I left a few comments.

f0, f0)
end

function flux_nonconservative_waruszewski(u_ll, u_rr, orientation::Integer,
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do we need to support StructuredMesh?

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For simple cases, it's nice to be able to test with the TreeMesh or StructuredMesh solver :)

# Calculate 2D flux for a single point
@inline function Trixi.flux(u, orientation::Integer,
equations::CompressibleEulerEquationsWithGravity2D)
rho, rho_v1, rho_v2, rho_etot, phi = u
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't remember if we already discussed about that. I'm sorry.

Since the equations already mentions that this is with gravity, I would naturally assume that the rho_e is indeed the total energy, without the need to specify.


Should be used together with [`UnstructuredMesh2D`](@ref).
"""
@inline function Trixi.boundary_condition_slip_wall(u_inner,
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In atmospheric flows I think we could avoid this a bit more expensive boundary condition and use the other more simple already implemented.

I'm note sure if actually there are differences and one has an advantage over the other.

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have added this BC because it's very robust and has been shown to be entropy stable. See, e.g.,

Hindenlang, F. J., Gassner, G. J., & Kopriva, D. A. (2020, August). Stability of wall boundary condition procedures for discontinuous Galerkin spectral element approximations of the compressible Euler equations. In Spectral and High Order Methods for Partial Differential Equations ICOSAHOM 2018: Selected Papers from the ICOSAHOM Conference, London, UK, July 9-13, 2018 (pp. 3-19). Cham: Springer International Publishing.

The classical reflective slip-wall BC that uses the numerical flux function is only provably ES for some Riemann solvers, such as LLF, HLL, HLLC (see paper). However, there's no guarantee that it will also be ES for other Riemann solvers. I don't know if LMARS provides entropy stability, for example.

You are right that this is probably not critical in most atmospheric science cases, but the additional robustness might be important in some cases. For instance, when simulating the blast from a volcano eruption ;)

I don't think this BC is intrinsically more expensive than the reflective boundary condition that uses the surface numerical flux function. In particular, the numerical flux function needs to compute the maximum wave speed and pressure, which is essentially what is needed for this BC.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants