API

Base Model Class

This is the base class from which all other models inherit. All of these initialization arguments are available to all of the other model types. This class is not called directly.

class pyqg.Model(nx=64, ny=None, L=1000000.0, W=None, dt=7200.0, twrite=1000.0, tmax=1576800000.0, tavestart=315360000.0, taveint=86400.0, f=10000.0, hb=None, useAB2=False, rek=5.787e-07, filterfac=23.6, diagnostics_list='all', ntd=1, quiet=False, logfile=None)

A generic pseudo-spectral inversion model.

Attributes

q (real array) Potential vorticity in real space
qh (complex array) Potential vorticity in spectral space
ph (complex array) Streamfunction in spectral space
u, v (real arrays) Velocity anomaly components in real space
ufull, vfull (real arrays) Full velocity components in real space
uh, vh (complex arrays) Velocity anomaly components in spectral space
nx, ny (int) Number of grid points in the x and y directions
L, W (float) Domain length in x and y directions
rek (float) Linear drag in lower layer
filterfac (float) Amplitdue of the spectral spherical filter
dt (float) Numerical timstep
twrite (int) Interval for cfl writeout (units: number of timesteps)
tmax (float) Total time of integration (units: model time)
tavestart (float) Start time for averaging (units: model time)
tsnapstart (float) Start time for snapshot writeout (units: model time)
taveint (float) Time interval for accumulation of diagnostic averages. (units: model time)
tsnapint (float) Time interval for snapshots (units: model time)
ntd (int) Number of threads to use. Should not exceed the number of cores on your machine.

Note

All of the test cases use nx==ny. Expect bugs if you choose these parameters to be different.

Note

All time intervals will be rounded to nearest dt interval.

Parameters:

nx : int

Number of grid points in the x direction.

ny : int

Number of grid points in the y direction (default: nx).

L : number

Domain length in x direction. Units: meters.

W :

Domain width in y direction. Units: meters (default: L).

rek : number

linear drag in lower layer. Units: seconds -1.

filterfac : number

amplitdue of the spectral spherical filter (originally 18.4, later changed to 23.6).

dt : number

Numerical timstep. Units: seconds.

twrite : int

Interval for cfl writeout. Units: number of timesteps.

tmax : number

Total time of integration. Units: seconds.

tavestart : number

Start time for averaging. Units: seconds.

tsnapstart : number

Start time for snapshot writeout. Units: seconds.

taveint : number

Time interval for accumulation of diagnostic averages. Units: seconds. (For performance purposes, averaging does not have to occur every timestep)

tsnapint : number

Time interval for snapshots. Units: seconds.

ntd : int

Number of threads to use. Should not exceed the number of cores on your machine.

describe_diagnostics()

Print a human-readable summary of the available diagnostics.

run()

Run the model forward without stopping until the end.

run_with_snapshots(tsnapstart=0.0, tsnapint=432000.0)

Run the model forward, yielding to user code at specified intervals.

Parameters:

tsnapstart : int

The timestep at which to begin yielding.

tstapint : int

The interval at which to yield.

spec_var(ph)

compute variance of p from Fourier coefficients ph

stability_analysis(bottom_friction=False)

Baroclinic instability analysis

vertical_modes()

Calculate standard vertical modes. Simply the eigenvectors of the stretching matrix S

Specific Model Types

These are the actual models which are run.

class pyqg.QGModel(beta=1.5e-11, rd=15000.0, delta=0.25, H1=500, U1=0.025, U2=0.0, V1=0.0, V2=0.0, **kwargs)

Two layer quasigeostrophic model.

This model is meant to representflows driven by baroclinic instabilty of a base-state shear \(U_1-U_2\). The upper and lower layer potential vorticity anomalies \(q_1\) and \(q_2\) are

\[\begin{split}q_1 &= \nabla^2\psi_1 + F_1(\psi_2 - \psi_1) \\ q_2 &= \nabla^2\psi_2 + F_2(\psi_1 - \psi_2)\end{split}\]

with

\[\begin{split}F_1 &\equiv \frac{k_d^2}{1 + \delta^2} \\ F_2 &\equiv \delta F_1 \ .\end{split}\]

The layer depth ratio is given by \(\delta = H_1 / H_2\). The total depth is \(H = H_1 + H_2\).

The background potential vorticity gradients are

\[\begin{split}\beta_1 &= \beta + F_1(U_1 - U_2) \\ \beta_2 &= \beta - F_2( U_1 - U_2) \ .\end{split}\]

The evolution equations for \(q_1\) and \(q_2\) are

\[\begin{split}\partial_t\,{q_1} + J(\psi_1\,, q_1) + \beta_1\, {\psi_1}_x &= \text{ssd} \\ \partial_t\,{q_2} + J(\psi_2\,, q_2)+ \beta_2\, {\psi_2}_x &= -r_{ek}\nabla^2 \psi_2 + \text{ssd}\,.\end{split}\]

where ssd represents small-scale dissipation and \(r_{ek}\) is the Ekman friction parameter.

Parameters:

beta : number

Gradient of coriolis parameter. Units: meters -1 seconds -1

rek : number

Linear drag in lower layer. Units: seconds -1

rd : number

Deformation radius. Units: meters.

delta : number

Layer thickness ratio (H1/H2)

U1 : number

Upper layer flow. Units: m/s

U2 : number

Lower layer flow. Units: m/s

layer2modal()

calculate modal streamfunction and PV

set_U1U2(U1, U2)

Set background zonal flow.

Parameters:

U1 : number

Upper layer flow. Units: m/s

U2 : number

Lower layer flow. Units: m/s

set_q1q2(q1, q2, check=False)

Set upper and lower layer PV anomalies.

Parameters:

q1 : array-like

Upper layer PV anomaly in spatial coordinates.

q1 : array-like

Lower layer PV anomaly in spatial coordinates.

class pyqg.BTModel(beta=0.0, rd=0.0, H=1.0, U=0.0, V=0.0, **kwargs)

Single-layer (barotropic) quasigeostrophic model. This class can represent both pure two-dimensional flow and also single reduced-gravity layers with deformation radius rd.

The equivalent-barotropic quasigeostrophic evolution equations is

\[\partial_t q + J(\psi, q ) + \beta \psi_x = \text{ssd}\]

The potential vorticity anomaly is

\[q = \nabla^2 \psi - \kappa_d^2 \psi\]
Parameters:

beta : number, optional

Gradient of coriolis parameter. Units: meters -1 seconds -1

rd : number, optional

Deformation radius. Units: meters.

U : number, optional

Upper layer flow. Units: meters.

set_UV(U, V)

Set background zonal flow.

Parameters:

U : number

Upper layer flow. Units meters.

class pyqg.SQGModel(beta=0.0, Nb=1.0, rd=0.0, H=1.0, U=0.0, V=0.0, **kwargs)

Surface quasigeostrophic model.

Parameters:

beta : number

Gradient of coriolis parameter. Units: meters -1 seconds -1

Nb : number

Buoyancy frequency. Units: seconds -1.

U : number

Background zonal flow. Units: meters.

set_UV(U, V)

Set background zonal flow

Lagrangian Particles

class pyqg.LagrangianParticleArray2D(x0, y0, periodic_in_x=False, periodic_in_y=False, xmin=-inf, xmax=inf, ymin=-inf, ymax=inf, particle_dtype='f8')

A class for keeping track of a set of lagrangian particles in two-dimensional space. Tries to be fast.

Parameters:

x0, y0 : array-like

Two arrays (same size) representing the particle initial positions.

periodic_in_x : bool

Whether the domain wraps in the x direction.

periodic_in_y : bool

Whether the domain ‘wraps’ in the y direction.

xmin, xmax : numbers

Maximum and minimum values of x coordinate

ymin, ymax : numbers

Maximum and minimum values of y coordinate

particle_dtype : dtype

Data type to use for particles

step_forward_with_function(uv0fun, uv1fun, dt)

Advance particles using a function to determine u and v.

Parameters:

uv0fun : function

Called like uv0fun(x,y). Should return the velocity field u, v at time t.

uv1fun(x,y) : function

Called like uv1fun(x,y). Should return the velocity field u, v at time t + dt.

dt : number

Timestep.

class pyqg.GriddedLagrangianParticleArray2D(x0, y0, Nx, Ny, grid_type='A', **kwargs)

Lagrangian particles with velocities given on a regular cartesian grid.

Parameters:

x0, y0 : array-like

Two arrays (same size) representing the particle initial positions.

Nx, Ny: int

Number of grid points in the x and y directions

grid_type: {‘A’}

Arakawa grid type specifying velocity positions.

interpolate_gridded_scalar(x, y, c, order=1, pad=1, offset=0)

Interpolate gridded scalar C to points x,y.

Parameters:

x, y : array-like

Points at which to interpolate

c : array-like

The scalar, assumed to be defined on the grid.

order : int

Order of interpolation

pad : int

Number of pad cells added

offset : int

???

Returns:

ci : array-like

The interpolated scalar

step_forward_with_gridded_uv(U0, V0, U1, V1, dt, order=1)

Advance particles using a gridded velocity field. Because of the Runga-Kutta timestepping, we need two velocity fields at different times.

Parameters:

U0, V0 : array-like

Gridded velocity fields at time t - dt.

U1, V1 : array-like

Gridded velocity fields at time t.

dt : number

Timestep.

order : int

Order of interpolation.

Diagnostic Tools

Utility functions for pyqg model data.

pyqg.diagnostic_tools.calc_ispec(model, ph)

Compute isotropic spectrum phr of ph from 2D spectrum.

Parameters:

model : pyqg.Model instance

The model object from which ph originates

ph : complex array

The field on which to compute the variance

Returns:

kr : array

isotropic wavenumber

phr : array

isotropic spectrum

pyqg.diagnostic_tools.spec_var(model, ph)

Compute variance of p from Fourier coefficients ph.

Parameters:

model : pyqg.Model instance

The model object from which ph originates

ph : complex array

The field on which to compute the variance

Returns:

var_dens : float

The variance of ph