Using an antidiffusive transport scheme in the vertical direction : a promising novelty for chemistry-transport models

The potential interest of the antidiffusive transport scheme proposed by Després and Lagoutière (1999) for resolving vertical transport in chemistry-transport models is investigated in an idealized framework with very encouraging results. We show that, compared to classical higher-order schemes, the Després and Lagoutière (1999) scheme reduces numerical diffusion and improves accuracy in idealized cases that are typical of atmospheric transport of tracers in chemistry-transport models. Increased accuracy and reduced diffusion is spectacular when, and only when vertical resolution is insufficient to properly 5 resolve vertical gradients, which is very frequent in chemistry-transport models. Therefore, we think that this scheme is an extremely promising solution for reducing numerical diffusion in chemistry-transport models.

) and its application to modelling the March 18, 2012 eruption of Mount Etna (Italy). They show that using 25 this transport scheme reduces numerical diffusion and permits a better representation of the volcanic plume after long-range advection, thereby showing that it is possible to strongly reduce numerical diffusion in CTMs without increasing the number of vertical levels. However, that study, set in the framefork of a full-fledge chemistry-transport model fed by real-life atmospheric fields makes it difficult to fully disentangle the effect of the transport scheme itself from other effects such as uncertainties on emission fluxes and mass-wind inconsistencies in the forcing meteorological fields. 30 Therefore, the present study aims at answering the questions on the use of Després and Lagoutière (1999) that could not be adressed in the realistic framework of Lachatre et al. (2020). For that purpose, we have designed three idealized test cases that permit to compare the performance of the Després and Lagoutière (1999) scheme with the classical schemes of Van Leer (1977) and Colella and Woodward (1984), not only in terms of diffusion but also in terms of accuracy compared to the exact solution, which could not be done in Lachatre et al. (2020). We also examine how the performance of Després and Lagoutière 35 (1999) compares to the two above-cited schemes and to the order-1 upwind Godunov scheme when resolution increases.
Section 2 describes the numerical methods that have been used, their implementation, the discretization strategies and the cases that have been designed for the study. Section 3 presents simulation outputs and the diagnostics that have been designed to compare the different transport strategies with each other. Section 4 discusses the results, and section 5 gives our conclusions.
2 Numerical methods and case description 40 Continuity equation for the motion of air is as follows: where u represents wind speed, C is air concentration (all species together) in molecules per unit volume and Φ = Cu represent the air flux vector.
The continuity equation for species s is as follows: where C s is concentration of species s in molecules per unit volume and Φ s = C s u is the flux vector for species s. Equivalently, Eq. 2 becomes: where α s is the mixing ratio for species s: 2 https://doi.org/10.5194/gmd-2020-304 Preprint. Discussion started: 19 November 2020 c Author(s) 2020. CC BY 4.0 License.
Chemistry-transport models try to solve Eq. 3 as accurately as possible on their discretized grids, while keeping the numerical cost of numerical resolution under control.

1d discretization of advection and advection schemes
In 1d, Eq. 2 becomes: Here we follow a semi-lagrangian approach by identifying the air parcels that have entered cell i through its left boundary, and conversely the air parcels that have left cell i through its right boundary (we suppose that the wind is positive). Let ∆ − x be so that the Lagrangian trajectory starting at x i− 1 2 − ∆ − x at time t passes through x i− 1 2 at time t + ∆t. ∆ − x is the distance travelled by the last air particle entering grid cell i at time t + ∆t. Let us define the wind speed representative of facet i − 1 2 60 between t and t + ∆t as u i− 1 2 = ∆ − x ∆t . If we define the wind speed representative of right facet u i+ 1 2 in a similar way as ∆ + x ∆t where ∆ + x is so that the Lagrangian trajectory starting at x i+ 1 2 −∆ + x at time t passes through x i+ 1 2 at time t+∆t, then Eq. 5 can be discretized over cell i as: and Eq. 6 is verified exactly with no particular hypothesis on the wind speed u (x, t) nor the concentration field C s (x, t). For air concentration, the continuity equation can be discretized in the same way: In order to permit monotonicity of the advection scheme in terms of mixing ratios (i.e. the mixing ratio for species s stays within its initial range) , we reformulate Eq. 6 by using fluxes and mixing ratios instead of winds and concentrations by introducing: Then equations 6 and 9 become: and Equations 12-13 are a flux-form reformulation of semi-lagrangian equations 6-9. The form of Eqs. 12-13 makes it straightforward to verify that if Eq. 13 is verified and if α s,i− 1 2 lies between α s,i−1 and α s,i (and if α s,i+ 1 2 lies between α s,i and α s,i+1 ) then the resulting advection scheme guarantees monotonicity of mixing ratios, which is physically desirable and is the reason why chemistry-transport models usually resolve advection of trace species using an approach based on Eq. 12 rather 85 than a straightforward resolution of Eq. 6.
Regarding chemistry-transport models, in practice, approximated values of F i− 1 2 are inferred from the wind and density fields provided by the forcing meteorological dataset. Here we will avoid the typical problem of mass-wind inconsistenciesdiscussed in, e.g., Jöckel et al. (2001); Emery et al. (2011);Lachatre et al. (2020), by working with analytically-defined non-divergent mass fluxes, and constant and uniform air density, so that Eq. 12 is verified exactly by construction.

90
This simplified framework will permit us to focus on the transport scheme of the chemistry-transport model whose task, in flux-form, is to estimate the values of α s,i± 1 2 that are needed for numerical resolution of Eq. 12.
This order-1 scheme is cheap, robust and mass-conservative but extremely diffusive. It is therefore important to find more accurate ways to estimate α s,i+ 1 2 .
The Van Leer (1977) scheme 100 The second-order slope-limited scheme of Van Leer (1977) brought to our notations and assuming uniform air concentration yields the following expression of α s,i+ 1 2 (for F i+ 1 2 > 0).
is the Courant number for the donor cell. If ν > 1, then more mass leaves the cell than the mass that was initially present and the Courant-Friedrichs-Lewy condition is violated, yielding numerical instability. If α s,i is a local 105 extremum of mixing ratio ((α s,k − α s,k−1 ) (α s,k+1 − α s,k ) ≤ 0), no interpolation is performed and α s,i+ 1 2 = α s,i is imposed: in this case, the scheme falls back to the simple Godunov donor-cell formulation (Eq. 14). This order-2 scheme has been used for decades in chemistry-transport modelling, being a good tradeoff between reasonably weak diffusion, at least compared to more simple schemes such as the Godunov donor-cell scheme, and small computational burden compared to higher-order schemes such as the Piecewise Parabolic Method (Colella and Woodward, 1984).

110
The Colella and Woodward (1984) Piecewise Parabolic Method The Colella and Woodward (1984) Piecewise Parabolic Method (PPM) consists in performing a parabolic reconstruction of the concentration field inside each model cell using information from three upwind cells and two downwind cells, and applying limiters to preserve the scheme's monotonicity and stability. The detailed procedure is described in the seminal Colella and Woodward (1984) paper. Our implementation of this method has been adapted from the CASTRO Compressible Astrophysical 115 Solver Almgren et al. (2010). While third-order by design, application of limiters in the vicinity of extrema introduces firstorder truncature errors in their vicinity so that third-order convergence is not expected with the PPM scheme as descibed in Colella and Woodward (1984) (Colella and Sekora, 2008). However, this limitation does not prevent the Colella and Woodward (1984) method to give much better results than simpler order-2 schemes such as Van Leer (1977), so that the Colella and Woodward (1984) PPM scheme has been used successfully for a wide range of applications including meteorological modelling 120 (Carpenter et al., 1990), chemistry-transport modelling (Vuolo et al., 2009), astrophysics (Almgren et al., 2010) etc. 5 https://doi.org/10.5194/gmd-2020-304 Preprint. Discussion started: 19 November 2020 c Author(s) 2020. CC BY 4.0 License.

The Després and Lagoutière (1999) scheme
The scheme of Després and Lagoutière (1999) is defined by their equations 2 to 4. If F i+ 1 2 > 0, these equations brought to the notations of Eq. 12, give: 125 with the same notations as for the Van Leer (1977) scheme (above).
interpolation is performed and the scheme falls back to the simple Godunov donor-cell formulation (Eq. 14). As stated by its authors, this scheme is antidiffusive. Unlike other schemes such as the Van Leer (1977) scheme described above, two unusual choices have been made by the authors in order to minimize diffusion by the advection scheme: -Their scheme is accurate only to the first order

130
-The scheme is linearly unstable, but non-linearly stable (their Theorem 1) The idea of the authors has been to make the interpolated value α s,i+ 1 2 as close as possible to the downstream value (α s,i+1 if the flux is positively oriented). This property is desirable because it is the key property in order to reduce numerical diffusion as much as mathematically possible while still maintaining the scheme stability. The authors present 1d case-studies with their scheme obtaining extremely interesting results: fields that are initially concentrated on one single cell do not occupy more 135 than 3 cells even after a long advection time (their Fig. 2), sharp gradients are very well preserved (their Fig. 1), and, more unexpectedly due to its antidiffusive character, the scheme also performs well in maintaining the shape of concentration fields with an initially smooth concentration gradient. After extensive testing, these authors also suggest (their Conjecture 1) that convergence of the simulated values towards exact values occur even if the time step is reduced before the space step: in simpler terms, this means that the scheme performs very well even at small CFL values, a property that is not shared by most 140 advection schemes. Comparison of Eq. 17 with 16 shows that the numerical cost of the Després and Lagoutière (1999) scheme is about the same as the Van Leer (1977) scheme.

2d discretization of advection and directional splitting
All the simulations in the present study are 2d x-z cases on a domain discretized over a regular cartesian mesh. In order to work with the usual units (concentrations per unit volume and not per unit area), a third, degenerate space dimension has been 145 introduced, with one grid cell in the y direction, δy = δx. Zero mass flux is prescribed in the y direction. Index i is attributed to the x direction, index k to the z direction, and no index is attributed to the degenerate y dimension. Eq. 12 becomes: Here, F i− 1 2 ,k is the time-averaged mass flux of air through the left boundary of cell i, k between t and t + ∆t, and similar definitions for F i+ 1 2 ,k and F i,k± 1 2 . α s,i− 1 2 ,k is the mixing ratio of species s in the air volume entering cell i, k through its left boundary between t and t + ∆t. If V is the geometric volume containing at time t all the air parcels that are going to cross the left boundary of cell i, k between t and t + ∆t then: From a practical point of view, it is extremely difficult to actually find the countours of V, reconstruct and integrate the concentrations of air and of tracer over this volume. This is why, in practice, Eulerian model in cartesian grids tend to split between the two (or three) space directions. Integrating separately direction x and then direction z over time ∆t (generally calles 'Lie splitting") gives order-1 error, and applying the so-called Strang splitting by integrating first in the x direction over ∆t 2 , then in the z direction over ∆t, and finally once again in the x direction over ∆t 2 reduces the splitting error to order-2.

Tested configurations
Since the present study is aimed at studying vertical transport only, we chose to test the Godunov, Van Leer (1977), Colella and Woodward (1984) and Després and Lagoutière (1999) schemes in the vertical direction with the same transport strategy over the x axis, namely the Colella and Woodward (1984), our best available option. Lie splitting has been applied for scheme compbinations with expected order-1 accuracy, while Strang splitting has been applied for scheme combinations with expected 165 order-2 accuracy. The configurations that have been tested are summarized in Table 1.

Test case definition
We have defined three test-cases designed to be representative of long-range tracer transport situations in the atmosphere. The clean air (zero tracer concentration) entering from these boundaries. We will use T = 86400 s, the length of a complete day on Earth, as time scale for the case studies, along with the corresponding pulsation ω = 2π T . The number density of the carrying fluid (air) will be assumed uniform. These simplifying assumptions are designed to ease the formalism and the formulation of exact solutions of the problem. Two situations of relevance for atmospheric tracer transport have been represented, along with another numerical experiment designed in order to investigate the properties of the tested transport configurations in terms of convergence rate. Case 1, presented in 2.2.1 aims to represent the formation of a thin plume from an initially thicker tracer volume through the action of zonal wind shear, a situation typical of long-range advection of polluted plumes in the free troposphere. Case 2, presented 2.2.2 represents long-range advection of a thin plume under the action of a strong zonal wind.
Except for tests of convergence rates for which increasingly fine discretizations have to be tested, the domain is discretized into 80 evenly-spaced cells from west to east (∆x = 25 km) and 24 evenly-spaced cells from bottom to top (∆z = 500 m). In this case, we consider the evolution of an inert tracer initially distributed as follows: α (t = 0, x, z) = 0 otherwise, Zonal wind is constant in time, zonally uniform and vertically sheared: 190 with U 0 = L 2T so that the horizontal motion of fluid at z = H 2 brings it back at its initial position after a time 2T , which will be the duration of the numerical experiment.
We add a vertical wind defined as: w(x, y, z, t) = w 0 cos (ωT ) The vertical wind speed scale is taken as 5 10 −2 m s −1 , a typical scale for synoptic-scale vertical motion in the troposphere.

195
Since ∂u ∂x = ∂w ∂z = 0 and since the density field is uniform, this mass flux is non-divergent. This case can describe a plume that is initially vertical, covering a 50-km wide column (two grid-cells), uniformly spanned over a 3km altitude range (4500 to 7500 m, corresponding to 6 grid cells, see Tab. 2). This initially thick vertical column later evolves under the action of wind shear into a thin layer.
Direct integration of Eq. 24 and then of Eq. 23 give access immediately to the position of a particle initially (t i = 0) at 200 position (x i ; z i ): Eqs. 25-26 describe the superposition of a horizontal motion forced by the horizontal wind speed defined in 23, and an elliptic motion with pulsation ω due to the oscillating vertical speed (Eq. 24) and its interaction with the horizontal wind shear.
(x (t) ; z (t)) being, for any given time t, affine functions of (x i ; z i ), the wind field of Eqs. 23-24 advects straight lines into straight lines. In particular, the initial rectangular zone containing the tracer will be advected, at any given time, into a parallelogram, which summits are readily given by applying Eqs. 25-26 to the summits of the initial rectangle. Inside this moving parallelogram, that is increasingly tilted with time, tracer mixing ratio is equal to 100 ppb, zero outside, giving access 210 to the exact solution of the case at any time.

Case 2: Long-range advection of thin layer
In this case, the initial tracer mixing ratio is as follows:

215
with h 2 = 500 m. These equations describe a zonally infinite layer contained between z = 5500 m and z = 6500 m. As in case 1, α m = 100 ppb.
Zonal wind is constant and uniform: so that the horizontal motion of the fluid brings it back at its initial x coordinate after time 2T , which will be the duration of 220 the experiment. Vertical wind is defined as: w(x, y, z, t) = w 0 cos 4πx L As for Case 1, w 0 = 5 × 10 −2 m s −1 . This vertical wind speed has two maxima and two minima over the horizontal domain.
Integration of Eqs. 29-30 is immediate and give the trajectory of a particle located at (x i ; z i ) at time t = 0: In particular, after a time kT , k being an integer, all the fluid particles are displaced by distance kL 2 in the horizontal direction and back to their initial altitude. At these times, since the initial tracer plume is zonally uniform and infinite, the field of mixing ration will be exactly back to its initial value everywhere.
This case can describe in a simplified way long-range advection of a 1 km thick layer of an inert tracer under the action of 230 a uniform zonal wind and submitted to the action of alternatively ascendant and subsident winds, representing for example in an extremely simplified way the advection of this layer through synoptic-scale structures. In the atmosphere, such fine layers of tracers are frequently formed by stretching of initially thicker polluted layers, as represented in Case 1.

Case 3: Fine layer advection and convergence rate test
This case has been designed to study numerical convergence rate of the various configurations that will be tested as a function 235 of space resolution. The case setup is the same as for Case 2, with wind speeds similar to Eqs. 29-30: 10 https://doi.org/10.5194/gmd-2020-304 Preprint. Discussion started: 19 November 2020 c Author(s) 2020. CC BY 4.0 License.
Due to the need of increasing resolution and therefore the numerical cost of simulations, we simulate only one spatial and temporal period of this case (instead of 2 spatial and temporal periods for Case 2), hence the differences between equations 33-240 34 and 29-30.
The initial tracer mixing ratio is prescribed as: with h 3 = 1500 m. This initial mixing ratio distribution has been designed to be C 2 (and in fact C 3 ) everywhere. It is 245 therefore smooth enaugh to permit a convergence experiment with all the transport schemes which rely, at most, on the existence and continuity of the second-order derivative of the transported field (for the PPM scheme).
For convergence rate tests, 5 different resolutions have been tested (Table 2).

Implementation
Implementation of the idealized experiments and transport scheme have been done within the under development ToyCTM   Table 3. Performance of simulations performed with the Upwind, VL, PPM and DL99 vertical advection schemes relative to the discretized exact solution for Case 1: percent relative error in · 1 and · 2 and percent of total tracer mass contained in the correct envelope diagnostics, by a wide margin: in this case study, the performance gain of DL99 relative to PPM is similar to the gain of PPM 265 relative to Godunov in terms of maximal mixing ratio, percentage of mass in the evvelope and accuracy.
12 https://doi.org/10.5194/gmd-2020-304 Preprint. Discussion started: 19 November 2020 c Author(s) 2020. CC BY 4.0 License. 3.2 Case 2: Long-range advection of thin layer Figure 2 shows the final state of Case 2 simulation for the four advection schemes that have been tested, as compared to the exact solution. Visually, the Després and Lagoutière (1999) scheme has performed best in bringing virtually all the tracer back into its original envelope after two complete vertical oscillations. Its performance in this case is best of all the tested schemes 270 (Table 4), with only a 6% reduction in the maximal value of tracer mixing ratio (49% with the PPM method, even more with the Godunov and Van Leer schemes), and very small error values in · 1 and · 2 compared to the other tested schemes. 93% of the mass is contained in the theoretical envelope where it should be after the end of the numerical experiment (50% only with the PPM scheme).

275
The numerical convergence rates of all the tested advection configurations for · 1 and · 2 based on the last segment in Figs. 3a-b (between n x = 160 and n x = 320) are shown in Tab. 5. In · 1 . These orders of convergence are around 0.8 for  Table 4. Performance of simulations performed with the Upwind, VL, PPM and DL99 vertical advection schemes relative to the exact solution for Case 2: percent relative error in · 1 and · 2 and percent of total tracer mass contained in the correct envelope (a) (b) Figure 3. Convergence rate results for the four tested vertical advection schemes as a function of the number of horizontal points nx in the horizontal direction, in · 1 (panel a) and · 2 (panel b) for Case 3. The black dashed lines shox the slope expected for order-1 (long-dash) and order-2 (short-dash) the DL99 and Upwind schemes because vertical resolution is still insufficient even at this high resolution to ensure theoretical convergence for these schemes. The Van Leer scheme yields a convergence rate of 1.80 in · 1 , 2.43 for the PPM scheme. A smoother test case has been performed to check that all schemes are able to obtain convergence up to their theoretical order 280 (see Appendix). Figure 3 shows that, when resolution is too coarse to appropriately resolve the thin layer, the Després and Lagoutière (1999) scheme strongly overperforms higher-order schemes and offers an accuracy that is substantially better than both the Van Leer (1977) scheme and the Colella and Woodward (1984) scheme. This is consistent with the results obtained on Cases 1 and 2, and explains why in these cases the Després and Lagoutière (1999) scheme permits to obtain excellent results in reproducing 285 the thin layer of tracer. On the other hand, when vertical resolution becomes sufficient to appropriately resolve the smooth tracer layer, higher-order schemes perform better than the Després and Lagoutière (1999) schemes, which in turns falls back towards an accuracy similar to the Godunov scheme, consistent with its expected order or accuracy.  Table 5. Convergence rates in · 1 and · 2

Discussion
The numerical experiments that have been presented confirm the interest of using the Després and Lagoutière (1999) scheme 290 for vertical transport in chemistry-transport models, as already claimed by Lachatre et al. (2020). The idealized framework set up here permits to examine some of the questions that were out of reach in the real case of Lachatre et al. (2020) due to uncertainties on the forcing meteorological fields, the volcanic emissions and to the lack of accurate measurement of plume structure.
First, we show not only that the Després and Lagoutière (1999) scheme is less diffusive than Van Leer (1977) and Colella 295 and Woodward (1984), as could be expected due to its antidiffusive design, but also that, in the presence of sharp vertical gradients that are not adequately resolved at model resolution, using this vertical transport scheme also increases model accuracy with the exact solution show that, at this resolution and for these cases, using Després and Lagoutière (1999) both reduces diffusion and increases accuracy compared to the schemes of Godunov, Van Leer (1977) and Colella and Woodward (1984).
While reduction of diffusion is in line with the results of Lachatre et al. (2020) and could be expected by design of the Després and Lagoutière (1999) schemes, improved accuracy in presence of sharp gradients is a strong argument in favor of using the Després and Lagoutière (1999) scheme for vertical transport in chemistry-transport models.

305
Improved accuracy of a low-order scheme compared to higher-order schemes for a given resolution is not impossible from a theoretical point of view but still counterintuitive since higher-order schemes are expected to reduce numerical error at any given resolution compared to lower-order schemes due to "smarter" reconstruction procedures. To understand this surprisingly good behaviour of Després and Lagoutière (1999) compared to teh higher-order Van Leer (1977) and Colella and Woodward (1984) schemes, we have performed a convergence test for advection of a 3000m-thick layer with a smooth vertical profile. This 310 convergence test (Fig. 3) show that the Després and Lagoutière (1999) performs better than these classical order-2 schemes if model vertical resolution ∆z is equal to 1000 m or 500 m, but that due to their faster convergence rate, order-2 schemes performs better if ∆z ≤ 250 m. In more general words, this result suggest that the Després and Lagoutière (1999) scheme may be expected to perform better than classical schemes in chemistry-transport models for advection of polluted plume thinner 6 ∆z. When vertical resolution becomes much too fine compared to the size of the modelled object (polluted plume thicker that 50 ∆z), accuracy of the simulations with the Després and Lagoutière (1999) scheme stops improving with resolution at order-1 rate (Fig. A1). Examination of the simulation outputs for these configuration reveal that progress in accuracy is hindered by the apparition of undesirable small-scale oscillations that degrade acccuracy (not shown).

320
The Després and Lagoutière (1999) scheme is shown to offer excellent performance in reproducing long-range advection of fine tracer layers when vertical resolution is so coarse that a correct representation of these layers could be expected to be impossible in practice. This improved performance of Després andLagoutière (1999) compared to Van Leer (1977) and Colella and Woodward (1984) translates into reduced numerical diffusion and in improved accuracy compared to these higher-order schemes. Convergence tests show that this improved performance exists only for tracer layers that are represented by less 325 than 6 grid cells in the vertical direction, while for finer resolutions higher-order schemes perform better due to their faster convergence towards the exact solution. In other terms, if model resolution is fine enaugh to represent properly the plume, then higher-order schemes are still a better choice.
We think that these results are important because they explain under which conditions the Després and Lagoutière (1999) is able to reduce excessive vertical diffusion in CTMs, as was observed in Lachatre et al. (2020), show that this reduction 330 in numerical diffusion is not obtained at the expanse of accuracy, which is improved as well. Since the vertical resolution of chemistry-transport model is usually coarse in the free troposphere, while the advected plumes tend to be extremely thin in this part of the atmosphere due to the action of wind shear, the present study along with Lachatre et al. (2020) advocates for using the Després and Lagoutière (1999) transport scheme for chemistry-transport modelling in the free troposphere, and probably even more in the stratosphere where vertical diffusion needs to be extremely reduced. More investigation is needed in real 335 and/or idealized cases to adress several questions: -Does the Després and Lagoutière (1999) perform better in the boundary layer, where CTM resolution is typically finer than in the free troposphere ?
-If not, is it possible to use a traditional transport scheme in the boundary layer and the Després and Lagoutière (1999) scheme in the free troposphere without intruducing numerical artefacts in the buffer zone ?

340
-Atmospheric chemistry being a non-linear process, how does reduction of excessive numerical diffusion in the troposphere affect representation of chemistry inside chemically active air masses such as volcanic or biomass burning plumes.
-Is it desirable, and under which conditions, to use antiffusive transport schemes in the horizontal directions as well ? This case has been designed to study numerical convergence rate of the various configurations that will be tested as a function 350 of space resolution. The case setup is the same as for Cases 2 and 3, with wind speeds from Eqs. 29-30, but the initial tracer mixing ratio is prescribed as: α i (x, z) = α m 4 1 + cos π(z − H/2) h 4 1 + cos π(x − L/2) δ (A1) with h 4 = H 2 = 6000 m and δ = L 2 = 5×10 5 m. This initial tracer distribution represents a cosine bell centered at x = L 2 ; z = H 2 , C ∞ everywhere with extremely smooth variations in space. Domain resolution and sizes are as shown in Tab. 2. Table A1 and Fig. A1 show that all the tested configurations exhibit the expected convergence order at least in · 1 , but accuracy with the DL99 scheme stops improving when vertical resolution becomes "too fine" compared to the thickness of the represented layer. In the present case, this "saturation" of convergence occurs between nx = 80 (nz = 48)and nx = 160 (nz=96):

360
Author contributions. All the authors have contributed to the design of the simulated cases, SM has performed and analyzed the simulations, SM has developed the software with RP, LM, MLand RP have contributed to writing and improving the manuscript.