Linear Semi-Lagrangian Vertical Transport

[James Kent (jameskent-metoffice)]

PR Summary

Sci/Tech Reviewer: Thomas Bendall (@tommbendall)
Code Reviewer: Oakley Brunt (@oakleybrunt)

This adds the semi-Lagrangian (SL) vertical transport to the linear model. The FFSL linear transport scheme has issues with potential temperature, partly because it is a Wtheta variable that FFSL puts on the shifted mesh. One way around this is to use the SL transport as this works on the native mesh.

The full SL vertical transport scheme solves the advection equation $dq/dt + w dq/dz = 0$ by interpolating the field to the departure point $q^{n+1}_a = q^{n}_d$, where $q$ is the field, $w$ is the vertical velocity, $n$ is the time index, and subscript $a$ and $d$ are the arrival and departure points respectively.

For the linear model we write the variables in terms of the linearisation state, $\bar{q}$ and $\bar{w}$, and the perturbation state, $q'$ and $w'$, and hence the advection equation now becomes

$$d\bar{q}/dt + dq'/dt + \bar{w} d\bar{q}/dz + \bar{w} dq'/dz + w' d\bar{q}/dz = 0$$

which, when removing the transport of $\bar{q}$ with $\bar{w}$, gives the advection equation for the perturbation field

$$dq'/dt + \bar{w} dq'/dz + w' d\bar{q}/dz = 0$$

The first terms can be computed using the standard SL scheme with the linearisation state wind, and the last term can be computed using the perturbation departure point (as this is the perturbation Courant number) multiplied by the gradient of $\bar{q}$ in the departure cell (note this is in the linearisation state departure cell)

$$q'^{n+1}_a = q'^{n}_d - \delta' \nabla \bar{q}$$

where $\delta' = w' \Delta{t} / \Delta{z}$ and $\nabla \bar{q}$ is computed using the difference of $\bar{q}$ values either side of the departure point.

Note that I did not start this branch from stable as I wanted to pick up the fixes to linear MoL in PR543.

Code Quality Checklist

• I have performed a self-review of my own code
• My code follows the project's style guidelines
• Comments have been included that aid understanding and enhance the readability of the code
• My changes generate no new warnings
• All automated checks in the CI pipeline have completed successfully

Testing

• I have tested this change locally, using the LFRic Apps rose-stem suite
• If any tests fail (rose-stem or CI) the reason is understood and acceptable (e.g. kgo changes)
• I have added tests to cover new functionality as appropriate (e.g. system tests, unit tests, etc.)
• Any new tests have been assigned an appropriate amount of compute resource and have been allocated to an appropriate testing group (i.e. the developer tests are for jobs which use a small amount of compute resource and complete in a matter of minutes)

There are kgo changes in the linear model app as I have changed the vertical transport for the FFSL tests to use the SL scheme (instead of MoL) for potential temperature.

Below shows the linear DCMIP gravity wave test with the vertical SL scheme for theta.

This can be compared with using MoL for theta (below).

trac.log

For some reason there are two lfricinputs failures. As my ticket does not touch any part of this code, I believe that it is unrelated to my PR.

Test Suite Results - lfric_apps - linear_vert_sl/run2

Suite Information

Item	Value
Suite Name	linear_vert_sl/run2
Suite User	james.kent
Workflow Start	2026-09-11T09:14:04
Groups Run	all

Dependency	Reference	Main Like
casim	MetOffice/casim@2026.07.1	True
jules	MetOffice/jules@2026.07.1	True
lfric_apps	jameskent-metoffice/lfric_apps@linear_vert_sl	False
lfric_core	MetOffice/lfric_core@60fc29a	True
moci	MetOffice/moci@2026.07.1	True
SimSys_Scripts	MetOffice/SimSys_Scripts@77a5166	True
socrates	MetOffice/socrates@2026.07.1	True
socrates-spectral	MetOffice/socrates-spectral@2026.07.1	True
ukca	MetOffice/ukca@9fc2b6d	True

Task Information

❌ failed tasks - 2

Task	State
run_lfricinputs_um2lfric-protogal_chem-N48L70_C48L70_ex1a_gnu_fast-debug-64bit	failed
run_lfricinputs_um2lfric-protogal_chem-N48L70_C48L70_ex1a_gnu_full-debug-64bit	failed

✅ succeeded tasks - 1597

⌛ waiting tasks - 1

Task	State
housekeep_ex1a	waiting

Security Considerations

• I have reviewed my changes for potential security issues
• Sensitive data is properly handled (if applicable)
• Authentication and authorisation are properly implemented (if applicable)

Performance Impact

• Performance of the code has been considered and, if applicable, suitable performance measurements have been conducted

AI Assistance and Attribution

• Some of the content of this change has been produced with the assistance of Generative AI tool name (e.g., Met Office Github Copilot Enterprise, Github Copilot Personal, ChatGPT GPT-4, etc) and I have followed the Simulation Systems AI policy (including attribution labels)

Documentation

• Where appropriate I have updated documentation related to this change and confirmed that it builds correctly

PSyclone Approval

• If you have edited any PSyclone-related code (e.g. PSyKAl-lite, Kernel interface, optimisation scripts, LFRic data structure code) then please contact the HPC Optimisation Team

Sci/Tech Review

• I understand this area of code and the changes being added
• The proposed changes correspond to the pull request description
• Documentation is sufficient (do documentation papers need updating)
• Sufficient testing has been completed

(Please alert the code reviewer via a tag when you have approved the SR)

Code Review

• All dependencies have been resolved
• Related Issues have been properly linked and addressed
• CLA compliance has been confirmed
• Code quality standards have been met
• Tests are adequate and have passed
• Documentation is complete and accurate
• Security considerations have been addressed
• Performance impact is acceptable

[github-actions]

⚠️ Hello James Kent (@jameskent-metoffice)!

Your CLA signature was found on the base branch, but you appear to have modified the CONTRIBUTORS.md file in this PR.

Please do not edit the CONTRIBUTORS.md file. If you have already signed the CLA, revert changes to the file and your signature will be picked up.

[DanStoneMO]

Testing with JEDI and it works fine. No linked PR needed for it

[cjohnson-pi]

Thanks James for adding this SL code. I just have a couple of comments from a linear code owner perspective.

[cjohnson-pi]

Please could you add ACTIVE/PASSIVE

[cjohnson-pi]

Note: In addition to this test, we also need to test the validity of the linearisation. I will allow this to be added on a separate ticket - and so I have added it to issue #507

[cjohnson-pi]

Please could the code be reordered so that PASSIVE fields (such as ls_field_local and grad_ls_field) are all computed first, followed by the ACTIVE fields. This will allow for easier implementation of the corresponding adjoint code.

[cjohnson-pi]

Please can you review the comments in the code to make sure that they actually reflect the code. For example, for a gradient I would expect a dividing by the distance between the points. But I think this must be incorporated into the departure distance number. I always get confused with the departure distance in LFRic - is it a distance, a velocity, a Courant number?

[cjohnson-pi]

Can you just note here why this has changed from true to false. Is it because the target is to have 5*false and this is a permanent change, or is it because the log-space code is still to be created and this is temporary.

[Thomas Bendall (tommbendall)]

I'm happy with the science and have just suggested fixing some typos and a couple of very trivial syntactical suggestions (which proves I've read the PR)

If you can address cjohnson-pi's comments, I'll be happy to approve

[Thomas Bendall (tommbendall)]

Suggested change

	!> @brief Kernel to compute the vertical cubic semi-Lagragian advection of a field
	!> @brief Kernel to compute the vertical cubic semi-Lagrangian advection of a field

[Thomas Bendall (tommbendall)]

Suggested change

	!! is solved using a cubic semi-Lagragian advection scheme. There are two
	!! is solved using a cubic semi-Lagrangian advection scheme. There are two

[Thomas Bendall (tommbendall)]

Suggested change

	! Compute pert dist based on whether this is W3 or W3theta field
	! Compute pert dist based on whether this is W3 or Wtheta field

[Thomas Bendall (tommbendall)]

Suggested change

	pert_dist(1:nlayers) = ( dep_dist_pert(w2_idx : w2_idx + nlayers - 1) &
	pert_dist(:) = ( dep_dist_pert(w2_idx : w2_idx + nlayers - 1) &

I think also works

[Thomas Bendall (tommbendall)]

Suggested change

	pert_dist(1:nl) = dep_dist_pert(w2_idx : w2_idx + nl - 1)
	pert_dist(:) = dep_dist_pert(w2_idx : w2_idx + nl - 1)

I think also works