Mesh generation and descriptions
This section describes the generation of the different meshes used in this study. In what follows, the following definitions are needed:
Definition 2.1.1 (Spherical polygon). It is a part of the spherical domain bounded by spherical arcs. A spherical arc, in turn, is the spherical equivalent of a straight segment, i.e. the shortest-length curve joining two points on the sphere.
Definition 2.1.2 (Mesh cell). Each polygon within a partition of the spherical domain into spherical polygons, not containing smaller polygons of the same mesh, is referred to as a cell of the mesh. Definition 2.1.3 (Neighbouring cells). Cells adjacent to each other by sharing a spherical arc are called neighbouring cells. The common spherical arc is an edge of both cells.
Quasi-uniform Voronoi-mesh generation
The following paragraph describes the construction of the quasi-uniform hexagonal-pentagonal Voronoi mesh used in this study. This mesh is inspired, in particular, by the mesh of Sadourny et al. , followed by the optimization procedure of Du et al. .
Definition 2.1.4 (Voronoi mesh). A Voronoi mesh is one generated from any distribution of points, each of which is the generator of a cell. The latter is composed of all the points of the domain which are closer to the considered generator than to any other generator.
One result of this definition is that the cells of a Voronoi mesh are polygons. In practice, we start with a regular icosahedron, which is a polyhedron with 20 flat triangular faces. The size of this icosahedron does not matter, as long as the triangular faces are equilateral, and equal to one another. The icosahedron is then projected, radially, onto a sphere. The result of this step is a sphere whose surface is divided into 20 equal spherical triangles, with curved edges and surfaces (Fig. 2.1). The resulting surface is referred to as an icosahedral mesh.
The following step consists of refining this icosahedral mesh. To do so, we follow the method of Sadourny et al. . Each spherical triangle of the spherical icosahedral mesh is considered. Two of its edges are divided into n arcs of equal length. The end of each arc is numbered from 1 to n increasingly from the point common to the two edges, and until the end of the considered edge. Each end of arc is joined with the corresponding end of arc on the second edge of the considered triangle. The lines joining these points together are in turn arcs. Each of the latter is divided into sub-arcs in turn, with an increasing number as we move away from the common triangle vertex. The first arc is kept as 1 segment, the second is divided into 2, and so forth until we divide the nth arc into n sub-arcs. This final arc is, in fact, the third edge of the initial triangle. The ends of these segments are joined together to result in n2 sub-triangles within the initial triangle. This step’s procedure falls into the category of triangular tessellation [Du et al., 1999]. The edges of the resultant spherical sub-triangles, which are no longer equal, have edges which vary to within 18% of each other’s lengths [Wang and Lee, 2011]. This resultant mesh is a Delaunay triangulated mesh.
Variable-resolution Voronoi-mesh generation
Similarly to the quasi-uniform Voronoi meshes, a variable resolution Voronoi mesh is built upon a Delaunay triangulation.
However, to account for the variable resolution, the starting point cannot, in this case, be a regular icosahedron. Instead, a set of points is randomly distributed on the surface of a unit sphere, with a higher probability in regions where higher resolution is sought. A Delaunay triangulation is formed using these points, and then the Voronoi dual is taken in the same way as in the previous case. Also, the target resolution is taken into account during the Lloyd iteration process [Du et al., 1999]. It is important to note that, with variable resolution, the Voronoi mesh is no longer hexagonal-pentagonal: at least some cells have seven edges or more.
Voronoi and Delaunay mesh elements
The following section describes the different parts of the primal and dual meshes, which are defined and annotated for use in the numerical schemes to come. The definitions below apply to the primal and dual meshes simultaneously. The indices used in this section are broadly inspired by the notations of Thuburn et al. .
Each cell is referred to using a fixed cell index, k. These indices are unstructured: there is no particular direction for the index increment (Fig. 2.2). A neighbouring cell may or may not have an index with a unit difference from the current cell index. The total number of indices is equal to the total number of cells, with no two cells sharing the same index. Dual mesh cells are referred to by an index, v.
The cell of index k has an area Ak. A dual cell with index v has an area Av.
Edge values of the mixing ratio
The mixing ratio at the midpoints of the cell edges is denoted by e. This edge value is needed in the numerical schemes detailed in Section 2.3. Its evaluation depends on the chosen numerical scheme, starting from the value at the respective cell center. This step of the numerical formulation of schemes is referred to as reconstruction [Dubey et al., 2015]. To transport the tracer field from one cell to one of its neighbouring cells, we reconstruct the tracer mixing ratio at the edge common to the two cells.
Discretized velocity field
The velocity field, u, on the other hand, is a vector field, and a different approach is used to discretize it. To each edge, of index e, we attach a real number, ue, which we define as the velocity vector—at the edge midpoint—projected onto the normal unit vector of the edge.
This quantity is useful for the second order numerical schemes requiring a discrete gradient (Section 2.3.3). This is the discrete counterpart of the continuous gradient of the continuous mixing ratio field. In our study, we go by the work of Dubey et al. , who demonstrate that a first order approximation of the gradient is sufficient to maintain second order scheme accuracy. Following their first order approximation, a first-order estimate of the gradient, rv is obtained at a each dual cell (with index v) as the dual-cell-averaged gradient: rv = 1 Av Z Vv rdA.
Table of contents :
1.1.1 Sharp plumes
1.1.2 The advection equation
1.2 Methods for the numerical resolution of the advection equation
1.2.1 Cartesian grid formulations on the sphere
1.2.2 Triangular and polygonal meshes
1.2.3 Numerical schemes
1.3 Research questions
1.4 General outline
2.1 Mesh generation and descriptions
2.1.1 Quasi-uniform Voronoi-mesh generation
2.1.2 Variable-resolution Voronoi-mesh generation
2.1.3 Voronoi and Delaunay mesh elements
2.2 Discrete representation of fields
2.2.1 Discrete degrees of freedom
2.2.2 Edge values of the mixing ratio
2.2.3 Discretized velocity field
2.2.4 Discretized gradient
2.3 Numerical solution of the advection equation
2.3.1 Time discretization
2.3.2 Method of lines: time-integration schemes
2.3.3 Method of lines: spatial schemes
2.3.4 Coupled time and space scheme
2.4 Slope limiting
2.5 Cartesian meshes
2.5.1 Cartesian-mesh generation and elements
2.5.2 Operator splitting
2.5.3 One-dimensional transport schemes
3 Two-Dimensional Test Cases
3.1 Tracer-concentration distributions
3.1.1 Uniform distribution
3.1.2 Double cosine-bell distribution
3.1.3 Single cosine-bell distribution
3.2 Wind fields
3.2.1 Solid-body rotation
3.2.2 Tilted solid-body-type periodic oscillation
3.2.3 Dual vortex with solid-body rotation (NL2010)
3.3 Limiting to the Cartesian band
3.4.4 Numerical diffusion
4 Findings on the Numerical Mesh-Scheme Pairs in 2D
4.3 Numerical diffusion
4.3.2 Preservation of non-linear relationships between tracers
4.4.1 Zonal solid-body-rotation convergence
4.4.2 Tilted periodic solid-body convergence
4.4.3 NL2010 low-shear convergence
4.4.4 NL2010 intermediate-shear convergence
4.4.5 NL2010 default-shear convergence
4.5 Partial conclusions
5 Puyehue-Cord´on Caulle 3D Simulation Results
5.1 Puyehue-Cord´on Caulle eruption
5.2 3D simulations
5.3 Trajectory of the plume
5.3.1 SO2 concentration at 330 hPa
5.3.2 Column-integrated SO2 concentration
5.3.3 Plume trajectory
6.1 Main results
6.2 Perspectives for chemistry-transport modeling
7 Sommaire en franc¸ais
7.1 R´ esultats principaux
7.2 Perspectives pour la mod´ elisation chimie-transport