7.2 Subgrid-scale modeling

For as long as we can reasonably speculate, even the most ambitious global simulations will only resolve a small fraction of the dynamical scales which are active in the solar interior. Thus, some type of model is necessary to account for the influence of motions on scales smaller than the grid spacing.

Current subgrid-scale (SGS) models assume that this influence is merely diffusive in nature, acting as an effective scalar viscosity, thermal diffusivity, and magnetic diffusivity which are many orders of magnitude larger than the corresponding molecular values. These scalar coefficients are allowed to vary with depth and are often assumed to be proportional to $r--1/2$ , as suggested by mixing-length arguments. Such parameterizations are very crude and do not accurately represent the complex dynamics known to occur in rotating, stratified, magnetized flows (e.g., Sections 7.1 and Section 8). More realistic models are necessary in order to make substantial further progress in global simulations.

The primary objectives of a subgrid-scale model may be outlined as follows:

to reduce the influence of dissipation on the largest scales,
to reliably account for cascade processes,
to model processes which are completely unresolved,
to minimize the number of free parameters.

We now proceed to elaborate on these objectives.

The extremely high Reynolds numbers characteristic of the solar convection zone suggest that global-scale motions must be essentially inviscid (Section 5.1). Thermal and magnetic diffusion are similarly expected to be insignificant on large scales. This is not the case for current global simulations in which diffusive transport still makes a substantial contribution to the net momentum and energy balance (e.g., Figure 14 ) and still influences the generation and evolution of the magnetic field. Thus, the first goal of any successful SGS model must be to reduce the influence of this artificial dissipation.

In a spectral model, the most straightforward way to accomplish this is by imposing hyperdiffusion wherein the Laplacian diffusion operator is replaced by or supplemented with a higher-order equivalent (e.g., $4 \~/$ or $8 \~/$ ). Thus, the effective diffusion on the largest scales can be greatly reduced while maintaining an efficient dissipation on the smallest scales, preventing a buildup of energy which would otherwise cause numerical instability.

Although hyperdiffusion has benefits, it also has drawbacks. It is a practical construct with little physical justification. Furthermore, higher-order radial derivatives require additional boundary conditions in order to make the problem well-posed, placing artificial constraints on the allowable solutions. Such constraints can be avoided if hyperdiffusion is only implemented on horizontal surfaces while keeping the radial diffusion second-order, an approach which has been used in geodynamo simulations (e.g., Glatzmaier, 2002 ). However, this introduces an unphysical and largely arbitrary anisotropy into the SGS transport. Hyperviscosity also can introduce spurious overshoot near sharp gradients (related to Gibbs ringing) and may have an adverse effect on dynamo simulations, fundamentally altering the field generation process (Zhang and Schubert, 2000 ; Busse, 2000 ). It is therefore important to consider alternatives.

Turbulent flows generally exhibit cascade processes, characterized by a self-similar (scale invariant) exchange of energy or some other ideal invariant between adjacent spectral modes. The most familiar example is the forward cascade of kinetic energy which occurs within the classical inertial range of 3D, homogeneous, isotropic, incompressible turbulence (e.g., Lesieur, 1997 ; Pope, 2000 ). Rotation, stratification, and magnetism can also give rise to forward and inverse cascades (e.g., Section 8.3). By narrowing the viscous dissipation range, hyperdiffusion can extend these cascade ranges and thereby better capture the essential dynamics of the largest scales. However, the dynamics within the dissipation range is not accurately represented. A better representation of the resolved flow on all scales might be achieved by assuming from the outset that it will be self-similar on scales comparable to the grid-spacing.

A variety of self-similarity methods have been developed, as reviewed by Meneveau and Katz (2000 ). These are all based on the Large-Eddy Simulation (LES) framework whereby a low-pass filter is applied to the equations of motion (e.g., Mason, 1994 ; Pope, 2000 ). One approach, known as a dynamic SGS model, is based on the Germano identity, which relates the turbulent stress tensor, $tij$ between two self-similar scales as follows (Germano et al., 1991 ):

where

and

In these equations, overbars and brackets denote two spatial filtering operations, characterized by two different cutoff wavenumbers, $k1$ and $k2$ . The first, $k1$ , corresponds to the grid scale and the associated velocity field $-- vi$ may be regarded as the resolved flow in the simulation. The second filter is applied at a larger scale, typically chosen such that $k1 = 2k2$ . The tensors $tij$ arise when filters are applied to the Navier-Stokes equations and are often referred to as the Leonard stress.

The left-hand-side of Equation (25 ) can be evaluated directly from the resolved velocity field. However, the right-hand-side involves the unknown correlations $---- vivj$ which must be modeled (this is essentially the Reynolds stress). If some parametric form is assumed for $tij$ , Equation (25 ) may then be used to compute the parameters. For example, if the turbulent transport is assumed to be diffusive, then $t = - 2n e ij t ij$ where $nt$ is the turbulentviscosity and $eij$ is the strain rate tensor. Equation (25 ) can then be used to derive $nt$ as a function of space and time. More commonly, $nt$ itself is assumed to be proportional to the trace of $eij$ as originally proposed by Smagorinsky (Smagorinsky , 1963 ; see also Pope , 2000 ). Equation (25 ) then yields the proportionality constant (Lesieur and Métais, 1996 ; Meneveau and Katz, 2000 ). The only remaining parameter is the ratio of filter scales $k /k 1 2$ , meeting objective 4 above.

Self-similarity models such as these may in principle be applied separately for velocity, thermal, and magnetic fields and they rank among the most promising SGS approaches for solar applications (other promising strategies are reviewed by Lesieur and Métais 1996 and Foias et al. 2001). However, they do not capture nonlocal spectral transfer between large and small scales. Furthermore, they do not account for distinct small-scale dynamics such as granulation which are entirely unresolved, possessing local energy maxima on scales below the grid resolution. For this, separate models must be developed as outlined in objective 3 above. Such models may be based on local-area simulations or on parameterizations and procedures developed in the context of mean-field theory (Section 5.3). In this respect, global solar convection and dynamo simulations may ultimately resemble global circulation models (GCMs) for the Earth’s atmosphere, where unresolved processes are parameterized and where a hierarchy of modeling efforts (macroscale, mesoscale, and microscale) may be used to devise more reliable parameterizations (e.g., Beniston, 1998 ).

The most straightforward way to evaluate whether an SGS model is reliable and robust is to compare simulations with different resolutions. An intermediate-resolution simulation which incorporates the SGS model should be able to reproduce results from a higher-resolution simulation with only Laplacian diffusion. Furthermore, the LES/SGS model should eventually converge on a statistically equivalent solution as the resolution is increased. Of course, these checks will only work if the assumptions of the model are met. For example, an SGS model which relies on scale invariance will only converge if the cutoff wavenumber corresponding to the grid spacing is well within the inertial range (or some equivalent cascade range). Furthermore, as the resolution is increased, the parameterizations for previously unresolved processes (objective 3) may need to be revised as their characteristic scales begin to overlap with the dynamical range captured by the simulation. This is occurring now in GCMs where increasing the resolution does not necessarily lead to better forecasts (e.g., Williamson, 2002 ).

Large-eddy simulations with subgrid-scale modeling generally perform well in fundamental turbulence applications (Mason, 1994 ; Meneveau and Katz, 2000 ; Pope, 2000 ). Results have been promising enough that the approach has become standard in engineering and atmospheric applications. It remains to be seen how reliable they will be for solar interior dynamics. Magnetism in particular poses difficult challenges for SGS modeling which have not yet been fully addressed. Rotation and stratification (i.e., buoyancy) must also be incorporated into a realistic model. Furthermore, LES/SGS approaches can run into problems near boundaries where the characteristic scales of the flow can decrease dramatically and where qualitatively different dynamics can occur (Mason, 1994 ). This is certainly an issue in solar applications where the boundaries of the convection zone are likely to be highly complex (Section 7.3). Still, the prospects are good that improved SGS modeling may lead to substantial advances in global solar convection simulations in the near future.