xinvert package

xinvert.apps module

Created on 2022.04.13

xinvert.apps.animate_iteration(app_name, F, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan}, loop_per_frame=5, max_frames=30)[source]

Animate the iteration process.

All the invert_xxx function can be animated here. This function will add a iter dimension similar to a time dimension that shows the results at different iteration steps.

Parameters:

app_name: {‘Poisson’, ‘PV2D’, ‘GillMatsuno’, ‘Eliassen’, ‘geostrophic’, ‘StommelMunk’, ‘RefState’, ‘Omega’}: Application name as a suffix of invert_xxx.
F: xarray.DataArray: A forcing function.
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’].
coords: str, optional: Coordinate combinations in which inversion is performed.
icbc: xarray.DataArray, optional: Prescribe inital condition/guess (IC) and boundary conditions (BC).
mParams: dict, optional: Model parameters. None for the Poisson model.
iParams: dict, optional: Iteration parameters.
loop_per_frame: int, optional: Iteration loop count per frame.
max_frames: int, optional: Max frames beyond which loop is stopped.

Returns:

xarray.DataArray: The result, in which an extra dimension called iter will be added.

xinvert.apps.cal_flow(S, dims, coords='lat-lon', BCs=['fixed', 'fixed'], vtype='streamfunction', mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027})[source]

Calculate flow vector using streamfunction or velocity potential.

Parameters:

S: xarray.DataArray: Streamfunction or velocity potential.
dims: list of str: Dimension combination for the inversion e.g., [‘lat’, ‘lon’].
coords: {‘lat-lon’, ‘z-lat’, ‘z-lon’, ‘cartesian’}, optional: Coordinate combinations in which inversion is performed.
BCs: dict: Boundary conditions e.g., [‘fixed’, ‘periodic’] for 2D case.
vtype: {‘streamfunction’, ‘velocitypotential’, ‘GillMatsuno’}, optional: Type of the given variable, which determins the returns.

Returns:

tuple: Flow vector components.

xinvert.apps.default_mParams = {'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}: Application functions

xinvert.apps.invert_3DOcean(F, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting 3D ocean flow.

The 3D ocean flow equation is given as:

\[c_3 \frac{\partial^2 \psi}{\partial z^2} + c_1 \frac{\partial^2 \psi}{\partial y^2} + c_1 \frac{\partial^2 \psi}{\partial x^2} + c_3 \frac{\partial \psi}{\partial z} + c_1 \frac{\partial \psi}{\partial y} + c_1 \frac{\partial \psi}{\partial x} = F\]

Invert this equation for the streamfunction \(\psi\) given the forcing function \(F\).

Parameters:

F: xarray.DataArray

A forcing function.

dims: list

Dimension combination for the inversion e.g., [‘lev’, ‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter at south BC on beta plane.

beta: Meridional derivative of Coriolis parameter.

epsilon: Linear damping coefficient for momentum.

N2: Buoyancy frequency = g/theta0 * dtheta/dz = -R*pi/p * dtheta/dp.
k:linear damping coefficient for buoyancy

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (streamfunction) of the SOR inversion.

xinvert.apps.invert_BrethertonHaidvogel(h, dims, coords='cartesian', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the Bretherton-Haiduogel model.

The Bretherton-Haiduogel model is given as:

\[\nabla^2\psi - \lambda D \psi=-f_0-\beta y-\frac{f_0}{D}h \Phi\]

Invert this equation for the geostrophic streamfunction \(\psi\) given the topography \(h\).

Parameters:

lapPhi: xarray.DataArray

Laplacian of geopotential.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter.

beta: Meridional derivative of Coriolis parameter.

D: Constant total depth of the fluid (>> h).

lambda: Lagrangian multiplier, determined by the total energy.

iParams: dict

Iteration parameters.

Returns:

xarray.DataArray: Results (geostrophic streamfunction) of the SOR inversion.

xinvert.apps.invert_Eliassen(F, dims, coords='z-lat', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the Eliassen balanced vortex model.

The Eliassen model is given as:

\[\frac{\partial}{\partial z}\left( A\frac{\partial \psi}{\partial z} + B\frac{\partial \psi}{\partial y} \right) + \frac{\partial}{\partial y}\left( B\frac{\partial \psi}{\partial z} + C\frac{\partial \psi}{\partial y} \right) = F\]

Invert this equation for the overturning streamfunction \(\psi\) given the forcing \(F\).

Parameters:

F: xarray.DataArray

Forcing function.

dims: list

Dimension combination for the inversion e.g., [‘lev’, ‘lat’] or [‘z’,’y’].

coords: {‘z-lat’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict

Parameters required for this model are:

A: Inertial stability.

B: baroclinic stability.

C: static stability.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results of the SOR inversion.

xinvert.apps.invert_Fofonoff(F, dims, coords='cartesian', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the Fofonoff (1954) model.

The equation is given as:

\[\nabla^2 \psi - c_0 \psi = c_1 - f\]

Invert the equation for \(\psi\) given \(f\).

Parameters:

F: xarray.DataArray

Forcing function. Note that Forcing is irrelevant, only coordinates are used here. The true forcing (Coriolis parameter f) will be calculated automatically depending on the geometry of the domain.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribed inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

c0: Linear coefficient (>0).

c1: A constant.

f0: Coriolis parameter.

beta: Meridional derivative of Coriolis parameter.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results of the SOR inversion.

xinvert.apps.invert_GeoAdjustment(PV0, dims, coords='lat', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

(PV) inversion for a geostrophic adjustment model.

The balanced free surface satisfies:

\[\frac{\partial}{\partial y}\left(A\frac{\partial h}{\partial y}\right) +B h &= F\]

where

\[A = \cos\phi / f B = - q_0 \ cos\phi / g F = - f \cos\phi / g\]

Invert this equation for geostrophically adjusted free surface \(h\) given the initial PV distribution \(q_0\).

Parameters:

PV0: xarray.DataArray: Initial PV distribution.
dims: list: Dimension name for the inversion e.g., [‘lat’].
coords: {‘lat-lon’, ‘cartesian’}, optional: Coordinate combinations in which inversion is performed.
icbc: xarray.DataArray, optional: Prescribe inital condition/guess (IC) and boundary conditions (BC).
mParams: dict: No explicit model parameters are required here.
iParams: dict, optional: Iteration parameters.

Returns:

xarray.DataArray: Results (free surface after adjustment) of the SOR inversion.

xinvert.apps.invert_GillMatsuno(Q, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Gill-Matsuno model.

The Gill-Matsuno model is given as:

\[\begin{split}\epsilon u = fv - \frac{\partial \phi}{\partial x}\\\\ \epsilon v = -fu - \frac{\partial \phi}{\partial y}\\\\ \epsilon \phi + \Phi\left(\frac{\partial u}{\partial x} +\frac{\partial v}{\partial y}\right) = -Q\end{split}\]

Invert this equation for the mass distribution \(\phi\) given the diabatic heating function \(Q\).

Parameters:

Q: xarray.DataArray

Diabatic heating function.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter at south BC if on beta plane.

beta: Meridional derivative of f.

epsilon: Linear damping coefficient.

Phi: Background geopotential.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (mass distribution) of the SOR inversion.

xinvert.apps.invert_GillMatsuno_test(Q, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Gill-Matsuno model (test use only).

The Gill-Matsuno model is given as:

\[\begin{split}\epsilon u &= fv - \frac{\partial \phi}{\partial x}\\\\ \epsilon v &= -fu - \frac{\partial \phi}{\partial y}\\\\ \epsilon \phi + \Phi\left(\frac{\partial u}{\partial x} +\frac{\partial v}{\partial y}\right) &= -Q\end{split}\]

Invert this equation for the mass distribution \(\phi\) given the diabatic heating function \(Q\).

Parameters:

Q: xarray.DataArray

Diabatic heating function.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter at south BC if on beta plane.

beta: Meridional derivative of f.

epsilon: Linear damping coefficient.

Phi: Background geopotential.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (mass distribution) of the SOR inversion.

xinvert.apps.invert_MultiGrid(invert_func, *args, ratio=3, gridNo=3, **kwargs)[source]

Using multi-grid method to do the inversion (test only now).

All the invert_xxx function can be solved using multi-grid method here.

Parameters:

app_name: {‘Poisson’, ‘PV2D’, ‘GillMatsuno’, ‘Eliassen’, ‘geostrophic’, ‘StommelMunk’, ‘RefState’, ‘Omega’}: Application name as a suffix of invert_xxx.
F: xarray.DataArray: A forcing function.
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’].
coords: str, optional: Coordinate combinations in which inversion is performed.
icbc: xarray.DataArray, optional: Prescribe inital condition/guess (IC) and boundary conditions (BC).
mParams: dict, optional: Model parameters. None for the Poisson model.
iParams: dict, optional: Iteration parameters.
ratio: int, optional: Ratio of multi grids.
gridNo: int, optional: Number of multi grids.

Returns:

xarray.DataArray: The result, in which an extra dimension called iter will be added.

xinvert.apps.invert_PV2D(PV, dims, coords='z-lat', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the QG PV equation.

The QG PV equation is given as:

\[\frac{\partial}{\partial p}\left(\frac{f_0}{N^2} \frac{\partial \psi}{\partial p}\right) + \frac{1}{f_0}\frac{\partial^2 \psi}{\partial y^2} + f = q\]

This is slightly changed to:

\[L(\psi) = \frac{\partial}{\partial p}\left(\frac{f_0^2}{N^2} \frac{\partial \psi}{\partial p}\right) + \frac{\partial^2 \psi}{\partial y^2} = (q - f)f_0\]

Invert this equation for QG streamfunction \(\psi\) given the PV distribution \(q\).

Parameters:

PV: xarray.DataArray

2D distribution of QGPV.

dims: list

Dimension combination for the inversion e.g., [‘lev’, ‘lat’].

coords: {‘z-lat’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter.

beta: Coriolis parameter.

N2: buoyancy frequency.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (QG streamfunction) of the SOR inversion.

xinvert.apps.invert_Poisson(F, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the Poisson equation.

The Poisson equation is given as:

\[\frac{\partial^2 \psi}{\partial y^2} + \frac{\partial^2 \psi}{\partial x^2} = F\]

Invert the Poisson equation for \(\psi\) given \(F\).

Parameters:

F: xarray.DataArray: Forcing function (e.g., vorticity).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’].
coords: {‘lat-lon’, ‘z-lat’, ‘z-lon’, ‘cartesian’}, optional: Coordinate combinations in which inversion is performed.
icbc: xarray.DataArray, optional: Prescribed inital condition/guess (IC) and boundary conditions (BC).
mParams: dict, optional: Model parameters. None for the Poisson model.
iParams: dict, optional: Iteration parameters.

Returns:

xarray.DataArray: Results of the SOR inversion.

xinvert.apps.invert_RefState(PV, dims, coords='z-lat', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

PV inversion for a balanced symmetric vortex.

The balanced symmetric vortex equation is given as:

\[\frac{\partial}{\partial \theta}\left(\frac{2\Lambda_0}{r^3} \frac{\partial \Lambda}{\partial \theta}\right) + \frac{\partial}{\partial r}\left(\frac{\Gamma g}{Q r} \frac{\partial \Lambda}{\partial r}\right) = 0\]

Invert this equation for absolute angular momentum \(\Lambda\) given the PV distribution \(Q\).

Parameters:

PV: xarray.DataArray

2D distribution of PV.

dims: list

Dimension combination for the inversion e.g., [‘lev’, ‘lat’].

coords: {‘z-lat’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict

Parameters required for this model are:

Ang0: Angular momentum Λ0 as the known coefficient.

Gamma: vertical function defined as Γ = Rd/p * (p/p0)^κ = κ * Π/p.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (angular momentum Λ) of the SOR inversion.

xinvert.apps.invert_RefStateSWM(Q, dims, coords='lat', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

(PV) inversion for a steady state of shallow-water model.

The balanced symmetric vortex equation is given as:

\[\frac{\partial}{\partial y}\left(A\frac{\partial\Delta M}{\partial y}\right) -B\Delta M &=F\]

where

\[A = 1 / r B = \frac{2C_0 Q_0 \sin\phi}{2\pi g r^3} F = \frac{C_0^2\sin\phi }{2\pi g r^3}-\frac{2\pi\Omega^2r\sin\phi}{g}-\frac{\partial }{\partial y}\left(\frac{1}{r}\frac{\partial M_0}{\partial y}\right)\]

Invert this equation for mass correction \(\Delta M\) given the PV distribution \(Q\).

Parameters:

Q: xarray.DataArray

A set of PV contours.

C0: xarray.DataArray

Initial circulation profile along the meridional direction.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict

Parameters required for this model are:

M0: initial guess of meridional mass profile.
- c0: Eulerian zonal-mean circulation.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (angular momentum Λ) of the SOR inversion.

xinvert.apps.invert_Stommel(curl, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Stommel model.

The Stommel model is given as:

\[- \frac{R}{D}\nabla^2 \psi - \beta\frac{\partial \psi}{\partial x} = - \frac{\hat\nabla \cdot \vec\tau}{\rho_0 D}\]

Invert this equation for the streamfunction \(\psi\) given wind-stress curl \(\hat\nabla \cdot \vec\tau\).

Parameters:

curl: xarray.DataArray

Wind stress curl.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

beta: Meridional derivative of Coriolis parameter.

R: Laplacian viscosity.

D: Depth of the ocean or mixed layer depth.

rho0: Density of the fluid.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (streamfunction) of the SOR inversion.

xinvert.apps.invert_StommelArons(Q, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Stommel-Arons model.

The Stommel-Arons model is given as:

\[\begin{split}\epsilon u = fv - \frac{\partial \phi}{\partial x}\\\\ \epsilon v = -fu - \frac{\partial \phi}{\partial y}\\\\ \left(\frac{\partial u}{\partial x} +\frac{\partial v}{\partial y}\right) = -Q\end{split}\]

Invert this equation for the mass distribution \(\phi\) given the diabatic heating function \(Q\).

Parameters:

Q: xarray.DataArray

Diabatic heating function.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter at south BC if on beta plane.

beta: Meridional derivative of f.

epsilon: Linear damping coefficient for velocity only.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (mass distribution) of the SOR inversion.

xinvert.apps.invert_StommelMunk(curl, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Stommel-Munk model.

The Stommel-Munk model is given as:

\[A_4\nabla^4\psi - \frac{R}{D}\nabla^2 \psi - \beta\frac{\partial \psi}{\partial x} = - \frac{\hat\nabla \cdot \vec\tau}{\rho_0 D}\]

Invert this equation for the streamfunction \(\psi\) given wind-stress curl \(\hat\nabla \cdot \vec\tau\).

Parameters:

curl: float or xarray.DataArray

Wind-stress curl (N/m^2).

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

beta: Meridional derivative of Coriolis parameter.

A4: Hyperviscosity.

R: Laplacian viscosity.

D: Depth of the ocean or mixed layer depth.

rho0: Density of the fluid.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results of the SOR inversion.

xinvert.apps.invert_Stommel_test(curl, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting Stommel model (test used only).

The Stommel model is given as:

\[- \frac{R}{D}\nabla^2 \psi - \beta\frac{\partial \psi}{\partial x} = - \frac{\hat\nabla \cdot \vec\tau}{\rho_0 D}\]

Invert this equation for the streamfunction \(\psi\) given wind-stress curl \(\hat\nabla \cdot \vec\tau\).

Parameters:

curl: xarray.DataArray

Wind stress curl.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

beta: Meridional derivative of Coriolis parameter.

R: Laplacian viscosity.

D: Depth of the ocean or mixed layer depth.

rho0: Density of the fluid.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (streamfunction) of the SOR inversion.

xinvert.apps.invert_geostrophic(lapPhi, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the geostrophic balance model.

The geostrophic balance model is given as:

\[\frac{1}{\partial y}\left(f\frac{\partial \psi}{\partial y}\right)+ \frac{1}{\partial x}\left(f\frac{\partial \psi}{\partial x}\right)=\nabla^2 \Phi\]

Invert this equation for the geostrophic streamfunction \(\psi\) given the Laplacian of geopotential field \(\nabla^2\Phi\).

Parameters:

lapPhi: xarray.DataArray

Laplacian of geopotential.

dims: list

Dimension combination for the inversion e.g., [‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter.

beta: Meridional derivative of Coriolis parameter.

iParams: dict

Iteration parameters.

Returns:

xarray.DataArray: Results (geostrophic streamfunction) of the SOR inversion.

xinvert.apps.invert_omega(F, dims, coords='lat-lon', icbc=None, mParams={'A': 100000.0, 'N2': 0.0002, 'Omega': 7.292e-05, 'Phi': 10000.0, 'R': 5e-05, 'Rearth': 6371200.0, 'ang0': 200000.0, 'beta': 2e-11, 'c0': 8e-09, 'c1': 8e-05, 'depth': 100, 'epsilon': 7e-06, 'f0': 1e-05, 'g': 9.80665, 'lambda': 1e-08, 'rho0': 1027}, iParams={'BCs': ['fixed', 'fixed'], 'debug': False, 'mxLoop': 5000, 'optArg': None, 'printInfo': True, 'tolerance': 1e-08, 'undef': nan})[source]

Inverting the omega equation.

The omega equation is given as:

\[\frac{f^2}{N^2}\frac{\partial^2 \omega}{\partial z^2} + \frac{\partial^2 \omega}{\partial y^2} + \frac{\partial^2 \omega}{\partial x^2} = \frac{F}{N^2}\]

This is slightly changed to:

\[\frac{1}{\partial z}\left(f^2\frac{\partial \omega}{\partial z}\right)+ \frac{1}{\partial y}\left(N^2\frac{\partial \omega}{\partial y}\right)+ \frac{1}{\partial x}\left(N^2\frac{\partial \omega}{\partial x}\right)=F\]

Invert this equation for the vertical velocity \(\omega\) given the forcing function \(F\).

Parameters:

F: xarray.DataArray

A forcing function.

dims: list

Dimension combination for the inversion e.g., [‘lev’, ‘lat’, ‘lon’].

coords: {‘lat-lon’, ‘cartesian’}, optional

Coordinate combinations in which inversion is performed.

icbc: xarray.DataArray, optional

Prescribe inital condition/guess (IC) and boundary conditions (BC).

mParams: dict, optional

Parameters required for this model are:

f0: Coriolis parameter at south BC on beta plane.

beta: Meridional derivative of Coriolis parameter.

N2: Buoyancy frequency = g/theta0 * dtheta/dz = -R*pi/p * dtheta/dp.

iParams: dict, optional

Iteration parameters.

Returns:

xarray.DataArray: Results (vertical velocity) of the SOR inversion.

xinvert.core module

Created on 2020.12.09

xinvert.core.inv_general2D(A, B, C, D, E, F, G, S, dims, iParams)[source]

Inverting a 2D slice of elliptic equation in general form.

\[A \frac{\partial^2 \psi}{\partial y^2} + B \frac{\partial^2 \psi}{\partial y \partial x} + C \frac{\partial^2 \psi}{\partial x^2} + D \frac{\partial \psi}{\partial y} + E \frac{\partial \psi}{\partial x} + F \psi = G\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lat’, ‘lon’], then for the horizontal slice, the 2nd dim is ‘lat’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
D: xr.DataArray: Coefficient D.
E: xr.DataArray: Coefficient E.
F: xr.DataArray: Forcing function F.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’]. Order is important, should be consistent with the order of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_general2D_bih(A, B, C, D, E, F, G, H, I, J, S, dims, iParams)[source]

Inverting a 2D slice of elliptic equation in the general form.

\[A \frac{\partial^4 \psi}{\partial y^4} + B \frac{\partial^4 \psi}{\partial y^2 \partial x^2} + C \frac{\partial^4 \psi}{\partial x^4} + D \frac{\partial^2 \psi}{\partial y^2} + E \frac{\partial^2 \psi}{\partial y \partial x} + F \frac{\partial^2 \psi}{\partial x^2} + G \frac{\partial \psi}{\partial y} + H \frac{\partial \psi}{\partial x} + I \psi = J\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lat’, ‘lon’], then for the horizontal slice, the 2nd dim is ‘lat’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
D: xr.DataArray: Coefficient D.
E: xr.DataArray: Coefficient E.
F: xr.DataArray: Coefficient F.
G: xr.DataArray: Coefficient G.
H: xr.DataArray: Coefficient H.
I: xr.DataArray: Coefficient I.
J: xr.DataArray: Forcing function J.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’]. Order is important, should be consistent with the order of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_general3D(A, B, C, D, E, F, G, H, S, dims, iParams)[source]

Inverting a 3D volume of elliptic equation in the general form.

\[A \frac{\partial^2 \psi}{\partial z^2} + B \frac{\partial^2 \psi}{\partial y^2} + C \frac{\partial^2 \psi}{\partial x^2} + D \frac{\partial \psi}{\partial z} + E \frac{\partial \psi}{\partial y} + F \frac{\partial \psi}{\partial x} + G \psi = H\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lev’, ‘lat’, ‘lon’], then for the 3D volume, the 3rd dim is ‘lev’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
D: xr.DataArray: Coefficient D.
E: xr.DataArray: Coefficient E.
F: xr.DataArray: Coefficient F.
G: xr.DataArray: Coefficient G.
H: xr.DataArray: Forcing function H.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’]. Order is important, should be consistent with the order of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_standard1D(A, B, F, S, dims, iParams)[source]

Inverting equations in 1D standard form.

\[\frac{1}{\partial x}\left( A\frac{\partial \psi}{\partial x} + B\psi= F\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lat’], then for the meridional series, the 1st dim is ‘lat’ .

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
F: xr.DataArray: Forcing function F.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list or str: Dimension combination for the inversion e.g., [‘lat’].
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_standard2D(A, B, C, F, S, dims, iParams)[source]

Inverting equations in 2D standard form.

\[\frac{1}{\partial y}\left( A\frac{\partial \psi}{\partial y} + B\frac{\partial \psi}{\partial x} \right) + \frac{1}{\partial x}\left( B\frac{\partial \psi}{\partial y} + C\frac{\partial \psi}{\partial x} \right) = F\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lat’, ‘lon’] then for the horizontal slice, the 2nd dim is ‘lat’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
F: xr.DataArray: Forcing function F.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’]. Order is important, should be consistent with the order of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_standard2D_test(A, B, C, D, E, F, S, dims, iParams)[source]

Inverting equations in 2D standard form (test only).

\[\frac{1}{\partial y}\left( A\frac{\partial \psi}{\partial y} + B\frac{\partial \psi}{\partial x} \right) + \frac{1}{\partial x}\left( B\frac{\partial \psi}{\partial y} + C\frac{\partial \psi}{\partial x} \right) + E\psi= F\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lat’, ‘lon’], then for the horizontal slice, the 2nd dim is ‘lat’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
D: xr.DataArray: Coefficient D.
E: xr.DataArray: Coefficient E.
F: xr.DataArray: Forcing function F.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lat’, ‘lon’]. Order is important, should be consistent with the order of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.core.inv_standard3D(A, B, C, F, S, dims, iParams)[source]

Inverting a 3D volume of elliptic equation in a standard form.

\[\frac{1}{\partial z}\left(A\frac{\partial \omega}{\partial z}\right)+ \frac{1}{\partial y}\left(B\frac{\partial \omega}{\partial y}\right)+ \frac{1}{\partial x}\left(C\frac{\partial \omega}{\partial x}\right)=F\]

Invert this equation using SOR iteration. If F = F[‘time’, ‘lev’, ‘lat’, ‘lon’] and we invert for the 3D spatial distribution, then 3rd dim is ‘lev’, 2nd dim is ‘lat’ and 1st dim is ‘lon’.

Parameters:

A: xr.DataArray: Coefficient A.
B: xr.DataArray: Coefficient B.
C: xr.DataArray: Coefficient C.
F: xr.DataArray: Forcing function F.
S: xr.DataArray: Initial guess of the solution (also the output).
dims: list: Dimension combination for the inversion e.g., [‘lev’, ‘lat’, ‘lon’]. Order is important, should be consistent with the dimensions of F.
iParams: dict: Parameters for inversion.

Returns:

xarray.DataArray: Solution \(\psi\).

xinvert.finitediffs module

Created on 2021.01.03

class xinvert.finitediffs.FiniteDiff(dim_mapping, BCs='extend', coords='lat-lon', fill=0, R=6371200.0)[source]

Bases: object

This class wrap some basic finite-difference operators supported for Cartesian coordinates (coords=’cartesian’) or latitude/longitude coordinates (coords=’lat/lon’), using centered different scheme.

This is designed particularly for Arakawa A grid (all the variables are defined on the same grid points). For grids of other types (variables are staggered), please use xgcm to calculate the finite difference in finite volumn fashion.

For derivative along a dimension, one may use xarray’s differentiate(). The problem with xarray’s differentiate() is that the boundary conditions are not flexible enough for our purpose. So we implement each operator here using xr.DataArray.pad() method to account for different BCs.

Attributes:

dmap: dict: Dimension mapping from those in xarray.DataArray to [‘T’, ‘Z’, ‘Y’, ‘X’].
BCs: dict: Default boundary conditions e.g., BCs={‘X’: ‘periodic’} for both end points along ‘X’ dimension; or BCs={‘Z’: (‘fixed’,’reflect’)} for fixed left BC and reflected right BC. Left indicates lower indices along ‘X’.
fill: float or dict: Fill value if BCs are fixed. A typical example can be: {‘Z’:(1,2), ‘Y’:(0,0), ‘X’:(1,0)}
coords: {‘lat-lon’, ‘cartesian’}: Types of coords. Should be one of [‘lat-lon’, ‘cartesian’].
R: float: Radius of Earth.

Methods

grad(scalar)	3D gradient.
divg(vector)	3D divergence.
vort(vector)	3D vorticity.
curl(vector)	vertical vorticity of vector.
laplacian(scalar)	Laplacian.
tension_strain(vector)	Tension strain.
shear_strain(vector)	Shear strain.
deformation_rate(vector)	Deformation rate.
Okubo_Weiss(vector)	Okubo Weiss parameter.

Laplacian(v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate the Laplacian of a scalar.

Parameters:

v: xarray.DataArray: A given scale variable.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: Laplacian of a scalar.

Okubo_Weiss(u, v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate Okubo-Weiss parameter.

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: Okubo-Weiss parameter.

curl(u, v, BCs=None, fill=None)[source]

Calculate vertical vorticity (k) component.

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

deformation_rate(u, v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate sqrt(tension^2+shear^2).

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: Deformation rate.

divg(vector, dims, BCs=None, fill=None)[source]

Calculate divergence as du/dx + dv/dy + dw/dz.

Parameters:

vector: xarray.DataArray or a list (tuple) of xarray.DataArray: Component(s) of a vector.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: Divergence.

Examples

>>> divX   = divg(u, 'X')            # du/dx
>>> divY   = divg(v, 'Y')            # dv/dy
>>> divZ   = divg(w, 'Z')            # dw/dz
>>> divXY  = divg((u,v), ['X','Y'])  # du/dx+dv/dy
>>> divVW  = divg((v,w), ['Y','Z'])  # dv/dy+dw/dz
>>> divXZ  = divg((u,w), ['X','Z'])  # du/dx+dw/dz
>>> divXYZ = divg((u,v,w), ['X','Y','Z'])  # du/dx+dv/dy+dw/dz

grad(v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate spatial gradient components along each dimension given.

Parameters:

v: xarray.DataArray: A scalar variable.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray or list: Gradient components.

Examples

>>> vx = grad(v, 'X', coords='cartesian')
>>> vx, vy = grad(v, ['lon', 'lat'])
>>> vz, vy, vx = grad(v, ['lev', 'lat', 'lon'])

shear_strain(u, v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate tension strain as dv/dx + du/dy.

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: shear strain.

tension_strain(u, v, dims=['X', 'Y'], BCs=None, fill=None)[source]

Calculate tension strain as du/dx - dv/dy.

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
dims: list of str: Dimensions for gradient. Order of dims is the same as that of the outputs: vx, vy, vz = grad(v, [‘X’, ‘Y’, ‘Z’]). Here use [‘X’, ‘Y’] in dim_mapping instead of true dimension names.
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray: tension strain.

vort(u=None, v=None, w=None, components='k', BCs=None, fill=None)[source]

Calculate vorticity component. All the three components satisfy the right-hand rule so that we only need one function and input different components accordingly.

Parameters:

u: xarray.DataArray: X-component velocity.
v: xarray.DataArray: Y-component velocity.
w: xarray.DataArray: Z-component velocity.
components: str or list of str: Component(s) of the vorticity. Order of component is the same as that of the outputs: vork, vorj, vori = vort(u,v,w, [‘k’,’j’,’i’])
BCs: dict: Boundary condition. If provided, overwrite the default one per call.
fill: tuple of floats: Fill values of fixed BC. If provided, overwrite the default one per call.

Returns:

xarray.DataArray or list: vorticity components.

Examples

>>> vori = vort(v=v, w=w, 'i')   # x-component (i) is: dw/dy - dv/dz
>>> vorj = vort(u=u, w=w, 'j')   # y-component (j) is: du/dz - dw/dx
>>> vork = vort(u=u, v=v, 'k')   # z-component (k) is: dv/dx - du/dy

>>> vori,vorj      = vort(u=u,v=v,w=w, ['i','j'])      # i,j components
>>> vori,vorj,vork = vort(u=u,v=v,w=w, ['i','j','k'])  # all components

xinvert.finitediffs.deriv(v, dim, BCs=('extend', 'extend'), fill=(0, 0), scale=1, scheme='center')[source]

First-order derivative along a given dimension.

The first-order derivative is calculated with proper boundary conditions (BCs) and finite difference scheme.

Parameters:

v: xarray.DataArray

A scalar variable.

dim: str

Dimension along which difference is taken.

BCs: tuple or list of str

Boundary conditions for the two end points. Types of BCs are:

‘fixed’: pad with fixed value.
‘extend’: pad with original edge value.
‘reflect: pad with first inner value. 1st-order derivative is exactly zero.
‘extrapolate’: pad with extrapolated value. (NOT implemented). 1st-order derivative equals the first inner point. 2nd-order derivative is exactly zero.
‘periodic’: pad with cyclic values.

fill: tuple of floats

Fill values as BCs if BC is fixed at two end points.

scale: float or xarray.DataArray

Scale of the result, usually the metric of the dimension.

scheme: str

Finite difference scheme in [‘center’, ‘forward’, ‘backward’]

Returns:

xarray.DataArray: First-order derivative along the dimension

xinvert.finitediffs.deriv2(v, dim, BCs=('extend', 'extend'), fill=(0, 0), scale=1)[source]

Second-order derivative along a given dimension

The second-order derivative is calculated with proper boundary conditions (BCs) and centered finite-difference scheme.

Parameters:

v: xarray.DataArray

A scalar variable.

dim: str

Dimension along which difference is taken.

BCs: tuple or list of str

Boundary conditions for the two end points. Types of BCs are:

‘fixed’: pad with fixed value.
‘extend’: pad with original edge value.
‘reflect: pad with first inner value. 1st-order derivative is exactly zero.
‘extrapolate’: pad with extrapolated value. (NOT implemented). 1st-order derivative equals the first inner point. 2nd-order derivative is exactly zero.
‘periodic’: pad with cyclic values.

fill: tuple of floats

Fill values as BCs if BC is fixed at two end points.

scale: float or xarray.DataArray

Scale of the result, usually the metric of the dimension.

scheme: {‘center’, ‘forward’, ‘backward’}

Finite difference scheme in [‘center’, ‘forward’, ‘backward’].

Returns:

xarray.DataArray: Second-order derivative along the dimension

xinvert.finitediffs.padBCs(v, dim, BCs, fill=(0, 0))[source]

Pad array with boundary conditions.

Pad (add two extra endpoints) original DataArray with BCs along a specific dimension. Types of boundary conditions are:

‘fixed’: pad with fixed value.

‘extend’: pad with original edge value.

‘reflect: pad with first inner value. 1st-order derivative is exactly zero.

‘extrapolate’: pad with extrapolated value. (NOT implemented). 1st-order derivative equals the first inner point. 2nd-order derivative is exactly zero.

‘periodic’: pad with cyclic values.

Parameters:

v: xarray.DataArray: A scalar variable.
dim: list of str: Dimension along which it is padded.
BCs: tuple or list of str: Boundary conditions for the two end points e.g., (‘fixed’,’fixed’).
fill: tuple or list of floats: Fill values as BCs if BC is fixed at two end points e.g., (0,0).

Returns:

p: xarray.DataArray: Padded array.

xinvert.numbas module

Created on 2020.12.09

xinvert.numbas.absNorm1D(S, undef)[source]: Sum up 1D absolute value S

xinvert.numbas.absNorm2D(S, undef)[source]: Sum up 2D absolute value S

xinvert.numbas.absNorm3D(S, undef)[source]: Sum up 3D absolute value S

xinvert.numbas.invert_general_2D(S, A, B, C, D, E, F, G, yc, xc, dely, delx, BCy, BCx, delxSqr, ratio, ratioQtr, ratioSqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 2D slice of elliptic equation in the general form.

\[A \frac{\partial^2 \psi}{\partial y^2} + B \frac{\partial^2 \psi}{\partial y \partial x} + C \frac{\partial^2 \psi}{\partial x^2} + D \frac{\partial \psi}{\partial y} + E \frac{\partial \psi}{\partial x} + F \psi = G\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first term.
B: numpy.array: Coefficient for the second term.
C: numpy.array: Coefficient for the third term.
D: numpy.array: Coefficient for the fourth term.
E: numpy.array: Coefficient for the fifth term.
F: numpy.array: Coefficient for the sixth term.
G: numpy.array: Known forcing function.
yc: int: Number of grid point in the y-dimension (e.g., Y or lat).
xc: int: Number of grid point in the x-dimension (e.g., X or lon).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’, ‘periodic’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension y (unit of m^2).
ratio: float: Ratio of delx to dely.
ratioQtr: float: Ratio of delx to dely, divided by 4.
ratioSqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

S: numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_general_3D(S, A, B, C, D, E, F, G, H, zc, yc, xc, delz, dely, delx, BCz, BCy, BCx, delxSqr, ratio2, ratio1, ratio2Sqr, ratio1Sqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 3D volume of elliptic equation in the general form.

\[A \frac{\partial^2 \psi}{\partial z^2} + B \frac{\partial^2 \psi}{\partial y^2} + C \frac{\partial^2 \psi}{\partial x^2} + D \frac{\partial \psi}{\partial z} + E \frac{\partial \psi}{\partial y} + F \frac{\partial \psi}{\partial x} + G \psi = H\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first term.
B: numpy.array: Coefficient for the second term.
C: numpy.array: Coefficient for the third term.
D: numpy.array: Coefficient for the fourth term.
E: numpy.array: Coefficient for the fifth term.
F: numpy.array: Coefficient for the sixth term.
G: numpy.array: Coefficient for the seventh term.
H: numpy.array: A known forcing function.
zc: int: Number of grid point in the z-dimension (e.g., Z or lev).
yc: int: Number of grid point in the y-dimension (e.g., Y or lat).
xc: int: Number of grid point in the x-dimension (e.g., X or lon).
delz: float: Increment (interval) in dimension z (unit of m, not degree).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCz: str: Boundary condition for dimension z in [‘fixed’, ‘extend’].
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension y (unit of m^2).
ratio2Sqr: float: Squared Ratio of delx to delz.
ratio1Sqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_general_bih_2D(S, A, B, C, D, E, F, G, H, I, J, yc, xc, dely, delx, BCy, BCx, delxSSr, delxTr, delxSqr, ratio, ratioSSr, ratioQtr, ratioSqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 2D slice of biharmonic equation in the general form.

\[A \frac{\partial^4 \psi}{\partial y^4} + B \frac{\partial^4 \psi}{\partial y^2 \partial x^2} + C \frac{\partial^4 \psi}{\partial x^4} + D \frac{\partial^2 \psi}{\partial y^2} + E \frac{\partial^2 \psi}{\partial y \partial x} + F \frac{\partial^2 \psi}{\partial x^2} + G \frac{\partial \psi}{\partial y} + H \frac{\partial \psi}{\partial x} + I \psi = J\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first term.
B: numpy.array: Coefficient for the second term.
C: numpy.array: Coefficient for the third term.
D: numpy.array: Coefficient for the fourth term.
E: numpy.array: Coefficient for the fifth term.
F: numpy.array: Coefficient for the sixth term.
G: numpy.array: Coefficient for the seventh term.
H: numpy.array: Coefficient for the eighth term.
I: numpy.array: Coefficient for the ninth term.
J: numpy.array: Known forcing function.
yc: int: Number of grid point in the y-dimension (e.g., Y or lat).
xc: int: Number of grid point in the x-dimension (e.g., X or lon).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’, ‘periodic’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSSr: float: Increment (interval) in dimension y (unit of m^2) to the power of 4.
delxTr: float: cubed increment (interval) in dimension y (unit of m^2).
delxSqr: float: Squared increment (interval) in dimension y (unit of m^2).
ratio: float: Ratio of delx to dely.
ratioSSr: float: Ratio of delx to dely to the power of 4.
ratioQtr: float: Ratio of delx to dely, divided by 4.
ratioSqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_standard_1D(S, A, B, F, xc, delx, BCx, delxSqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 1D series of elliptic equation in standard form.

\[\frac{1}{\partial x}\left( A\frac{\partial \psi}{\partial x}\right) + B\psi = F\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the 2nd-order derivative.
B: numpy.array: Coefficient for the linear term.
F: numpy.array: Forcing function.
xc: int: Number of grid point in x-dimension (e.g., X or lon).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension x (unit of m^2).
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

S: numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_standard_2D(S, A, B, C, F, yc, xc, dely, delx, BCy, BCx, delxSqr, ratioQtr, ratioSqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 2D slice of elliptic equation in standard form.

\[\frac{1}{\partial y}\left( A\frac{\partial \psi}{\partial y} + B\frac{\partial \psi}{\partial x} \right) + \frac{1}{\partial x}\left( B\frac{\partial \psi}{\partial y} + C\frac{\partial \psi}{\partial x} \right) = F\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first dimensional derivative.
B: numpy.array: Coefficient for the cross derivatives.
C: numpy.array: Coefficient for the second dimensional derivative.
F: numpy.array: Forcing function.
yc: int: Number of grid point in y-dimension (e.g., Y or lat).
xc: int: Number of grid point in x-dimension (e.g., X or lon).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’, ‘periodic’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension x (unit of m^2).
ratioQtr: float: Ratio of delx to dely, divided by 4.
ratioSqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_standard_2D_test(S, A, B, C, D, E, F, yc, xc, dely, delx, BCy, BCx, delxSqr, ratioQtr, ratioSqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 2D slice of elliptic equation in standard form.

\[\frac{1}{\partial y}\left( A\frac{\partial \psi}{\partial y} + B\frac{\partial \psi}{\partial x} \right) + \frac{1}{\partial x}\left( C\frac{\partial \psi}{\partial y} + D\frac{\partial \psi}{\partial x} \right) + E\psi = F\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first dimensional derivative.
B: numpy.array: Coefficient for the cross derivatives.
C: numpy.array: Coefficient for the cross derivatives.
D: numpy.array: Coefficient for the second dimensional derivative.
E: numpy.array: Coefficient for the linear term.
F: numpy.array: Forcing function.
yc: int: Number of grid point in y-dimension (e.g., Y or lat).
xc: int: Number of grid point in x-dimension (e.g., X or lon).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’, ‘periodic’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension x (unit of m^2).
ratioQtr: float: Ratio of delx to dely, divided by 4.
ratioSqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

S: numpy.array: Results of the SOR inversion.

xinvert.numbas.invert_standard_3D(S, A, B, C, F, zc, yc, xc, delz, dely, delx, BCz, BCy, BCx, delxSqr, ratio2Sqr, ratio1Sqr, optArg, undef, flags, mxLoop, tolerance)[source]

Inverting a 3D volume of elliptic equation in standard form.

\[\frac{1}{\partial z}\left(A\frac{\partial \omega}{\partial z}\right)+ \frac{1}{\partial y}\left(B\frac{\partial \omega}{\partial y}\right)+ \frac{1}{\partial x}\left(C\frac{\partial \omega}{\partial x}\right)=F\]

Parameters:

S: numpy.array (output): Results of the SOR inversion.
A: numpy.array: Coefficient for the first dimensional derivative.
B: numpy.array: Coefficient for the cross derivatives.
C: numpy.array: Coefficient for the second dimensional derivative.
F: numpy.array: Forcing function.
zc: int: Number of grid point in z-dimension (e.g., Z or lev).
yc: int: Number of grid point in y-dimension (e.g., Y or lat).
xc: int: Number of grid point in x-dimension (e.g., X or lon).
delz: float: Increment (interval) in dimension z (unit of m or Pa).
dely: float: Increment (interval) in dimension y (unit of m, not degree).
delx: float: Increment (interval) in dimension x (unit of m, not degree).
BCz: str: Boundary condition for dimension z in [‘fixed’, ‘extend’].
BCy: str: Boundary condition for dimension y in [‘fixed’, ‘extend’, ‘periodic’].
BCx: str: Boundary condition for dimension x in [‘fixed’, ‘extend’, ‘periodic’].
delxSqr: float: Squared increment (interval) in dimension x (unit of m^2).
ratio2Sqr: float: Squared Ratio of delx to delz.
ratio1Sqr: float: Squared Ratio of delx to dely.
optArg: float: Optimal argument ‘omega’ (relaxation factor between 1 and 2) for SOR.
undef: float: Undefined value.
flags: numpy.array: Length of 3 array, [0] is flag for overflow, [1] for converge speed and [2] for how many loops used for iteration.
mxLoop: int: Maximum loop count, larger than this will break the iteration.
tolerance: float: Tolerance for iteraction, smaller than this will break the iteraction.

Returns:

numpy.array: Results of the SOR inversion.

xinvert.numbas.trace(a, b, c, d)[source]

Trace method for solving tri-diagonal equation set.

Parameters:

a: numpy.array: Lower coefficients of the matrix (N-1).
b: numpy.array: Diagonal coefficients of the matrix (N).
c: numpy.array: Upper coefficients of the matrix (N-1).
d: numpy.array: Vector on the right-hand side of the equation (N).

Returns:

numpy.array: Results of the unknown (N).

xinvert.numbas.traceCyclic(a, b, c, d, a0, cn)[source]

Trace method for solving tri-diagonal equation set with periodic BCs.

Parameters:

a: numpy.array: Lower coefficients of the matrix (N-1).
b: numpy.array: Diagonal coefficients of the matrix (N).
c: numpy.array: Upper coefficients of the matrix (N-1).
d: numpy.array: Vector on the right-hand side of the equation (N).
a0: float: Cyclic coefficient for a.
cn: float: Cyclic coefficient for c.

Returns:

numpy.array: Results of the unknown (N).

xinvert.utils module

Created on 2021.01.03

xinvert.utils.loop_noncore(data, dims=None)[source]

Loop over the non-core dimensions using generator.

The non-core dimensions are given outside the list in dims.

Parameters:

data: xarray.DataArray: A given multidimensional data.
dims: list of str: Core dimensions. The remaining dimensions are non-core dimension

Yields:

dict: dict indicates the portion of data can be extracted