2. Model and Method of Numerical Calculation

A 2-dimensional non-divergent flow on a rotating sphere governed by the following vorticity equation will be examined.

Here, ψ (λ, μ, t): stream functions; ω (λ, μ, t): vertical component of vorticity (≡ Δ ψ); λ: longitude; μ: sine latitude; t: time; Ω: angular speed of rotation of sphere; ν2p: hyper-viscosity coefficient, and p: order of hyper-viscosity. Furthermore, Δ is the horizontal Laplacian and is defined as follows.

Note that,

The dependent variable ψ is expanded using spherical harmonics as


Based on the above expansion, the simultaneous ordinary differential equations for the expansion coefficient ψnm(t) is yielded from the vorticity equation using the spectral Galerkin method. The non-linear terms are calculated in real space, and a technique is adopted for returning the results back into the wavenumber space by Fourier-Legendre transformation. To eliminate aliasing, global surface domain is covered by a transform grid of 2,048 (longitudinal) × 1,024 (latitudinal) when the truncation wavenumber N = 682, and a grid of 1,024×512 when N = 341.

In the experiment for decaying turbulence, an initial turbulent field will be given to observe time evolution. The energy density is given by

The initial energy spectrum is constrained to be the distribution function employed by Cho and Polvani (1996)

and each band wavenumber component m will be given random amplitude and phase. Here, coefficient A will be set so that the total dimensionless kinetic energy Etotal ≡ Σn=2NE(n,t=0) will be equal to 1. n0 is a parameter that gives the central wavenumber, and will take the values of either 10, 50, or 100. γ is a parameter that determines the spectral width, and a larger value for γ will confine the initial energy to a narrower band. In the present study, γ was set to 1,000.

The angular velocity of rotation Ω will be used as the experimental parameter and given six values, 0, 25, 50, 100, 200, and 400. When Ω ≠ 0, the planet will make Ω/(2π) rotations during unit time.

Note here that the wavenumber at which the non-linear term of the vorticity equation and the linear rotational effect term are nearly equal (wavenumber corresponding to the Rhines scale), nβ ≡ √ (π Ω/(4U))= √ (π Ω/(4√(2× Etotal)) ) =√(π Ω/(4√ 2)), will be an important indicator. Actually, the manner of time evolution of the decaying turbulence is strongly dependent on the magnitude relation between nβ and the initial energy peak wavenumber n0. When nβ > n0, then the linear term will dominate over the non-linear term from the initial stages, and so the inverse energy cascade predicted by the turbulence theory will not always appear.

For the numerical viscosity term, p is set to 8, and ν2p = 1×10-43 (for N = 682), and 1×10-38 (for N = 341). Time integration is done using the 4th-order Runge-Kutta method, with the time increment being set to Δt = 5×10-4 (for N = 682) and 1×10-3 (for N = 341), and time evolution will be computed to t = 5.

For the numerical calculations, the efficient transform code of ISPACK (ishioka, 1999) has been employed.

<<index <previous next>