Initial conditions
This module contains the initialization functions.
- SWAMPE.initial_conditions.ABCDE_init(Uic, Vic, etaic0, Phiic0, mus, I, J)
Initializes the auxiliary nonlinear components.
- Parameters:
Uic (array of float or complex) – zonal velocity component
Vic (array of float or complex) – meridional velocity component
etaic0 (array of float or complex) – initial eta
Phiic0 (array of float or complex) – initial Phi
mustile (array of float or complex) – reshaped mu array to fit the dimensions
- Returns:
J by I data arrays for eta0, eta1, delta0, delta1, and phi0, phi1
- Return type:
array of float or complex
- SWAMPE.initial_conditions.coriolismn(M, omega)
Initializes the Coriolis parameter in spectral space.
Parameters
- param M:
Spectral dimension.
- type M:
int
- param omega:
Planetary rotation rate, in radians per second.
- type omega:
float
Returns
- return:
fmn, The Coriolis parameter in spectral space.
- rtype:
array of float, size (M+1, M+1)
- SWAMPE.initial_conditions.spectral_params(M)
Generates the resolution parameters according to Table 1 and 2 from Jakob et al. (1993). Note the timesteps are appropriate for Earth-like forcing. More strongly forced planets will need shorter timesteps.
- Parameters:
M (int) – spectral resolution
Returns
- return:
N - int - the highest degree of the Legendre functions for m = 0
I - int - number of longitudes
J - int - number of Gaussian latitudes
dt - float - timestep length, in seconds
lambdas - array of float of length I - evenly spaced longitudes
mus - array of float of length J - Gaussian latitudes
w - array of float of length J - Gaussian weights
- SWAMPE.initial_conditions.state_var_init(I, J, mus, lambdas, test, etaamp, *args)
Initializes state variables.
Parameters
- param I:
number of longitudes.
- type I:
int
- param J:
number of latitudes.
- type J:
int
- param mus:
Array of Gaussian latitudes of length J.
- type mus:
array of float
- param lambdas:
Uniformly spaced longitudes of length I.
- type lambdas:
array of float
- param test:
The number of the regime being tested from Williamson et al. (1992)
- type test:
int
- param etaamp:
Amplitude of absolute vorticity.
- type etaamp:
float
- param *args:
Additional initialization parameters for tests from Williamson et al. (1992)
- type *args:
arrays of float
Returns
- return:
etaic0 - Initial condition for absolute vorticity, (J,I).
etaic1 - Second initial condition for absolute vorticity, (J,I).
deltaic0 - Initial condition for divergence, (J,I).
deltaic1 - Second initial condition for divergence, (J,I).
Phiic0 - Initial condition for geopotential, (J,I).
Phiic1 - Second initial condition for geopotential, (J,I).
- rtype:
tuple of arrays of float
- SWAMPE.initial_conditions.test1_init(a, omega, a1)
Initializes the parameters from Test 1 in Williamson et al. (1992), Advection of Cosine Bell over the Pole.
Parameters
- param a:
Planetary radius, in meters.
- type a:
float
- param omega:
Planetary rotation rate, in radians per second.
- type omega:
float
- param a1:
Angle of advection, in radians.
- type a1:
float
Returns
- return:
- SU0
Amplitude parameter from Test 1 in Williamson et al. (1992)
- sina
sine of the angle of advection.
- cosa
cosine of the angle of advection.
- etaamp
Amplitude of absolute vorticity.
- Phiamp
Amplitude of geopotential.
- rtype:
float
- SWAMPE.initial_conditions.velocity_init(I, J, SU0, cosa, sina, mus, lambdas, test)
Initializes the zonal and meridional components of the wind vector field.
Parameters
- param I:
number of latitudes.
- type I:
int
- param J:
number of longitudes.
- type J:
int
- param SU0:
Amplitude parameter from Test 1 in Williamson et al. (1992)
- type SU0:
float
- param cosa:
cosine of the angle of advection.
- type cosa:
float
- param sina:
sine of the angle of advection.
- type sina:
float
- param mus:
Array of Gaussian latitudes of length J
- type mus:
array of float
- param lambdas:
Array of uniformly spaces longitudes of length I.
- type lambdas:
array of float
- param test:
when applicable, number of test from Williamson et al. (1992).
- type test:
int
Returns
- return:
Uic - (J,I) array - the initial condition for the latitudinal velocity component,
Vic - (J,I) array - the initial condition for the meridional velocity component.
- rtype:
array of float