New numerical approaches for modeling thermochemical convection in a compositionally stratified fluid

Abstract

Seismic imaging of the mantle has revealed large and small scale heterogeneities in the lower mantle; specifically structures known as large low shear velocity provinces (LLSVP) below Africa and the South Pacific. Most interpretations propose that the heterogeneities are compositional in nature, differing in composition from the overlying mantle, an interpretation that would be consistent with chemical geodynamic models. The LLSVP's are thought to be very old, meaning they have persisted throughout much of Earth's history. Numerical modeling of persistent compositional interfaces presents challenges, even to state-of-the-art numerical methodology. For example, some numerical algorithms for advecting the compositional interface cannot maintain a sharp compositional boundary as the fluid migrates and distorts with time dependent fingering due to the numerical diffusion that has been added in order to maintain the upper and lower bounds on the composition variable and the stability of the advection method. In this work we present two new algorithms for maintaining a sharper computational boundary than the advection methods that are currently openly available to the computational mantle convection community; namely, a Discontinuous Galerkin method with a Bound Preserving limiter (DGBP) and a Volume-of-Fluid (VOF) interface tracking algorithm. We compare these two new methods with two approaches commonly used for modeling the advection of two distinct, thermally driven, compositional fields in mantle convection problems; namely, an approach based on a high-order accurate finite element method (FEM) advection algorithm that employs an artificial viscosity technique known as 'Entropy Viscosity' (FEM-EV) to maintain the upper and lower bounds on the composition variable as well as the stability of the advection algorithm and the advection of particles that carry a scalar quantity representing the location of each compositional field. All four of these algorithms are implemented in the open source FEM code ASPECT, which we use to compute the velocity, pressure, and temperature fields associated with the underlying flow field. We compare the performance of these four algorithms by computing the solution to an initially compositionally stratified fluid that is subject to a thermal gradient at a Rayleigh number of Ra = 10 5 with a buoyancy ratio of B = 1.0. For B = 1.0, a value for which the initial stratification of the compositional fields persists indefinitely, our computations demonstrate that the entropy viscosity-based method has far too much numerical diffusion to yield meaningful results. On the other hand, the DGBP method yields good results, although small amounts of each compositional field are numerically entrained within the other compositional field. The particle method yields yet better results than this, but some particles representing the denser fluid are entrained in the upper, less dense fluid and are advected to the top of the computational domain and, similarly, particles representing the less d...

Implementation of a Volume-of-Fluid method in a finite element code with applications to thermochemical convection in a density stratified fluid in the Earth’s mantle

“…In other words, it is unstable [12] and perturbations due to different numerical methods can yield vastly differing results. In fact, in [58] we demonstrated that it suffices to change only the algorithm with which the composition variable C is advected in order to obtain clearly different results at the same output time. For example, see Figure 11 of [58] or compare Figures 5(c)-(d) of [59] to our results here or in [58] or to the results in [73].…”

“…See [27] for a description of this algorithm in ASPECT and a comparison with the advection method with entropy viscosity described above. Also see [58] for a comparison of these two methods with the VOF method described here and a particle method for modeling the solution of the composition equation (22).…”

“…we obtain from each of the six slopes in (55) and (56) in the 3 × 3 block B ij centered on the cell of interest Ω ij and use the line that minimizes the difference between the given volume fractions and the volume fractions due to the lines in (58). We now explain this procedure in a bit more detail.…”

3-D numerical modelling of crustal polydiapirs with volume-of-fluid methods

“…At this point need is to mention the two main classes of VOF methods: methods that try to reconstruct exactly the interface between fluids (e.g. Puckett et al 2018), which requires significant computational time, and methods that do not, such as in the present approaches with JADIM and OpenFOAM. With no interface reconstruction, the thickness of the interfacial region is defined by 0 < C i < 1, and typically occupies two to three grid cells.…”

“…This motivates the use of a volume-of-fluid (VOF) method, dedicated to the conservation of chemical interfaces. Recently, Puckett et al (2018) compared a VOF method (using interface reconstruction) against other standard methods in models of mantle convection. They qualitatively studied two-layer Rayleigh-Taylor and Rayleigh-Bénard systems in 2-D, and concluded that their VOF method may be the most appropriate method for modelling interfaces separating chemical compositions.…”