Thermally perfect gas as an `AbstractEquationOfState` by DanielDoehring · Pull Request #3079 · trixi-framework/Trixi.jl

DanielDoehring · 2026-06-13T08:59:08Z

For high temperature flows, $c_p$ and $\gamma$ should no longer be treated as constant. In turn, they can modeled as functions of only temperature ("thermally perfect") while preserving the ideal gas relation $p = R_\text{spec} \rho T$ which makes them very attractive.
The go-to data that relates $c_p$ to $T$ is given by the NASA piecewise 9-degree polynomials.

This is an experimental implementation that implements thermally perfect gases as an equation of state.
However, I would like to not restrict this to NonIdealCompressibleEulerEquations since the original CompressibleEulerEquations could readily be extended with a non-constant gamma = gamma(T) . In fact, we already have this behaviour for the viscosity:

Trixi.jl/src/equations/compressible_navier_stokes_2d.jl

Lines 189 to 194 in 7091f3b

    
           # In the simplest cases, the user passed in `mu` or `mu()` 
        
           # (which returns just a constant) but 
        
           # more complex functions like Sutherland's law are possible. 
        
           # `dynamic_viscosity` is a helper function that handles both cases 
        
           # by dispatching on the type of `equations.mu`. 
        
           mu = dynamic_viscosity(u, equations)

Trixi.jl/src/equations/compressible_navier_stokes.jl

Lines 125 to 137 in 7091f3b

    
           """ 
        
               dynamic_viscosity(u, equations) 
        
           Wrapper for the dynamic viscosity that calls 
        
           `dynamic_viscosity(u, equations.mu, equations)`, which dispatches on the type of 
        
           `equations.mu`. 
        
           For constant `equations.mu`, i.e., `equations.mu` is of `Real`-type it is returned directly. 
        
           In all other cases, `equations.mu` is assumed to be a function with arguments 
        
           `u` and `equations` and is called with these arguments. 
        
           """ 
        
           dynamic_viscosity(u, equations) = dynamic_viscosity(u, equations.mu, equations) 
        
           dynamic_viscosity(u, mu::Real, equations) = mu 
        
           dynamic_viscosity(u, mu::T, equations) where {T} = mu(u, equations)

I also added varnames for cons2thermo (this is probably something for a different PR) to be able to save the thermal variables to get e.g. plots of temperature conveniently:

Comparison to ideal gas:
One can see that the shock is slightly further away from the cylinder compared to the real gas case.

github-actions · 2026-06-13T08:59:18Z

codecov · 2026-06-13T12:35:32Z

Codecov Report

❌ Patch coverage is 63.63636% with 44 lines in your changes missing coverage. Please review.
✅ Project coverage is 96.81%. Comparing base (1d411ba) to head (8832995).

Files with missing lines	Patch %	Lines
...elixir_euler_therm_perf_cylinder_bowshock_mach6.jl	0.00%	39 Missing ⚠️
...uations/equation_of_state_thermally_perfect_gas.jl	93.55%	4 Missing ⚠️
...e_1d_dgsem/elixir_euler_therm_perf_density_wave.jl	95.00%	1 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main    #3079      +/-   ##
==========================================
- Coverage   96.89%   96.81%   -0.08%     
==========================================
  Files         647      650       +3     
  Lines       50024    50145     +121     
==========================================
+ Hits        48466    48543      +77     
- Misses       1558     1602      +44

Flag	Coverage Δ
unittests	`96.81% <63.64%> (-0.08%)`	⬇️

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

☔ View full report in Codecov by Harness.
📢 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.

jlchan · 2026-06-13T13:42:31Z

This is great to see, @DanielDoehring; I was hoping to add this myself eventually.

Agreed that standard CompressibleEulerEquations could be extended with gamma(T), but the associated routines (EC fluxes, etc) may still assume gamma is constant.

What about adding type parameters to NonIdealCompressibleEulerEquations which identify whether the pressure EOS is non-ideal? This could differentiate between CPG, TPG, and fully non-ideal EOS?

… into ThermallyPerfGas

DanielDoehring · 2026-06-14T08:15:38Z

Agreed that standard CompressibleEulerEquations could be extended with gamma(T), but the associated routines (EC fluxes, etc) may still assume gamma is constant.

I noticed since there are flux_terashima_etal and flux_terashima_etal_central this might not be required since we can reuse a lot of machinery. 👍

What about adding type parameters to NonIdealCompressibleEulerEquations which identify whether the pressure EOS is non-ideal? This could differentiate between CPG, TPG, and fully non-ideal EOS?

We could think about that, does this change something in the implementation/allow for simplifications? I did not go through the entire equation of state and nonideal Euler equations, so maybe there are some spots where things could be more efficient for an ideal gas pressure/density/temperature relation.

jlchan · 2026-06-15T02:51:58Z

I think of thermally perfect as modifying the internal energy/temperature relationship, while a nonideal EOS modifies the pressure equation.

Specializing on thermally perfect this way might not make too much of a difference computationally, however.

jlchan · 2026-06-15T02:56:50Z

+    tol = 100 * eps(Float64)
+    dp = pressure(V, T, eos()) - p
+    iter = 1
+    while abs(dp) / abs(p) > tol && iter < 100
+        dp = pressure(V, T, eos()) - p
+        dpdT_V = ForwardDiff.derivative(T -> pressure(V, T, eos()), T)
+        T = max(tol, T - dp / dpdT_V)
+        iter += 1
+    end
+    if iter == 100
+        @warn "Solver for temperature(V, p) did not converge"
+    end


We use this fairly often for non ideal fluids; maybe it makes sense to factor this into a helper function? Happy to also do it if you file an issue

Agreed👍 - then we should also query eos as equations.eos or unpack it (this is why I added eos() as a function up there (albeit unnecessary convoluted)

Should we do the following: You @jlchan prepare this change across the repository, and I dedicate some time to adding documentation for the thermally perfect gas and the specific 9degree piecewise polynomial (following your suggestion below?)

Sure. To clarify, are you thinking of making a separate PR for AbstractThermallyPerfectGas <: AbstractEquationOfState first, or making a PR to this PR?

I would change this PR I think 👍

jlchan · 2026-06-15T03:00:04Z

+    # Fetch temperature interval index for a given temperature to select the correct polynomial coefficients.
+    # `clamp` ensures that even for temperatures outside the provided bounds,
+    # a valid index is returned to prevent a bounds error.
+    return clamp(searchsortedlast(eos.temperature_bounds, T), 1, N)


How efficient is searchsortedlast here?

knstmrd · 2026-06-15T14:05:04Z

A bit of an outsider's comment: the naming is a bit confusing (to me, at least), as the equations implement a specific (NASA 9th degree polynomial) model; so maybe something like ThermallyPerfectGas9PolyFit?

jlchan · 2026-06-15T14:14:26Z

A bit of an outsider's comment: the naming is a bit confusing (to me, at least), as the equations implement a specific (NASA 9th degree polynomial) model; so maybe something like ThermallyPerfectGas9PolyFit?

I drafted a comment on this last night but looks like it didn't send.

I think something like an AbstractThermallyPerfectGas <: AbstractEquationOfState with ThermallyPerfectGas9PolyFit <: AbstractThermallyPerfectGas would make sense, and allow one to just specialize on a smaller number of non-ideal EOS subroutines.

…ock_mach6.jl

DanielDoehring · 2026-06-15T19:29:59Z

A bit of an outsider's comment: the naming is a bit confusing (to me, at least), as the equations implement a specific (NASA 9th degree polynomial) model; so maybe something like ThermallyPerfectGas9PolyFit?

I drafted a comment on this last night but looks like it didn't send.

I think something like an AbstractThermallyPerfectGas <: AbstractEquationOfState with ThermallyPerfectGas9PolyFit <: AbstractThermallyPerfectGas would make sense, and allow one to just specialize on a smaller number of non-ideal EOS subroutines.

Fine with me 👍

… into ThermallyPerfGas

DanielDoehring · 2026-06-20T08:44:53Z

Should wait for #3093

… into ThermallyPerfGas

DanielDoehring added 4 commits May 28, 2026 22:36

start thermally perfect with other EoS

a6230ef

Merge branch 'main' into ThermallyPerfGas

2eeb97a

cont

f3de329

test and elixir

7dfd67b

DanielDoehring requested a review from jlchan June 13, 2026 08:59

DanielDoehring added enhancement New feature or request discussion labels Jun 13, 2026

DanielDoehring changed the title ~~Thermally perfect gas as Equation of state~~ Thermally perfect gas as AbstractEquationOfState Jun 13, 2026

DanielDoehring changed the title ~~Thermally perfect gas as AbstractEquationOfState~~ Thermally perfect gas as an AbstractEquationOfState Jun 13, 2026

DanielDoehring added 2 commits June 13, 2026 11:03

adaps

6be7088

Merge branch 'main' into ThermallyPerfGas

3e518b4

DanielDoehring added 7 commits June 13, 2026 16:49

coverage

c4c1482

Merge branch 'ThermallyPerfGas' of github.com:DanielDoehring/Trixi.jl…

7ca45cd

… into ThermallyPerfGas

Merge branch 'main' into ThermallyPerfGas

66e6d62

fmt

022450b

Merge branch 'ThermallyPerfGas' of github.com:DanielDoehring/Trixi.jl…

f50e4da

… into ThermallyPerfGas

fix

ca94fe3

bowshock

4439d16

DanielDoehring mentioned this pull request Jun 13, 2026

varnames cons2thermo #3084

Merged

DanielDoehring added 4 commits June 14, 2026 08:23

Merge branch 'main' into ThermallyPerfGas

7d80b08

tols

0d0b681

Merge branch 'ThermallyPerfGas' of github.com:DanielDoehring/Trixi.jl…

79289e1

… into ThermallyPerfGas

hll

28d40ed

Merge branch 'main' into ThermallyPerfGas

a93195b

jlchan reviewed Jun 15, 2026

View reviewed changes

DanielDoehring commented Jun 15, 2026

View reviewed changes

Comment thread examples/p4est_2d_dgsem/elixir_euler_therm_perf_cylinder_bowshock_mach6.jl

Update examples/p4est_2d_dgsem/elixir_euler_therm_perf_cylinder_bowsh…

6e5d189

…ock_mach6.jl

DanielDoehring added 3 commits June 15, 2026 21:30

Merge branch 'main' into ThermallyPerfGas

8832995

adapt structure

489cfff

Merge branch 'ThermallyPerfGas' of github.com:DanielDoehring/Trixi.jl…

a3646e6

… into ThermallyPerfGas

DanielDoehring mentioned this pull request Jun 20, 2026

temperature as function of V and p for nonideal Eq of State #3093

Open

DanielDoehring added 2 commits June 20, 2026 09:07

fix examples

0e51188

Merge branch 'main' into ThermallyPerfGas

7f97480

DanielDoehring added 4 commits June 20, 2026 10:59

store gamma

ef3a35d

Merge branch 'ThermallyPerfGas' of github.com:DanielDoehring/Trixi.jl…

773c0fe

… into ThermallyPerfGas

rm copy paste

b3981a8

rm

c7f96bc

	# In the simplest cases, the user passed in `mu` or `mu()`
	# (which returns just a constant) but
	# more complex functions like Sutherland's law are possible.
	# `dynamic_viscosity` is a helper function that handles both cases
	# by dispatching on the type of `equations.mu`.
	mu = dynamic_viscosity(u, equations)

	"""
	dynamic_viscosity(u, equations)

	Wrapper for the dynamic viscosity that calls
	`dynamic_viscosity(u, equations.mu, equations)`, which dispatches on the type of
	`equations.mu`.
	For constant `equations.mu`, i.e., `equations.mu` is of `Real`-type it is returned directly.
	In all other cases, `equations.mu` is assumed to be a function with arguments
	`u` and `equations` and is called with these arguments.
	"""
	dynamic_viscosity(u, equations) = dynamic_viscosity(u, equations.mu, equations)
	dynamic_viscosity(u, mu::Real, equations) = mu
	dynamic_viscosity(u, mu::T, equations) where {T} = mu(u, equations)

Conversation

DanielDoehring commented Jun 13, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

github-actions Bot commented Jun 13, 2026

Review checklist

Purpose and scope

Code quality

Documentation

Testing

Performance

Verification

Uh oh!

codecov Bot commented Jun 13, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Codecov Report

Uh oh!

jlchan commented Jun 13, 2026

Uh oh!

DanielDoehring commented Jun 14, 2026

Uh oh!

jlchan commented Jun 15, 2026

Uh oh!

Uh oh!

jlchan Jun 15, 2026

Choose a reason for hiding this comment

Uh oh!

DanielDoehring Jun 15, 2026

Choose a reason for hiding this comment

Uh oh!

DanielDoehring Jun 16, 2026 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jlchan Jun 16, 2026

Choose a reason for hiding this comment

Uh oh!

DanielDoehring Jun 17, 2026

Choose a reason for hiding this comment

Uh oh!

jlchan Jun 15, 2026

Choose a reason for hiding this comment

Uh oh!

knstmrd commented Jun 15, 2026

Uh oh!

jlchan commented Jun 15, 2026

Uh oh!

Uh oh!

DanielDoehring commented Jun 15, 2026

Uh oh!

DanielDoehring commented Jun 20, 2026

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

DanielDoehring commented Jun 13, 2026 •

edited

Loading

codecov Bot commented Jun 13, 2026 •

edited

Loading

DanielDoehring Jun 16, 2026 •

edited

Loading