The spectral element method on variable resolution grids : evaluating grid sensitivity and resolution-aware numerical viscosity

Introduction Conclusions References


Introduction
In climate and weather forecast applications, there is an increased demand for variable-resolution capabilities.This demand is motivated by the need to resolve various temporal and spatial scales in forecast and regional climate studies with limited computational resources.Several approaches can be employed to this end, including nesting techniques, multi-physics modeling, and multi-resolution simulations, recently overviewed in Ringler et al. (2011).
Here, we focus on the multi-resolution approach made possible by the spectral element method (SEM).We use global spherical grids, constructed to have uniform high resolution over a region of interest, uniform low resolution over the rest of the globe, and a transition region between them.The SEM has a long history of running on highly unstructured grids, including spherical geometry with climate applications (Fournier et al., 2004;St.-Cyr et al., 2008;Baer et al., 2006;Marras et al., 2014).Here, our goal is to evaluate and improve the multi-resolution capabilities of the SEM formulation from the High-Order Method Modeling Environment (HOMME), recently adopted as the default dynamical core by the Community Atmosphere Model (CAM) (Dennis et al., 2012).HOMME uses a locally conservative/mimetic formulation from Taylor and Fournier (2010) and relies on a constant-coefficient hyperviscosity term to both dissipate energy near the grid scale and to damp grid scale modes with spurious propagation (Ainsworth and Wajid, 2009).This hyperviscosity operator is not suitable for variable-resolution grids, and thus we consider two resolution-aware extensions.The first is the straightforward extension of allowing the hyperviscosity coefficient ν to depend on the local element length scale.This approach was used in Zarzycki et al. (2014a, b) for realistic CAM simulations.Here, we use the shallow-water equations to show some deficiencies with this approach and then that better results are obtained with a tensor-based hyperviscosity operator, which can better represent both length scales within non-square spectral elements.We evaluate this approach using the shallow-water equations on the sphere with the two-dimensional version of HOMME's spectral element dynamical core.There have been other modifications of the viscosity and hyperviscosity operators.For example, Yu et al. (2014) uses a flowdependent coefficient for the Laplacian.In Dobrev et al. (2012), several tensor coefficient viscosity operators were considered, including a formulation in which directional viscosity coefficients were chosen to depend on directional length scales of a deformed Lagrangian volume.Here, we follow a similar approach, but the length scales come from the deformation of the elements, especially those in the grid transition region.
The mimetic SEM requires conforming quadrilateral grids.To generate variable-resolution grids, we employ a well-known grid-generating toolkit, CUBIT (https://cubit.sandia.gov).In the grid-transition region, the CUBIT grids have many valence 6 nodes (corner nodes shared by six elements).For grids in spherical geometry, it is possible to construct transition regions with mostly valence 5 nodes, and such elements will have less acute angles.To generate grids with low-valence nodes, we use the Spherical Quadrilateral Grid Generator (SQuadGen, http://climate.ucdavis.edu/squadgen.php).This toolkit uses a paving technique (Blacker and Stephenson, 1991) in combination with a set of lowvalence tiles to generate smooth quadrilateral grids based on cubed-sphere geometry.Regions of enhancement are determined via a user-specified image file, which is mapped onto a cubed-sphere grid.Grid smoothing is performed via straightforward application of spring dynamics in threedimensional geometry (Persson and Strang, 2004).Grids obtained via this technique exhibit several improved characteristics, including greater uniformity in the transition region, and elements with angles that are closer to 90 • .
In this study, we use multi-resolution grids with a single region of quasi-uniform high resolution ( x high ), which transitions to a quasi-uniform grid of low resolution ( x low ), covering most of the globe.
While evaluating the model's performance, it is natural to compare the multi-resolution simulation with the corresponding x low and x high uniform grid simulations.Motivated by climate applications, and following Weller et al. (2009); Ringler et al. (2011), we evaluate x low / x high variable-resolution simulations with two criteria: 1. Refinement does no harm to the global scales.For the shallow-water equation initial value problem, at short times, it is reasonable to expect many features in the refined region would not be contaminated by information from the low-resolution region.In contrast, for longer times, we expect all features to be influenced by information from both the low-and high-resolution regions.So, in general, global errors should depend on x low and may not be improved by the presence of a x high resolution region.However, they should not be worsened by the presence of refinement and we thus expect global errors in a multi-resolution simulation to be as low as those obtained by a uniform x low simulation.
2. Local scales are resolved in the refined region.The purpose of the multi-resolution approach in climate modeling is not to reduce the initial value problem error, but to resolve features of interest, such as hurricanes, eddies, or topographically driven features in select regions at a lower cost.We thus expect that, in the refined region, the multi-resolution simulation can resolve the same scales as the uniform x high resolution simulation.With hyperviscosity, this requires a resolutionaware formulation, which locally matches what would be used in a uniform resolution simulation.
In the study, we use the popular set of shallow-water test cases on the sphere compiled by Williamson et al. (1992) to show that the SEM satisfies both requirements when tensor hyperviscosity is used.We show that tensor hyperviscosity is both more accurate and more robust with respect to grid quality.The rest of the paper is organized as follows: in Sect.2, we introduce two dissipation mechanisms -scalar and tensor hyperviscosity; in Sect.3, we discuss grid refinement techniques; in Sect.4, we describe shallow-water test cases; in Sect.5, we present numerical results.

Hyperviscosity formulations
In many climate models, a hyperviscosity term is added to the right hand side of the dynamical equations for both physical and numerical reasons.Hyperviscosity is preferred over regular viscosity because it is more scale selective.It strongly damps grid scale modes while having less of an impact on the resolved large-scale modes.To introduce the various types of hyperviscosity that will be considered here, we first consider a model equation containing only the hyperviscosity operator: acting on scalar fields with a constant coefficient ν.Here, subscript t denotes differentiation with respect to time.The tensor formulation replaces ν 2 by (∇ •τ ∇) for a symmetric positive definite matrix τ .That is, instead of Eq. (1), we consider Our intention is to derive a formulation for τ that is suitable for uniform and quasi-uniform grids and that can be extended to non-uniform grids.
We start by writing equation Eq. (2) in weak form.Since Eq. ( 2) is equivalent to the set of equations Q t = −(∇ •τ ∇)q, q = Q, we rewrite the set as a system of integral equations: Here, is the problem domain, which, in our case, will be the sphere of the radius R.This system of equations is discretized by the standard SEM and solved for all SEM test functions φ 1 and φ 2 .We first decompose the domain into a set of quadrilateral elements on the surface of the sphere ( m , m = 1, . .., M, as in Fig. 1), and then write Each term in this sum is then written as an integral over the reference element ref = [−1, 1] × [−1, 1].We define r(x; m) as the map from ref to m , with r ∈ m as a point on the sphere and x = (x 1 , x 2 ) ∈ ref .We require this map to be differentiable and invertible, and further define where D is a 2 × 2 matrix whose columns are the covariant basis vectors expressed in spherical coordinates.The map and analytic expressions for D are given in the Appendix.
The integral over each spherical element m can then be written with respect to ref , using derivatives with respect to the reference element coordinates: where det(D) dx 1 dx 2 is the transformed area measure, and τ is the tensor expressed in spherical coordinates.Note that the discrete operator will conserve mass ( Q), since the change in mass is obtained by taking test function φ 1 = 1, and then the right hand side of Eq. ( 5) will be 0 if we ensure ∇1 = 0.
We now consider the eigenvalues of the metric tensor and its inverse: with the orthonormal matrix E whose columns are the basis vectors, which diagonalize the Laplace operator.For any practical grid, both D and D −1 are well defined, and hence the symmetric metric tensor is guaranteed to have such an eigenvalue decomposition with positive eigenvalues.We note that these two eigenvalues can be used to define the two length scales associated with each element m .For the special case of a grid of rectangular elements in Cartesian geometry of size l x × l y , we have E as the identity matrix and For a general, possibly distorted element, we define its two length scales by 2 √ λ 1 and 2

Constant-coefficient hyperviscosity
The traditional constant-coefficient hyperviscosity is obtained by taking τ = νI, with the identify matrix I.For uniform resolutions with an average grid spacing of x, often ν = c 0 ( x) 3.2 , for some constant c 0 > 0. This scaling is obtained by experimentation and is found to be effective for several different dynamical cores over a wide range of resolutions (Boville, 1991;Takahashi et al., 2006 , 2012).We take a slightly more general form and allow ν = c 0 ( x) s for a scaling parameter s.The constantcoefficient hyperviscosity is used for quasi-uniform grids, where we follow the convention of defining x by the average number of degrees of freedom on the equator.For square elements, l x = p x, where p is the polynomial order of the basis functions in the SEM.
In order to motivate how we generalize this operator to a full tensor, we first express τ = νI in the basis, which diagonalizes the Laplace operator (the local element basis defined by E).Some algebra shows that Below, we ensure that the more general tensor formulations retain this scaling in regions in which the grid has a uniform resolution of λ 1 λ 2 .

Scalar hyperviscosity
For scalar hyperviscosity, we once again take τ = νI and now allow ν to vary in space.The natural choice is to use the same scaling as with the constant-coefficient operator, ν = c 0 ( x) s , but with x now chosen locally for each element.To preserve the scaling for the constant-coefficient operator, but also to ensure that the coefficient does not become too small (and thus provide insufficient dissipation), we use Eq. ( 10) and approximate the resolution locally by taking l x = 2(max{λ 1 , λ 2 }) 1/2 and x = l x /p.For scalar hyperviscosity, this tensor scales according to On a Cartesian grid with square elements, λ 1 = λ 2 = (l x /2) 2 , and so the scalar and constant-coefficient operators are identical.For a quasi-uniform grid, where λ 1 λ 2 (l x /2) 2 , the scalar and constant-coefficient operators will have the same scaling with resolution.

Tensor hyperviscosity
Now consider a grid with only rectangles of size l x l y .Based on our expected scaling hyperviscosity with resolution, we note that the scalar hyperviscosity above would give us the desired amount of dissipation in the x direction, but would have excessive dissipation acting in the y direction.The natural choice for a grid of rectangles is a tensor coefficient: For a grid of pure rectangles, D is diagonal and E = I, so τ expressed in the E basis is given by We use this same formulation for unstructured grids by defining the two locally varying element length scales as in Eq. ( 10) and taking On a Cartesian grid of pure rectangles, the scalar and tensor operators are identical.For a quasi-uniform grid of rectangles, where λ 1 (l x /2) 2 and λ 2 (l y /2) 2 , the scalar and tensor operators will have the same scaling with resolution.
For direct comparison, we summarize the three different approaches: -Constant-coefficient: for quasi-uniform grids with average grid spacing x; -Scalar: ν = ν(r) depends on local element length scales; -Tensor: τ depends on local element length scales; For a smoothly deformed grid, the matrix entries of τ will be smooth functions over the domain m .For our discrete grids, we ensure τ is continuous across element edges by applying the standard SEM projecting operation to each entry of τ .We further reduce variations in τ by computing it only at element corner points, forming a bilinear fit to these corner values, and using this bilinear approximation at all nodes within the element.

Hyperviscosity acting on vector fields
In our method, we represent vector fields u in spherical coordinates, but care must be taken not to differentiate individual vector components when represented in spherical coordinates, since they are multiply valued at the poles.Instead, we transform u from spherical to Cartesian coordinates and solve Eqs.(5-6) for each Cartesian component of the velocity field, and then transform the result back to spherical coordinates.
Figure 2. Schematic idea of constructing refined grids for the conformation of quadrilaterals on a sphere: we start from uniform grid with x low .Next, the region of desired refinement is replaced with uniform elements of size x high .The grey area approximately defines a transition region, which is constructed by substituting quadrilaterals with x low by transition templates.After the transition region is assembled, spring dynamics can be used to smooth the grid.

High-and low-connectivity conforming quadrilateral grids on the sphere
The mimetic formulation of the SEM that we are using requires the conformation of quadrilateral grids.The cubed sphere is a popular way to construct these grids on the sphere with quasi-uniform resolution.An inscribed cube is projected onto the surface of the sphere and each panel is further subdivided into a grid of elements, as shown in Fig. 1.
For multi-resolution, we consider grids with a single refined region over an area of interest.We define a coarse resolution x low and fine resolution x high .We restrict ourselves to choices x high = x low /N , N = 2, 4 and 8. Starting from a cubed-sphere grid with a resolution of x low , the region under refinement is substituted by uniform elements with x high , as shown in Fig. 2. The approximate placement of the transition region is colored grey.For each N , we generate a family of grids with different low-resolution regions ( x low ).Following Ringler et al. (2011), we refer to these family of grids as ×1, ×2, ×4, and ×8.The ×1 family is the set of uniform cubed-sphere grids, with resolutions ranging from 3 to 0.5 • .The ×2 family (shown in Fig. 3) is similar, but each ×2 grid has a refined region with twice the resolution (N = 2).The ×4 family has a refined region with four times the resolution (N = 4) and the ×8 family has N = 8.
It is nontrivial to construct the transition region.We need to avoid hanging nodes and prefer the elements to be as close to squares as possible.In Fig. 4, we provide two example ×8 grids with the same x low = 3 • .Here, we represent two approaches to construct the transition region.Both are based on periodic templates, as seen in Fig. 5.The transition region in Fig. 5a is constructed by CUBIT, a grid-generating software for complex geometries in two and three dimensions (https://cubit.sandia.gov).Figure 5b contains the transition generated by SQuadGen.SQuadGen was developed to generate two-dimensional refined spherical grids based on a cubed sphere (http://climate.ucdavis.edu/squadgen.php).
As seen in Fig. 5, the transition region in Fig. 5a contains nodes of higher valence comparing to the similar region in Fig. 5b.In this context, valence of a node is a number of edges to which it is connected.Node valence greater than 4 results in quadrilaterals with more acute angles and more distorted elements, and thus lower valence grids are usually preferred.In Fig. 5a, most nodes are of valence 3-6, with a few of valence 7.In Fig. 5b, most nodes are of valence 3-5 with a few nodes of valence 6.For the approach used in Fig. 5a, it is possible to avoid valence-7 nodes altogether with a less automated, more user-dependent procedure, but it is not possible with valence-6 nodes.
After the transition has been constructed, it is standard procedure to apply a smoothing algorithm to further improve grid quality.SQuadGen employs the algorithm of Pers-O.Guba et al.: The spectral element method (SEM) on variable-resolution grids (a) Grid with a narrow, more distorted transition region from CUBIT mesh generation software with grid smoothing turned off (b) Grid with a wider, more uniform transition region from SQuadGen and grid smoothing with spring dynamics Moreover, the user can choose a halo size around the inner and outer boundaries of the transition region, in terms of graph distances.This halo is then used to define a region in which smoothing will be applied.To investigate the performance of our resolution-aware hyperviscosity operators, we take two extremes: non-smooth grids with higher valence nodes generated by CUBIT, and smoothed grids with lower valence nodes generated by SQuadGen.We call the former high-connectivity (or highly distorted) grids and latter lowconnectivity grids.

Shallow-water test cases
The shallow-water equations on a rotating sphere are given by Here, h is the fluid thickness, u represents velocity, ζ = k • ∇ × u the vorticity, f is the Coriolis parameter, g is gravity, and b denotes bottom topography.The equations are discretized following Taylor and Fournier (2010) with the hyperviscosity operator discretized, as per Sect.2. We take p = 3 for a fourth-order accurate spatial discretization and use the second-order accurate leapfrog-trapezoidal timestepping method.
For our studies, we choose two standard shallow-water test cases from Williamson et al. (1992): test case 2 (TC2) and test case 5 (TC5).TC2 represents a global steady state zonal geostrophic flow.Since the analytical solution is known, TC2 is often used to investigate convergence rates.The solution is very smooth, resulting in small errors, but the errors are very sensitive to local fluctuations in truncation error, such as those caused by grid irregularities.Following Ringler et al. (2011), we run TC2 for 12 simulation days, instead of the originally proposed 5 days, in order to allow longer time for the error growth to disrupt the steady state solution.
TC5 consists of a more realistic zonal flow over an isolated mountain run for 15 days.Error measures are obtained from a high-resolution reference simulation.Establishing convergence rates is difficult as the rate decreases to 0 as the errors approach the uncertainty in the reference solution.Instead, global errors are used primarily to measure the impact of the refined region.TC5 has much larger errors than those in TC2, which are less sensitive to small fluctuations in the local truncation error.In TC5, we examine the vorticity, which contains small-scale structures that are only captured at high resolution.We use the vorticity results to ensure that these structures can also be captured within the high-resolution region of a variable-resolution grid.
One of the purposes of this study is to confirm that, in the SEM with hyperviscosity, the large-scale errors are not harmed by the presence of refinement and thus are primarily controlled by coarse resolution ( x low ).For this, we collect a series of grids with ×2, ×4, and ×8 refinements.For both TC2 and TC5, the refined region covers a circle with coordinates λ = 3π/2, θ = π/6 and radius π/9.This placement is centered over the TC5 mountain.We summarize charac-  teristics of the grids in Table 1.In Table 2, we summarize some parameters for simulations in Figs.6-10.Resolutions x low and x high are computed according to the formula for an equatorial uniform resolution, considering that the whole sphere is covered by corresponding large or small elements.

Grid and hyperviscosity sensitivity in TC2
We present error plots for the TC2 height field h after 12 days in Fig. 6.The error is contoured for several different meshes all with a relatively low x low = 3 • resolution.Part a contains a plot for a uniform resolution and constant-coefficient hyperviscosity.Plots b and d are simulations with the scalar hyperviscosity on ×8 grids.Plots c and e are simulations with the tensor hyperviscosity on ×8 grids.We first note that the errors are quite small relative to the height field (which ranges from 1000 to 3000 m).The height field is not plotted since it would be identical to the analytic solution in Williamson et al. (1992).For the uniform grid with constant-coefficient hyperviscosity, the error is quite uniform with no indication of any grid sensitivity.There is no visible m = 4 mode that might be expected because of the cubed-sphere grid.
To investigate performance and robustness of two dissipation mechanisms, we chose the two ×8 grids shown in Fig. 4: a highly distorted unsmoothed grid generated by CUBIT, deg3×8highconn, and a grid with low-connectivity nodes and selectively applied smoothing generated by SQuadGen, deg3×8lowconn.Both the ×8 grids have the same lowresolution region as the uniform grid, and thus the errors in the ×8 grids should be equal or lower than the errors in panel a.This is obviously not the case when scalar hyperviscosity is used, as seen in panels b and d.Those results are contaminated with significant numerical noise with larger errors than panel a.This is even true in the Southern Hemisphere, away from the region of local refinement.Comparing the distorted grid (panel b) with the low-connectivity grid (panel d), we see that the scalar hyperviscosity has some grid sensitivity, as the better quality grid in panel d shows more In addition, the quality of the underlying grid significantly improves the outcome around the refined region when using scalar hyperviscosity, as seen by comparing (b) and (d).Contrary to scalar hyperviscosity, tensor hyperviscosity is more robust with respect to mesh quality, as follows from comparing (c) and (e).Simulations (c) and (e) with tensor hyperviscosity also exhibit a substantial error reduction in the vicinity of the refinement compared to the coarse uniform resolution in (a).
zonal contours similar to panel a and somewhat less noise in the refined mesh region, although panel d does have larger minimum and maximum errors (given in the figure caption) than panel b.Both panels b and d have minimum and maximum errors significantly larger than those in panel a.Contrary to this, results using the tensor coefficient hyperviscosity are very close to panel a and much less sensitive to the different types of refinement.The minimum and maximum errors with either the distorted grid (panel c) or the lowconnectivity grid (panel e) are slightly less than the values obtained on the uniform grid (panel a).In both panels c and e, the error contours in the Southern Hemisphere are almost identical to panel a.In the Northern Hemisphere, the errors are sensitive to the presence of refinement, but are actually lower in this region than they are with panel a.Thus, with the tensor hyperviscosity, the presence of mesh refinement does no harm to the solution and actually results in a minor local reduction in the error.

Vorticity in TC5
We now examine the vorticity field for TC5.The vorticity after 15 days is plotted in Fig. 7 high-resolution grid.Note the sharp gradient in the flow that is well resolved in our reference solution but not present in the low-resolution result.
We show results computed using locally refined grids in Figs. 9 and 10.Based on the results presented for TC2, here we compare the worst and best performing extremes: scalar hyperviscosity running on the highly distorted ×8 mesh in Fig. 4a, and tensor hyperviscosity on the low-connectivity ×8 mesh in Fig. 4b.
The vorticity is plotted over the mesh refinement region for these two simulations in Fig. 9.The computation with scalar hyperviscosity develops unphysical oscillations, which are not present in the tensor hyperviscosity result or reference solution.Figure 9b shows a very smooth field across the highly non-uniform transition region, and the sharp gradient that is present in Fig. 8b is resolved without numerical noise.Note that the exact matching of Figs.8b and 9b is not expected because the reference solution used a grid with 3 times finer resolution than the finest resolution used in the ×8 grids.
To quantify these observations, we plot the error in the vorticity field in Fig. 10.We show the error for the two ×8 simulations, as well as a uniform low-resolution simulation.The error is computed using our high-resolution reference solution as an approximation to the exact solution.The noise seen in the scalar hyperviscosity vorticity field (Fig. 9) is more evident throughout the refined region in the vorticity error plot (Fig. 10b).With the tensor viscosity on the low-connectivity grid, Fig. 10c shows very little noise in the refinement region and the mesh transition region.In addition, the error is substantially reduced in the refinement region as compared to the low-resolution uniform grid solution (Fig. 10b).The fact that the error can be reduced by local mesh refinement after 15 days in this region suggests that the solution contains standing features induced by the mountain, which benefit from mesh refinement and are not sensitive to the solution in the rest of the domain in which both grids have the same x low resolution.In fact, the improved resolution of these standing features leads to slightly less errors than are obtained by the global x low resolution grid.
To investigate effects of tensor hyperviscosity, one can compare panels Fig. 10b (scalar hyperviscosity and the highly distorted grid) and Fig. 10d (tensor hyperviscosity and the highly distorted grid).The numerical noise in the transition region present in the simulation with scalar hyperviscosity is practically eliminated when tensor formulation is used.If the simulation in Fig. 10c with a better-quality grid (low-connectivity grid) is considered optimal, then one can conclude that the tensor hyperviscosity simulation provides a very close-to-optimal result, even if a low-quality grid is used.

Convergence under grid refinement
We now present mesh convergence results for several choices of local refinement.We use our best configuration: tensor hyperviscosity running with low-connectivity grids.We compare the convergence properties of the method with uniform grids using constant-coefficient hyperviscosity, uniform grids using tensor coefficient hyperviscosity, the ×2 family of grids, the ×4 family of grids, and the ×8 family of grids.The ×2, ×4, and ×8 simulations all use tensor hyperviscosity.
We want to focus only on the spatial error, and, in all cases, use time steps in which the time truncation error is negligible, as compared to the spatial error.While running simulations for convergence, we made sure that temporal errors did not dominate.For TC2, we obtain fourth-order convergence when using time steps near the CFL limit.For TC5, we reduced time steps so that time truncation errors are of fourth order.For example, for the resolution with spatial scales x,   the time is t.For the refined grid with spatial scales x/2, we execute a simulation with t/4.The global l 2 errors for the TC2 height field for all these families of grids are shown in Fig. 11.We use the relative l 2 , defined in Williamson et al. (1992).As noted in Sect.2.1, the hyperviscosity scaling with resolution is typically chosen as s = 3.2.For TC2, we instead choose s = 4, so that the hyperviscosity term goes to 0 at a fourth-order rate and we can confirm the fourth-order accuracy of our p = 3 SEM spatial discretization.For TC2, mesh refinement adds no value to the simulation, and, as we saw in Sect.5.1, it can lead to a slight reduction in the local error, but no reduction in the error away from the refinement region.We thus examine convergence with respect to x low for each family of grids.As expected, the uniform resolution simulations with fourth-order constant-coefficient hyperviscosity demonstrate fourth-order convergence, with the tensor and scalar hyperviscosity results being nearly identical.Similar results are obtained for the ×2, ×4, and ×8 family of grids.The error is completely determined by x low for all grids, and all grids show fourth-order convergence under mesh refinement with respect to x low .Thus, the presence of mesh refinement, with refinements as much as 8×, does no harm to the global errors.
The global l 2 errors for the TC5 height field for all these families of grids are shown in Fig. 12, again plotted as a Figure 12.As in Fig. 11, except for TC5.In case of tensor-based hyperviscosity, three error curves for grid refinement with ×2, ×4 or ×8 local refinement are practically indistinguishable and obtain the same convergence rates as seen with uniform grids.For a given x low , the presence of some local refinement (×2) does decrease the global error, but ×4 or ×8 levels of local refinement produce no further benefit.2 and 5, the scalar approach had noticeable noise and oscillations near regions of local mesh refinement, which was not present with the tensor formulation.Results for both formulations were sensitive to the grid quality, as shown by comparing results on a highly distorted grid with sharp mesh transitions and a smooth grid with less acute angles due to its lower valence nodes.But, in TC2, the tensor formulation showed less grid sensitivity and obtained excellent results on both grids.
When running with tensor hyperviscosity in the SEM, the presence of local mesh refinement in TC2 had no impact on the global errors.The SEM obtained its formal order of accuracy for all grids tested (up to 8× regional refinement).In TC5, with refinement over the mountain, the presence of refinement, again, did no harm to the global errors and actually resulted in a small improvement.Asymptotically, the global errors were controlled by the coarse resolution, and the locally refined meshes obtained the same convergence rates as the global uniform meshes.
The tensor and scalar hyperviscosity operators were constructed to be resolution aware and to preserve the resolution scaling often used with uniform grids.In a high-resolution region of quasi-uniform elements, the resolution-aware operators are very similar to the SEM's well-proven constantcoefficient operator.Hence, we expect that, within the highresolution region of a locally refined mesh, the SEM can locally resolve the same types of structures that the uniform high-resolution grid can resolve.This was verified by looking at steep gradients in TC5 and comparing them to coarse and high uniform resolution solutions.
Altogether, we have

Figure 5 .
Figure 5. Different types of refined conforming quadrilateral grids.Plots (a) and (b) show close-ups of the transition regions from plots in Figs.4a and b, respectively.

O.Figure 6 .
Figure 6.Error plots for TC2.A contour spacing of 0.25 m is the same for all plots.Maximum and minimum values of error Eq. (15) are given in captions.Tensor hyperviscosity produces smoother fields comparing to scalar hyperviscosity, as follows from comparing pairs (b), (c), (d), and (e).In addition, the quality of the underlying grid significantly improves the outcome around the refined region when using scalar hyperviscosity, as seen by comparing (b) and (d).Contrary to scalar hyperviscosity, tensor hyperviscosity is more robust with respect to mesh quality, as follows from comparing (c) and (e).Simulations (c) and (e) with tensor hyperviscosity also exhibit a substantial error reduction in the vicinity of the refinement compared to the coarse uniform resolution in (a).

O.
Guba et al.: The spectral element method (SEM) on variable-resolution grids (a) Coarse resolution solution with x = 3 .(b) High resolution reference solution, using a uniform grid with x = 0.125 .

Figure 8 .
Figure 8.As in Fig. 7 but plotted in a subregion of the global domain.

( a )Figure 9 .
Figure 9. Contour plots of TC5 vorticity are shown in the top panels, while the grid is shown in the bottom panel.The contour interval is 5.0 × 10 −6 s −1 .All panels show a subregion of the global domain, which contains most of the refined region and is identical to the subregion used in Fig. 8. Panel (a) shows results using scalar hyperviscosity on the highly distorted grid shown in (c).Panel (b) shows results using tensor hyperviscosity on the low-connectivity grid shown in (d).The improved hyperviscosity and mesh quality result in significantly improved results.

Figure 10 .
Figure 10.The error in the TC5 vorticity field is plotted for the global domain.The color scheme given in (a) is the same for all plots.The tensor hyperviscosity again produces the best results with very little noise.From comparing panel (c) (tensor hyperviscosity and the lowconnectivity grid) and panel (d) (tensor hyperviscosity and the highly distorted grid), we conclude that tensor hyperviscosity is not sensitive to grid quality.

Table 2 .
Summary for simulations.TC stands for a test case, and the numbers in c 0 and s are parameters in the hyperviscosity coefficient ν = c 0 ( x) s .