OceanTransportMatrixBuilder

A Julia package to build ocean transport matrices from CMIP model output.

Warning

This is work in progress. Breaking changes expected.

The purpose of this package is to build transport matrices from standard CMIP model output as part of the Transport Matrix Intercomparison Project (TMIP).

By rearranging the 3D grid of the ocean into a vector, the divergence of the flow of any tracer can be conveniently expressed in matrix form. That is, for the vector x representing the 3D tracer concentrations x(r) at locations r, the divergence of the flow is linear in x and can be represented by T x as a transport matrix T acting on x. The code of this package is an extension of the work of Matt Chamberlain, who built transport matrices from the ACCESS1.3 model¹, which has been successfully used in multiple following projects²³⁴⁵.

The main application driving this project is for the validation of marine Carbon Dioxide Removal (mCDR) by computing the timescales and pathways for water in the deep ocean to reemerge to the surface⁶⁷. However, these matrices are useful in a number of contexts, e.g., for avoiding spin ups, optimization, or novel diagnostics⁸. The motivation for sharing this package is thus to facilitate the use of novel diagnostics across CMIP models.

Example use

The intended use would be simply to feed the ocean transport output from a given CMIP model to the functions in this package to build the desired transport matrix. After the required NetCDF files have been created (see some preliminary but working Python code for creating such files in TMIP-code/notebooks)

using NetCDF
using YAXArrays
using OceanTransportMatrixBuilder

inputdir = "/Users/benoitpasquier/Data/TMIP/data/ACCESS-ESM1-5/historical/r1i1p1f1/Jan1990-Dec1999"
# The path above is for my local machine. See test/onlinebuild.jl for an example using data from the cloud.

# Load datasets
umo_ds = open_dataset(joinpath(inputdir, "umo.nc"))
vmo_ds = open_dataset(joinpath(inputdir, "vmo.nc"))
mlotst_ds = open_dataset(joinpath(inputdir, "mlotst.nc"))
volcello_ds = open_dataset(joinpath(inputdir, "volcello.nc"))
areacello_ds = open_dataset(joinpath(inputdir, "areacello.nc"))

# Load the required variables
# Note: mass transport is expected to live on the (u,v) coordinates of an Arakawa C-grid
umo = umo_ds.umo # mass transport across "east" cell face
vmo = vmo_ds.vmo # mass transport across "north" cell face
mlotst = mlotst_ds.mlotst # MLD
areacello = areacello_ds.areacello # top area of cells
volcello = volcello_ds.volcello # volume of cells
lon = volcello_ds.lon # longitude of cell centers
lat = volcello_ds.lat # latitude of cell centers
lev = volcello_ds.lev # depth of cell centers
lon_vertices = volcello_ds.lon_verticies # cell vertices
lat_vertices = volcello_ds.lat_verticies # cell vertices
# "verticies" name is an xmip bug: https://github.com/jbusecke/xMIP/issues/369
# (my local data was preprocessed with xmip)

# Make the required data from grid geometry
gridmetrics = makegridmetrics(; areacello, volcello, lon, lat, lev, lon_vertices, lat_vertices)

# Make the indices for going back and forth between 3D and 1D
indices = makeindices(gridmetrics.v3D)

# Make arrays of the flux on each face for each grid cell
ϕ = facefluxesfrommasstransport(; umo, vmo, gridmetrics, indices)

# Some parameter values
ρ = 1035.0    # density (kg/m^3)
κH = 500.0    # horizontal diffusivity (m^2/s)
κVML = 0.1    # mixed-layer vertical diffusivity (m^2/s)
κVdeep = 1e-5 # background vertical diffusivity (m^2/s)

# Make the transport matrix (should take a few seconds)
(; T) = transportmatrix(; ϕ, mlotst, gridmetrics, indices, ρ, κH, κVML, κVdeep)

That's it! You've got yourself the transport matrix of your dreams!

Warning

This does not work for all CMIP models! See below for a list of models that pass the rudimentary tests.

List of models tested

ACCESS1-0
ACCESS1-3
ACCESS-ESM1-5
ACCESS-CM2

Tip

The test/ directory contains up-to-date examples for building the transport matrix and doing some simple calculations.

Citation

This code is © Benoît Pasquier (2024) and contributors, and it is made available under the MIT license enclosed with the software.

Over and above the legal restrictions imposed by this license, if you use this software for an academic publication then you are obliged to provide proper attribution. This can be to this code directly,

Benoît Pasquier (2024) “OceanTransportMatrixBuilder.jl: A Julia package to build ocean transport matrices from CMIP model output”. Zenodo. doi: 10.5281/zenodo.13864524.

or to the paper (currently in preparation) that describes it, or (ideally) both. You can also find the citation(s) in BibTeX format in the CITATION.bib file.

Acknowledgements

This package simply implements published work¹, and the Julia code available here was essentially translated from code by Matt Chamberlain at CSIRO. Benoît Pasquier (@briochemc) is funded through CSIRO's CarbonLock Future Science Platform, supervised by Richard Matear at CSIRO.

Footnotes

1. Chamberlain et al. (2019) Transport matrices from standard ocean-model output and quantifying circulation response to climate change. ↩ ↩²

2. Holzer et al. (2020) Climate-driven changes in the ocean's ventilation pathways and time scales diagnosed from transport matrices. ↩

3. Pasquier et al. (2023) Optimal parameters for the ocean's nutrient, carbon, and oxygen cycles compensate for circulation biases but replumb the biological pump. ↩

4. Pasquier et al. (2024a) The biological and preformed carbon pumps in perpetually slower and warmer oceans. ↩

5. Pasquier et al. (2024b) Deoxygenation and Its Drivers Analyzed in Steady State for Perpetually Slower and Warmer Oceans. ↩

6. DeVries et al. (2012) The sequestration efficiency of the biological pump. ↩

7. Siegel et al. (2021) Assessing the sequestration time scales of some ocean-based carbon dioxide reduction strategies. ↩

8. John et al. (2020) AWESOME OCIM: A simple, flexible, and powerful tool for modeling elemental cycling in the oceans. ↩