the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
| 31 Mar 2025
The tracer nudging method for correcting and preventing uneven tracer distributions in geodynamical models
Abstract. Tracers/markers/particles are commonly used in geodynamical models to track composition and sometimes temperature throughout the domain. A common problem is that over time, gaps in the tracer distribution can develop, often resulting in cells with no tracers as well as bunching of tracers. Here a correction method that perturbs or “nudges” the positions of tracers in such a way as to close gaps and eliminate bunching is presented. Test results show that this tracer nudging method is highly effective. Starting from an extremely heterogeneous tracer distribution with large regions of the domain devoid of tracers, it can produce an even distribution in only a few nudge iterations. In a time-stepping situation with a nudge every time-step, the amplitudes of the nudges are small yet sufficient to prevent gaps and bunches, allowing a low-order tracer advection method to be used while maintaining a tracer distribution that is more even than that obtained using higher-order advection methods alone. The nudge essentially corrects any non-conservation error inherent in an advection method. The computational cost is small because the method simply requires solving a Poisson equation.
- Preprint
(1447 KB) - Metadata XML
- BibTeX
- EndNote
Status: final response (author comments only)
-
RC1: 'Comment on egusphere-2025-1354', Anonymous Referee #1, 05 May 2025
This manuscript tackles an important numerical difficulty in geodynamic modelling: the development of uneven tracer distributions during Lagrangian particle advection. The authors propose a tracer-nudging algorithm, derived from the requirement that material density remain constant, which iteratively redistributes tracers until a uniform spatial density is achieved. The idea is elegant and, if widely adopted, could mitigate one of the longest-standing practical problems in high-resolution mantle-convection and lithosphere-deformation studies. While I acknowledge the novelty and potential impact of the work, several aspects of the presentation and validation need to be strengthened before the paper is suitable for publication.
Major comments
- Demonstrate the method on realistic geodynamic problems
The manuscript shows simple circulation tests only. Please include at least one geologically meaningful application—e.g. a high-viscosity-contrast convection benchmark or a 2-D subduction experiment—to illustrate how tracer nudging behaves in complex, time-dependent flow fields and when and how often the nudging is needed.
2. Compare with established schemes
The paper states that the computational cost of tracer nudging is small, yet no comparison is provided. As far as I can see from the manuscript, the additional computational cost is probably higher than the existing remedies that only add a correction items to the velocity interpolation. It would be nice if the author make a comparison with other method and emphasize the advantages (and limitations) of this method3. Clarify applicability to compressible versus incompressible flow
Lines 34–35 imply the method corrects non-divergence-free advection errors. However, many geodynamic models employ compressible Stokes flow, where ∇·v ≠ 0 by definition.-
- State explicitly that the current formulation targets incompressible Stokes problems.
- Discuss whether, and how, the nudging algorithm could be adapted to compressible flow, and what errors might arise if it were used without modification. 4. Streamline the figures
Figures 2 and 3 convey nearly the same information. I suggest keeping the one that best illustrates the tracer-density evolution and moving the other to the supplement or removing it entirely.
Citation: https://doi.org/10.5194/egusphere-2025-1354-RC1 - Demonstrate the method on realistic geodynamic problems
-
RC2: 'Comment on egusphere-2025-1354', Anonymous Referee #2, 19 May 2025
The study by P. Tackley presents a new algorithm for tracer advection for the marker-in-cell method. The algorithm corrects the displacement of tracers by solving for mass conservation of tracers. It is a neat, physics-based approach; however, I find that the choice of the potential function requires more discussion, and, if this algorithm is to be useful to the geodynamics community—or more broadly, to anyone using the marker-in-cell method, the writing can be improved and more tests would be useful. I detail these points below. Considering that the revision work I propose is considerable, I recommend a major revision of the manuscript.1. Choice of the potential function for displacement vector in Eq. 3.It is not clear why the author chose Eq. 4 (instead of Eq. 7) besides convenience to solve a linear Poisson equation. The justification is not there. I agree that it is more convenient to define dx=grad(phi)/rho_t, but what is the physical meaning of this equation? What would be the physical units of such a displacement? Maybe the author has thought about this, but the displacement should not depend on tracer density, and only Eq. 7 is a physical choice. Neither Eq. 4 nor Eq. 9 is balanced in terms of units. Has the author tried solving the non-linear Poisson equation in Eq. 8?L77, Eqs. 7-8: Why is the non-linear equation problematic in areas where rho_t is zero? Not sure this is an issue because in the marker-in-cell method, zero tracer density in control volume is a violation, so this is avoided (which makes the models in Fig. 2 somewhat artificial).Minor comments regarding mathematical theory:L46, Eq. 1: Need definition (equation) for tracer density in the main text. It was obtained from the code but was not obvious. Density = #nlocal_mark * nx*ny*nz / #ntracersEqually, for the mean density, which is taken equal to 1 here.Mean density = nmarkx_cell*nmarkz_cell * nx*ny*nz / #ntracers = 1L47: tracer density is "defined as the mass (or number) of tracers per unit volume". This statement is correct only if it is assumed that all tracers have the same mass, and the density relates to the number of tracers. This needs to be clarified.Eq. 2 to Eq. 3: the discrete version is missing for v=dx/dt.L55-57, Eq. 1: rho_c is mean density of tracers in this study. To generalise, better to use 'initial' instead of 'correct' (i.e., rho_0). Also for the purpose of this study, rho (without the t subscript) can be used for the tracer density.2. Generalisation and stopping criterion for nudgingL111: would the equations change for variable grid spacing? Also, is the scheme generalisable for other boundary conditions, other than impermeable (L72)? Can the boundary conditions for Poisson be generalised from the velocity boundary conditions?The number of 'nudges' required seems arbitrary and it is not clear what is the computational cost of each nudge. I wonder if a stopping criterion can be derived? Example, norm(density error)<=tolerance, where tolerance = 1e-1 or 1e-2 from Fig. 4. Is that an overkill for using the algorithm?3. Geodynamic and performance testsThe tests presented are highly simplified, and maybe artificial (Fig. 2). To demonstrate better the impact of the algorithm, I recommend adding another test showing a geodynamic problem—with sharp velocity interfaces or rotational flows near corners. It would be interesting to see, for example, the problem cited in L27 on eruption/intrusion, which was stated, but not explained, and how this algorithm deals with it.Then, there is the issue of the computational cost of this algorithm (with nudges). Simple Matlab scripts run fast, but the marker-in-cell method becomes expensive for 2-D and 3-D problems. What is the cost of this algorithm? How does it vary with grid size (resolution)?L201-202: says the computational cost of solving Poisson is small, but need to demonstrate. For example, using a geodynamics test case, show time to solve for Stokes and time for tracer correction (Poisson*number_nudges).4. Focus and style of writingThe algorithm corrects the displacement of markers by solving for conservation of tracer mass. This does not transpire from the nudging and bunching described. The writing is too focused on the algorithm (which comes out almost as a trick to deal with marker dispersion/clustering), and not so much on the physics or why is it better than other methods. I think it is not doing a great job advertising why this could be a useful technique for geodynamics, and the informal style of writing does not help (e.g. usage of regular+, Euler+nudge+random without clear definition).Other examples,L40: the goal shouldn't be to nudge the tracers, but to correct their displacement in order to keep uniform density.L46: the goal is to preserve the initial tracer distribution.The Abstract should also be revised because I had the following questions:L9: is the correction method physics-based? A sentence describing the method is lacking. Only results are described.L14: what does the author refer to as non-conservation errors? How do they occur? They are introduced somewhat late, even in the abstract.The details on the content of the Matlab routines in Section 3 (L112-118) and Fig 1. should be in an appendix. Instead, the author can outline a pseudo-algorithm such as
- Initialise tracer distribution, rho = rho_0
- Time loop:
- advect tracers with velocity field v,- calculate rho and density error,- correct tracer location by solving for displacement dx (Poisson equation).Minor commentsL6-7, 17: tracers may track other fields besides composition and temperature throughout the domain. Their advantage is that they may perform better at advecting these fields compared to grid-based methods.L7: why does the problem occur? i.e., rotating fields in a box domain or at sharp interfaces.L8: Suggest to use 'clustering' or 'accumulation' instead of 'bunching'.L20: (e.g., \cite not \citep)L21: why all the major codes use marker-in-cell? It is better to justify advantages rather than state a common practice.L24: (e.g. 5-50) not necessary. The accuracy of the marker in cell technique relies on having a high density of tracers in each cell.L29: why the errors occur in the first place?L30: physics-based remedy is requiredL31: one previous solution is to create new tracers, but why is this not a good remedy?L44, L96: just because it is already implemented in StagYY, it is not a statement of validation or verification. StagYY should cite the work in the current manuscript to demonstrate the method, not the other way around.L109: first time density bounds are given; this should occur earlier when density is defined.L149: unclear statement. Tracers at the boundaries were nudged only a fraction of the calculated displacement to avoid crossing the boundaries.Eqs. 11-12: could plot either the stream function or plot the velocity field to understand why the center particles might get dispersed.Section 4.2 notation: need partial differentials (\partial).Eq. 12: mistake Vx = 1/Lz sin() cos()L172: define Euler+nudge+randomL174: calculated as in Eq. 12 (not above)Figures 2 and 3 do not have colorbars or state that blue points represent tracers.Fig 5: what is the computational cost or runtime between the cases? Add colorbars.Fig 6: all cases should start with the same initial condition, and then additionally a case with completely random distribution.Matlab: Add in readme how to run the script. I.e. change ntest = 1; choice between 1-11.Citation: https://doi.org/10.5194/egusphere-2025-1354-RC2 - AC1: 'Comment on egusphere-2025-1354', Paul Tackley, 23 Jun 2025
Model code and software
Software accompanying manuscript "The tracer nudging method for correcting and preventing uneven tracer distributions in geodynamical models" ETH Zurich and Paul J. Tackley https://doi.org/10.5281/zenodo.15065273
Viewed
| HTML | XML | Total | BibTeX | EndNote | |
|---|---|---|---|---|---|
| 608 | 35 | 16 | 659 | 15 | 25 |
- HTML: 608
- PDF: 35
- XML: 16
- Total: 659
- BibTeX: 15
- EndNote: 25
Viewed (geographical distribution)
| Country | # | Views | % |
|---|
| Total: | 0 |
| HTML: | 0 |
| PDF: | 0 |
| XML: | 0 |
- 1