Skip to content

Conversation

@nhartney
Copy link
Contributor

@nhartney nhartney commented Dec 4, 2025

This PR adds the capability to evaluate physics schemes in the inner loop of the semi-implicit quasi-Newton time stepping scheme. Physics is evaluated both explicitly and implicitly, in a similar way to forcing. The parameter physics_beta controls the semi-implicit off-centring, where a physics_beta of 1 is fully implicit and a physics_beta of 0 is fully explicit.

The PR also changes the names of the physics schemes lists, replacing the pre-existing physics_schemes with final_physics_schemes.

Finally, the pull request introduces a moist shallow water linear solver, MoistThermalSWSolver. This includes moisture in the linear solve stage and is necessary when doing physics in the inner loop. The solver is only suitable for use in a saturated atmosphere where the saturation function has the usual dependence on depth and buoyancy.

These changes are tested with a new inner loop physics test that uses the moist solver.

@tommbendall tommbendall added enhancement Involves adding a new capability equation Adding or enhancing a new equation labels Dec 12, 2025
Copy link
Contributor

@tommbendall tommbendall left a comment

Choose a reason for hiding this comment

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

Thanks for all of this Nell. I like the tests. There's a few things I'm suggesting tidying-up:

  • some tweaks to the solver so that every term has a trial function
  • removing the copy_active_tracers list

Note that the new solver will be reorganised by Alex in #648 !

x_out: The output set of fields
"""

for name in self.non_solver_prognostics:
Copy link
Contributor

Choose a reason for hiding this comment

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

I think we want to keep this as non_solver_prognostics

- beta_u * 0.5 * Dbar * bbar * div(w) * dx
- beta_u * 0.5 * Dbar * b * div(w) * dx
- beta_u * 0.5 * bbar * div(w*(D-Dbar)) * dx
+ beta_u * 0.5 * jump((D-Dbar)*w, n) * avg(bbar) * dS
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
+ beta_u * 0.5 * jump((D-Dbar)*w, n) * avg(bbar) * dS
+ beta_u * 0.5 * jump((D)*w, n) * avg(bbar) * dS

Alex pointed out to me that terms only involve reference values and no trial functions shouldn't be here

inner(w, (u - u_in)) * dx
- beta_u * (D - Dbar) * div(w*bbar) * dx
+ beta_u * jump(w*bbar, n) * avg(D-Dbar) * dS
- beta_u * 0.5 * Dbar * bbar * div(w) * dx
Copy link
Contributor

Choose a reason for hiding this comment

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

Alex pointed out to me that terms only involve reference values and no trial functions shouldn't be here

Suggested change
- beta_u * 0.5 * Dbar * bbar * div(w) * dx

eqn = (
# u equation
inner(w, (u - u_in)) * dx
- beta_u * (D - Dbar) * div(w*bbar) * dx
Copy link
Contributor

Choose a reason for hiding this comment

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

Alex pointed out to me that terms only involve reference values and no trial functions shouldn't be here

Suggested change
- beta_u * (D - Dbar) * div(w*bbar) * dx
- beta_u * (D) * div(w*bbar) * dx

# u equation
inner(w, (u - u_in)) * dx
- beta_u * (D - Dbar) * div(w*bbar) * dx
+ beta_u * jump(w*bbar, n) * avg(D-Dbar) * dS
Copy link
Contributor

Choose a reason for hiding this comment

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

Alex pointed out to me that terms only involve reference values and no trial functions shouldn't be here

Suggested change
+ beta_u * jump(w*bbar, n) * avg(D-Dbar) * dS
+ beta_u * jump(w*bbar, n) * avg(D) * dS

+ beta_u * jump(w*bbar, n) * avg(D-Dbar) * dS
- beta_u * 0.5 * Dbar * bbar * div(w) * dx
- beta_u * 0.5 * Dbar * b * div(w) * dx
- beta_u * 0.5 * bbar * div(w*(D-Dbar)) * dx
Copy link
Contributor

Choose a reason for hiding this comment

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

Alex pointed out to me that terms only involve reference values and no trial functions shouldn't be here

Suggested change
- beta_u * 0.5 * bbar * div(w*(D-Dbar)) * dx


# Check that the third field is buoyancy
if not equation.field_names[2] == 'b':
raise NotImplementedError("Field 'b'must exist to use the moist thermal linear solver in the SIQN scheme")
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
raise NotImplementedError("Field 'b'must exist to use the moist thermal linear solver in the SIQN scheme")
raise NotImplementedError("Field 'b' must exist to use the moist thermal linear solver in the SIQN scheme")

self.setup_transporting_velocity(aux_scheme)

# Set up tracers to copy (those not included in the linear solve)
self.tracers_to_copy = []
Copy link
Contributor

Choose a reason for hiding this comment

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

So I think I removed this list at some point -- could we remove all of these changes to add it back in?

I added a solver_prognostics argument to the SIQN stepper -- I think when we use the MoistSWLinearSolver, we should be using that argument and passing all of our variables

@tommbendall tommbendall merged commit f3a89b1 into main Dec 17, 2025
5 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancement Involves adding a new capability equation Adding or enhancing a new equation

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants