Skip to content

Conversation

@benegee
Copy link
Collaborator

@benegee benegee commented Oct 23, 2025

Implemented following Souza et al..

@benegee
Copy link
Collaborator Author

benegee commented Oct 23, 2025

@MarcoArtiano Maaybe some things still need to be adapted to your equations.

@MarcoArtiano
Copy link
Collaborator

@MarcoArtiano Maaybe some things still need to be adapted to your equations.

Thanks! I think the source term might need an adjustment. I'll take a look at that.

Copy link
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 Benedict! I've added a few comments/suggestions.

du3 = -k_v * (u[3] - dotprod * x[2])
du4 = -k_v * (u[4] - dotprod * x[3])

du5 = -k_T * u[1] * c_v * (temperature - T_equi)
Copy link
Collaborator

Choose a reason for hiding this comment

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

To make this source term consistent with the potential temperature evolution we can use:

Suggested change
du5 = -k_T * u[1] * c_v * (temperature - T_equi)
du5 = -k_T * u[5] / temperature * (temperature - T_equi)

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I will blindly trust you, but I could not retrace this step immediately.

Copy link
Collaborator

@MarcoArtiano MarcoArtiano Nov 3, 2025

Choose a reason for hiding this comment

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

I can see two ways to treat that term in terms of potential temperature evolution. The first one would be to derive the evolution of the potential temperature and see how the source term would get transformed. If you do that you would get:
$$s_{\varrho \theta} = \frac{\partial \varrho \theta}{\partial p} (\gamma - 1) \left (V \cdot s_{\varrho V} + s_{\varrho E} \right ) $$.

The other approach would be simply to keep the source term as $$s_{\varrho \theta} = K \Delta T$$, where $K$ has to be consistent with the dimensions of [\varrho \theta]/(sec [T]), and that's why I choose $$s_{\varrho \theta} = -k_T \frac{\varrho \theta}{T} \Delta T$$.

I will run both versions to be sure they are consistent and stable.


lat_lon_trees_per_dim = 10
layers = 8
mesh = Trixi.T8codeMeshCubedSphere(lat_lon_trees_per_dim, layers, EARTH_RADIUS, 30000.0,
Copy link
Collaborator

Choose a reason for hiding this comment

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

This is just out of curiosity, since I've never used T8code, why this one instead of P4est?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Just because t8code was used in my current project (because others have linked t8code to their codes).

solver = DGSEM(polydeg = polydeg, surface_flux = surface_flux,
volume_integral = VolumeIntegralFluxDifferencing(volume_flux))

lat_lon_trees_per_dim = 10
Copy link
Collaborator

Choose a reason for hiding this comment

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

Since this is already extremely costly, could we use (8,4) elements? What do you think?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Sure!
I just used this setting in the initial version because this is what is used in the paper by Souza et al..


###############################################################################
# ODE solvers, callbacks etc.
T = 365 # 1 year
Copy link
Collaborator

@MarcoArtiano MarcoArtiano Oct 28, 2025

Choose a reason for hiding this comment

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

The test case seems to be really expensive to make it run for a full year. Could we lower that to 10 days as we do for the baroclinic instability or anything that you might think is appropriate, yet not that long as 1 year?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Yes, of course.
Again, this was because of the reference (although I think it is even longer in the paper). AFAIK this test case is usually run for a very long time and then evaluated using averages in time and space.

333.2487183585763
],
tspan=(0.0, 0.01 * SECONDS_PER_DAY),
lat_lon_trees_per_dim=2, layers=2)
Copy link
Collaborator

Choose a reason for hiding this comment

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

probably we would need to decrease abstol and reltol of the time integrator to avoid CI failing on MacOS and/or Windows. That worked perfectly fine for me here: test Euler Energy

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

The current values were just copied and are not useful. Once we have decided on the resolution and time span to be used in the CI run, we can adjust them.

@benegee benegee marked this pull request as ready for review October 28, 2025 16:17
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.

3 participants