Research into the use of unstructured mesh methods in oceanography has been growing steadily over the past decade. The advantages of this approach for domain representation and non-uniform resolution are clear. However, a number of issues remain, in particular those related to the computational cost of models produced using unstructured mesh methods compared with their structured mesh counterparts. Mesh adaptivity represents an important means to improve the competitiveness of unstructured mesh models, where high resolution is only used when and where necessary. In this paper, an optimization-based approach to mesh adaptivity is described where emphasis is placed on capturing anisotropic solution characteristics. Comparisons are made between the results obtained with uniform isotropic resolution, isotropic adaptive resolution and fully anisotropic adaptive resolution.
Abstract The determination of the temperature in and above the slab in subduction zones, using models where the top of the slab is precisely known, is important to test hypotheses regarding the causes of arc volcanism and intermediate-depth seismicity. While 2D and 3D models can predict the thermal structure with high precision for fixed slab geometries, a number of regions are characterized by relatively large geometrical changes over time. Examples include the flat slab segments in South America that evolved from more steeply dipping geometries to the present day flat slab geometry. We devise, implement, and test a numerical approach to model the thermal evolution of a subduction zone with prescribed changes in slab geometry over time. Our numerical model approximates the subduction zone geometry by employing time dependent deformation of a Bézier spline that is used as the slab interface in a finite element discretization of the Stokes and heat equations. We implement the numerical model using the FEniCS open source finite element suite and describe the means by which we compute approximations of the subduction zone velocity, temperature, and pressure fields. We compute and compare the 3D time evolving numerical model with its 2D analogy at cross-sections for slabs that evolve to the present-day structure of a flat segment of the subducting Nazca plate.
Abstract The thermal structure of subduction zones is fundamental to our understanding of the physical and chemical processes that occur at active convergent plate margins. These include magma generation and related arc volcanism, shallow and deep seismicity, and metamorphic reactions that can release fluids. Computational models can predict the thermal structure to great numerical precision when models are fully described but this does not guarantee accuracy or applicability. In a trio of companion papers, the construction of thermal subduction zone models, their use in subduction zone studies, and their link to geophysical and geochemical observations are explored. In this part II, the finite element techniques that can be used to predict thermal structure are discussed in an introductory fashion along with their verification and validation. Steady-state thermal structure for the updated subduction zone benchmark. a) Temperature predicted by TF for case 1; b) temperature difference between TF and Sepran using the penalty function (PF) method for case 1 at f m =1 where f m represents the smallest element sizes in the finite element grids near the coupling point; c) slab top temperature comparison for case 1; and d)–f) as a)–c) but now for case 2. The star indicates the position or temperature conditions at the coupling point.
Small-scale experiments of volcanic ash particle settling in water have demonstrated that ash particles can either settle slowly and individually, or rapidly and collectively as a gravitationally unstable ash-laden plume. This has important implications for the emplacement of tephra deposits on the seabed. Numerical modelling has the potential to extend the results of laboratory experiments to larger scales and explore the conditions under which plumes may form and persist, but many existing models are computationally restricted by the fixed mesh approaches that they employ. In contrast, this paper presents a new multiphase flow model that uses an adaptive unstructured mesh approach. As a simulation progresses, the mesh is optimized to focus numerical resolution in areas important to the dynamics and decrease it where it is not needed, thereby potentially reducing computational requirements. Model verification is performed using the method of manufactured solutions, which shows the correct solution convergence rates. Model validation and application considers 2-D simulations of plume formation in a water tank which replicate published laboratory experiments. The numerically predicted settling velocities for both individual particles and plumes, as well as instability behaviour, agree well with experimental data and observations. Plume settling is clearly hindered by the presence of a salinity gradient, and its influence must therefore be taken into account when considering particles in bodies of saline water. Furthermore, individual particles settle in the laminar flow regime while plume settling is shown (by plume Reynolds numbers greater than unity) to be in the turbulent flow regime, which has a significant impact on entrainment and settling rates. Mesh adaptivity maintains solution accuracy while providing a substantial reduction in computational requirements when compared to the same simulation performed using a fixed mesh, highlighting the benefits of an adaptive unstructured mesh approach.
Abstract The thermal structure of subduction zones is fundamental to our understanding of physical and chemical processes that occur at active convergent plate margins. These include magma generation and related arc volcanism, shallow and deep seismicity, and metamorphic reactions that can release fluids. Computational models can predict the thermal structure to great numerical precision when models are fully described but this does not guarantee accuracy or applicability. In a trio of companion papers, the construction of thermal subduction zone models, their use in subduction zone studies, and their link to geophysical and geochemical observations are explored. In part I, the motivation to understand the thermal structure is presented based on experimental and observational studies. This is followed by a description of a selection of thermal models for the Japanese subduction zones.
Abstract Mineral grain size in the mantle affects fluid migration by controlling mantle permeability; the smaller the grain size, the less permeable the mantle is. Mantle shear viscosity also affects fluid migration by controlling compaction pressure; high mantle shear viscosity can act as a barrier to fluid flow. Here we investigate for the first time their combined effects on fluid migration in the mantle wedge of subduction zones over ranges of subduction parameters and patterns of fluid influx using a 2‐D numerical fluid migration model. Our results show that fluids introduced into the mantle wedge beneath the forearc are first dragged downdip by the mantle flow due to small grain size (<1 mm) and high mantle shear viscosity that develop along the base of the mantle wedge. Increasing grain size with depth allows upward fluid migration out of the high shear viscosity layer at subarc depths. Fluids introduced into the mantle wedge at postarc depths migrate upward due to relatively large grain size in the deep mantle wedge, forming secondary fluid pathways behind the arc. Fluids that reach the shallow part of the mantle wedge spread trench‐ward due to the combined effect of high mantle shear viscosity and advection by the inflowing mantle and eventually pond at 55–65 km depths. These results show that grain size and mantle shear viscosity together play an important role in focusing fluids beneath the arc.
Abstract Numerical simulations of thermal convection in the Earth's mantle often employ a pseudoplastic rheology in order to mimic the plate‐like behavior of the lithosphere. Yet the benchmark tests available in the literature are largely based on simple linear rheologies in which the viscosity is either assumed to be constant or weakly dependent on temperature. Here we present a suite of simple tests based on nonlinear rheologies featuring temperature, pressure, and strain rate‐dependent viscosity. Eleven different codes based on the finite volume, finite element, or spectral methods have been used to run five benchmark cases leading to stagnant lid, mobile lid, and periodic convection in a 2‐D square box. For two of these cases, we also show resolution tests from all contributing codes. In addition, we present a bifurcation analysis, describing the transition from a mobile lid regime to a periodic regime, and from a periodic regime to a stagnant lid regime, as a function of the yield stress. At a resolution of around 100 cells or elements in both vertical and horizontal directions, all codes reproduce the required diagnostic quantities with a discrepancy of at most in the presence of both linear and nonlinear rheologies. Furthermore, they consistently predict the critical value of the yield stress at which the transition between different regimes occurs. As the most recent mantle convection codes can handle a number of different geometries within a single solution framework, this benchmark will also prove useful when validating viscoplastic thermal convection simulations in such geometries.