the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 10 Jan 2024
Split-explicit external mode solver in finite volume sea ice ocean model FESOM2
Abstract. A novel split-explicit (SE) external mode solver for the Finite volumE Sea ice–Ocean Model (FESOM2) and its sub-versions (example 2.5) is presented. It is compared with the semi-implicit (SI) solver currently used in FESOM2. The split-explicit solver utilizes a dissipative asynchronous (forward–backward) time-stepping scheme. Its implementation with Arbitrary Lagrangian-Eulerian vertical coordinates like Z-star (Z ∗) and Z-tilde (Z˜) is explored. The comparisons are performed through multiple test cases involving idealized and realistic global simulations. The SE solver demonstrates lower phase errors and dissipation, but maintain a simulated mean ocean state very similar to the SI solver. The SE solver is also shown to possess better run-time performance and parallel scalability across all tested workloads.
- Preprint
(1711 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
-
RC1: 'Comment on gmd-2023-208', Mark R. Petersen, 27 Feb 2024
Thank you, it is helpful for the ocean model development community to see your success with a split explicit method. This paper takes great care in the detailed presentation of the numerical methods for the new split explicit solver in FESOM2.
What order of convergence do you expect for your time-stepping method? Please test the convergence rate in time. The test case would need to have both slow and fast dynamics to be a useful test, but you could test them first separately.
You could test JUST the barotropic part with the Kelvin Wave or Inertial-Gravity wave tests here in Bishnu et al. (2024). Just use redundant layers in the vertical. Then you can compare against an exact solution for the convergence test.
You could then use the baroclinic channel test case from Ilicak et al. (2012) and Petersen et al (2015), refine in time, and compute convergence by comparing to a short-timestep case. The duration of the simulation would need to be short (30 minutes, say) and sufficiently laminar that the results do not differ due to small differences in the turbulent flow. This test case definitely has baroclinic dynamics. It is not designed for barotropic dynamics, but should include some surface gravity waves as the SSH and barotropic eastward flow adjust to a geostrophic balance. Thanks to Ange Ishimwe (Université catholique de Louvain) for pointing out this test case during his talk at AGU Ocean Sciences. He was able to show second-order convergence for his split explicit scheme using this method in his recent paper published in December, in Ishimwe et al. (2023) figure 6. This would be a good paper to reference for comparison in the current work.
Smaller comments:
- You subcycle the external mode to exactly one baroclinic time step, and do not use a filter, as explained in lines 30-34. I think this is important enough to put in the abstract, as it is a potential 2x speed-up for the barotropic stage compared to models that subcycle to n+2
- On equation 8, first line, I believe the sign of gH grad eta should be positive.
- The surface flux W was dropped in equation 8. It would be good to just comment that W was added to the baroclinic (eqn 14), and not the barotropic mode, and your reasoning for that. Because you have a function at the top with the delta in equation 14, it makes sense that this is a baroclinic addition.
- Fig 3 caption m2 needs superscript
I appreciate your description of the reasoning behind your choices. For example, the discussion of whether to add the bottom drag term to the barotropic dynamics at line 95, and the discussion on abandoning filtering at line 110. I have considered these exact issues and it is good to hear your thoughts on these.
Thank you,
Mark R. Petersen, Los Alamos National LaboratoryBishnu, S., Petersen, M.R., Quaife, B., Schoonover, J.A. A Verification Suite of Test Cases for the Barotropic Solver of Ocean Models. Submitted to JAMES. https://essopenarchive.org/users/566154/articles/612943-a-verification-suite-of-test-cases-for-the-barotropic-solver-of-ocean-models
Ishimwe, A.P., Deleersnijder, E., Legat, V., Lambrechts, J., 2023. A split-explicit second order Runge–Kutta method for solving 3D hydrodynamic equations. Ocean Modelling 186, 102273. https://doi.org/10.1016/j.ocemod.2023.102273
Ilicak, A.J. Adcroft, S.M. Griffies, R.W. Hallberg. Spurious dianeutral mixing and the role of momentum closure Ocean Modell., 45 (2012), pp. 37-58, 10.1016/j.ocemod.2011.10.003
Petersen, M.R., D.W. Jacobsen, T.D. Ringler, M.W. Hecht, M.E. Maltrud, Evaluation of the arbitrary Lagrangian–Eulerian vertical coordinate method in the MPAS-Ocean model, Ocean Modelling, Volume 86, February 2015, Pages 93-113, ISSN 1463-5003, http://dx.doi.org/10.1016/j.ocemod.2014.12.004.
Citation: https://doi.org/10.5194/gmd-2023-208-RC1 -
CC1: 'Comment on gmd-2023-208', Ange Ishimwe, 29 Feb 2024
I have two main remarks.
At line 79, it is stated that the Coriolis term is not included in the Baroclinic equations and is computed in the Barotropic equations, as shown in Equation 6. However, my concern is that the Coriolis term is calculated only using the Barotropic velocity (U_bar in Equation 6). There is no correction or adjustment in the Baroclinic part. This implies that the Coriolis force is considered uniform vertically, which is inaccurate. For instance, in an Ekman spiral, the Coriolis force changes direction vertically.
My second concern focus to bottom stress. At line 81, it is mentioned that frictional force is treated explicitly. Typically, this significantly impacts the time step, reducing a lot its maximum value. Why not consider a semi-implicit/Patankar approach, where the bottom drag term is treated as C * norm(u_explicit) * u_implicit? This approach is quite common.
In the same vein, how is the bottom drag coefficient C calculated? Is it constant or does it vary with bathymetry? What are its values in the two presented test cases? I am impressed by the baroclinic time steps, especially with explicit friction considered. If there is no bottom drag, would it be possible to show a vertical profile of velocity at certain locations?
Thank you.
Citation: https://doi.org/10.5194/gmd-2023-208-CC1 -
AC1: 'Reply on CC1', Tridib Banerjee, 01 Mar 2024
Hi Ange,
Thank you for your questions. I would like to provide the following clarifications.
About the coriolis, we also account for the coriolis in the baroclinic part of our solution. If you look at equation 3 (line 70), we already have coriolis taken into account for our baroclinic estimate (U*). The only reason equation 5 (line 75) doesnt include it is because it is the RHS of the barotropic subsystem where we are explicitly solving coriolis at every sub-step.
About the bottom drag, yes we can follow as you suggested also. In reality, since for oceans we have much larger bottom layers, unlike coastal applications, it should not be a big factor. Same with the coefficient of drag. In FESOM, it is a constant C_d=0.0025. We can do other sophisticated formulations but typically, this suffices for ocean applications.
Thanks again for asking and hope this helps with clarifying your doubts. Please dont hesitate to ask if anything else comes to mind.
Citation: https://doi.org/10.5194/gmd-2023-208-AC1 -
CC2: 'Reply on AC1', Ange P. Ishimwe, 01 Mar 2024
Thank you very much for your respond.
I have no more question qbout the drag coefficient.
Can I ask why you have to compute the Coriolis term in the barotropic mode ? Putting this term only in the baroclinic mode should be enough, at least for the stability of the scheme. What is the reasoning for that ?
Ange
Citation: https://doi.org/10.5194/gmd-2023-208-CC2 -
AC2: 'Reply on CC2', Tridib Banerjee, 04 Mar 2024
Hi,
It is from our desire to have accurate dispersion relationship for inertia gravity waves and simulate barotropic response without phase errors. Indeed, Coriolis will not be so important stability wise as the courant number will rather be limited by (sqrt(gHK^2)dt) than (fdt). You can see similar discussion in A.F. Shchepetkin, J.C. McWilliams / Ocean Modelling 9 (2005) 347–404, around equation 3.48.
Hope this helps. Please ask again if anything else is unclear.
Best, Tridib
Citation: https://doi.org/10.5194/gmd-2023-208-AC2
-
AC2: 'Reply on CC2', Tridib Banerjee, 04 Mar 2024
-
CC2: 'Reply on AC1', Ange P. Ishimwe, 01 Mar 2024
-
AC1: 'Reply on CC1', Tridib Banerjee, 01 Mar 2024
-
RC2: 'Comment on gmd-2023-208', Anonymous Referee #2, 26 Mar 2024
Review
The authors describe the implementation of a split-explicit algorithm for the FESOM2 ocean model, using for time integration of the external mode the modified FB scheme or the modified AB3-AM4 scheme proposed by Demange et al. The implementation is described for the z-star and z-tilde vertical coordinates. Numerical integration tests are carried out for the z-star split-explicit model, and compared with the current semi-implicit algorithm of FESOM2. These tests show numerical solutions that are qualitatively close to those obtained with the semi-implicit version. They also show better scalability, particularly for highly parallelized workflows.
General comments :
The paper addresses the interesting question of how to construct a split-explicit algorithm based on the FESOM time-stepping scheme, which has the peculiarity of staggering the prognostic variables in time. Original questions arise about how to phase the time integration of 3D variables, momentum and tracers, and the time integration of the barotropic mode. However, it seems to me that the paper needs to be clarified on several points and requires major revisions.
Section 2 details the proposed split-explicit algorithm. Some clarifications seem necessary to me concerning (i) the choice made for the order of integration of the variables Ubar and eta with the modified FB scheme and the modified AB3-AM4 scheme and the implications of this choice, (ii) the initialisation of the variable Ubar for the barotropic integration, (iii) the correction of the 3d variables after completion of the barotropic integration.
Section 3 details the amplitude and phase errors of the modified FB scheme, the modified AB3-AM4 scheme and the semi-implicit FESOM scheme in the context of the SW equations. This analysis, although largely based on published results by Demange et al., could nicely illustrate the contrast between the errors produced by explicit schemes and the semi-implicit FESOM scheme. The presentation and the figure should however be improved. It seems to me that sections 2 and 3 should then be swapped, so that the SW-based results of the current section 3 introduce/motivate the implementation of the split-explicit algorithm of current section 2.
Section 4 reports the results of numerical experiments. It is somewhat disturbing that the solutions for the idealised case show little difference between the split-explicit algorithm and the semi-implicit algorithm. This is probably due to the fact that the dynamics of this idealised case is essentially baroclinic and little affected by the representation of the external dynamics. It would have been desirable to consider an idealised case where barotropic dynamics play a more important role (for example, think of an idealised case of internal tide generation on a bathymetry by barotropic current oscillation). However, the results show a certain viability of the split-explicit algorithm and a gain in efficiency compared to the semi-implicit algorithm.
Specific comments :
l37 : I suppose ‘temporal interpolations’ should be replaced by ‘temporal integrations’?
l47 : A schematic describing the time-staggering of variables and the structure of the algorithm could be referred to here, upstream of the equations, to help the reader.
l53 : Should the variable h in the expression for the pressure gradient be indexed with k?
l62 : Should the transport U be indexed with k?
l100 : It seems to me that the abbreviation SE isn't very well chosen. It would be clearer to use another which suggests that the basic time-stepping is the modified FB. For the same reason,
l101-102 : Here, with the SE scheme, the choice is made to integrate U_bar before eta. With the SESM scheme, the choice will be reversed. What are the reasons for these different choices? What are the implications? In the SE case, this choice is associated with a semi-implicit discretisation of the Coriolis term, which probably requires specific numerical processing ; if so, these should be explained. On the other hand, still in the SE case, this choice leads to the expression of the mean barotropic transport (12), which would be different with a reverse order of integration.
l103 : The quoted value of theta has been obtained under rather special conditions: it is the value that, for a constant stratification of N=10-2 s-1, a water column of H=4000 m and a splitting ratio of 20, allows the stabilisation of split-explicit algorithms whose stability is limited by large-scale barotropic waves. It is not clear that the algorithm proposed in this paper is constrained by these waves and that this value is relevant here. This is plausible, however, because the FESOM scheme is a variant of the FB scheme, which is used by Demange et al. in their study.
l110 : It might be clearer to change ‘at the end of the barotropic step’ into ‘during the barotropic step’.
l114 : The reference to Shchepetkin and McWilliams on the line 112 could be moved to just after the ‘traditional notation’ line 114.
l120 : It is mentioned that it is not clear how to initialise Ubar (due, I suppose, to the time-staggering between U and Ubar), and it is noted that "We return to this topic below". But it's not clear to me where in the paper this is explained.
l125 : In this expression we see the correction by the mean barotropic transport <> of the 3D transports U_k^{n+1/2}, which is then used to transport tracers. There doesn't seem to be any place in your implementation for the other transport correction that classically appears in split-explicit algorithms based on synchronous schemes, i.e the correction by the barotropic transport of the 3D transports U_k^{n+1} (see for example p 384 in Shchepetkin and McWilliams). I'm a little embarrassed by the lack of an equivalent to this correction. Isn't there a cog missing to ensure consistency between the 3D and 2D integrations?
l156 : Should ‘temporal interpolations’ be replaced by ‘temporal integrations’?
l157-159 : It seems to me that this sentence is confusing. What is done in this section is not really the analysis of the external mode in the context of split-explicit algorithms (which is done, for example, in Demange et al. section 4.2). Rather, it is the analysis of the modified FB, modified AB3-AM4 and semi-implicit FESOM schemes in the context of SW equations.
l190 : The solutions of the continuous problem are e^{+i c} and e^{-i c}.
l 214 : The sentence ‘This CFL number …mesh size.’ is unnecessary here because the analysis is exact in space and does not depend on Delta x.
fig. 1: For the two panels on the left, why not show the amplitude error of the explicit schemes over larger ranges on each of the two axes? This would show the damping for large values of c, and the loss of stability of the schemes. This would make it a little easier to compare the amplitude error of the explicit schemes with that of the semi-implicit scheme. For the two panels on the right, the phase error could be plotted as the ratio of the discrete phase velocity to the exact phase velocity. why not use the same horizontal axis as for the left panels, with c values varying between 0 and 30? This would make the figure easier to read.
fig. 3 : the caption refers to diffusivities that are not traced.
Citation: https://doi.org/10.5194/gmd-2023-208-RC2
Model code and software
FESOM2.5 with preliminary Split-Explicit Subcycling Banerjee Tridib, Danilov Sergey, Scholz Patrick, Klingbeil Knut, and Sidorenko Dimitry https://zenodo.org/doi/10.5281/zenodo.10040943
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 268 | 66 | 29 | 363 | 8 | 12 |
- HTML: 268
- PDF: 66
- XML: 29
- Total: 363
- BibTeX: 8
- EndNote: 12
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1