Public API

NumericalEarth.NumericalEarthModule

NumericalEarth.jl

🌎 A framework for simulating the Earth system at all scales with prescribed or prognostic atmosphere, ocean, sea ice and land components, based on Oceananigans.

Overview

NumericalEarth.jl provides infrastructure for running Earth system model components—ocean, atmosphere, sea ice, and others—coupled together or driven by prescribed datasets. The coupling interface is generic: plug in Oceananigans for ocean dynamics, ClimaSeaIce for sea ice, SpeedyWeather or other atmospheric models, or use reanalysis products like JRA55 and ERA5 as prescribed forcing.

The package handles the complexity of component communication: interpolating between grids, computing air-sea fluxes via similarity theory, managing radiative transfer, and synchronizing time-stepping across components with different temporal resolutions.

NumericalEarth.jl also serves as a sandbox for developing and testing interface parameterizations—bulk flux formulations, roughness length models, albedo schemes, and other boundary layer physics—in a modular environment where they can be validated against observations before deployment in production climate models.

Installation instructions

NumericalEarth is a registered Julia package. So to install it,

  1. Download Julia (version 1.10 or later).

  2. Launch Julia and type

julia> using Pkg

julia> Pkg.add("NumericalEarth")

This installs the latest version that's compatible with your current environment.

Data Wrangling

Running realistic Earth system simulations requires wrangling gigabytes of observational and reanalysis data into formats your model can ingest. NumericalEarth.jl abstracts away this pain. Point the package at a dataset and a target grid, and it handles the downloading, caching, and regridding automatically.

The Metadatum abstraction provides a unified interface: whether you're initializing ocean temperature from reanalysis or prescribing atmospheric boundary conditions, the workflow is the same.

using NumericalEarth
using Dates

# Load temperature from reanalysis on a specific date
T_init = Metadatum(:temperature; date=DateTime(1993, 1, 1), dataset=ECCO2Daily())

# Build a prescribed atmosphere from JRA55
atmosphere = JRA55PrescribedAtmosphere(arch)

A core abstraction: EarthSystemModel

The coupling infrastructure is anchored by EarthSystemModel, which encapsulates the component models—ocean, sea ice, atmosphere—and specifies how they communicate. Each component can be either prognostic (time-stepped by its own dynamics) or prescribed (interpolated from data). The model handles flux computations at interfaces, grid interpolation between components, and synchronized time-stepping.

We conceive of EarthSystemModel as a model in its own right, not just a container for components. This means it works with all the Oceananigans tools you'd use for any other model—run!(simulation), Callback, Checkpointer, output writers, and the rest.

To illustrate, here's a global ocean simulation driven by prescribed atmospheric reanalysis:

using Oceananigans
using Oceananigans.Units
using Dates
using CUDA
import NumericalEarth

arch = GPU()
grid = LatitudeLongitudeGrid(arch,
                             size = (1440, 560, 10),
                             halo = (7, 7, 7),
                             longitude = (0, 360),
                             latitude = (-70, 70),
                             z = (-3000, 0))

bathymetry = NumericalEarth.regrid_bathymetry(grid)
grid = ImmersedBoundaryGrid(grid, GridFittedBottom(bathymetry))

# Build an ocean simulation initialized to the ECCO state estimate
ocean = NumericalEarth.ocean_simulation(grid)
start_date = DateTime(1993, 1, 1)
set!(ocean.model,
     T=NumericalEarth.Metadatum(:temperature; date=start_date, dataset=NumericalEarth.ECCO2Daily()),
     S=NumericalEarth.Metadatum(:salinity;    date=start_date, dataset=NumericalEarth.ECCO2Daily()))

# Couple the ocean to JRA55 atmospheric forcing
atmosphere = NumericalEarth.JRA55PrescribedAtmosphere(arch)
coupled_model = NumericalEarth.OceanOnlyModel(ocean; atmosphere)
simulation = Simulation(coupled_model, Δt=20minutes, stop_time=30days)
run!(simulation)

This simulation achieves approximately 8 simulated years per day of wall time on an Nvidia H100 GPU.

Since ocean.model is an Oceananigans.HydrostaticFreeSurfaceModel, we can leverage Oceananigans features directly. For example, to plot the surface speed at the end of the simulation:

u, v, w = ocean.model.velocities
speed = Field(sqrt(u^2 + v^2))

using GLMakie
heatmap(view(speed, :, :, ocean.model.grid.Nz), colorrange=(0, 0.5), colormap=:magma, nan_color=:lightgray)

which produces

image

Installation

NumericalEarth is a registered Julia package. To install from a Julia REPL:

julia> using Pkg

julia> Pkg.add("NumericalEarth")

Use Pkg.add(url="https://github.com/NumericalEarth/NumericalEarth.jl.git", rev="main") to install the latest development version.

Citing

If you use NumericalEarth for your research, teaching, or fun 🤩, everyone in our community will be grateful if you give credit by citing the corresponding Zenodo record, e.g.,

Wagner, G. L. et al. (2026). NumericalEarth/NumericalEarth.jl: v0.1.1 (v0.1.1). Zenodo. https://doi.org/10.5281/zenodo.18665405

and also the recent preprint submitted to the Journal of Advances in Modeling Earth Systems that presents an overview of all the things that make Oceananigans unique:

"High-level, high-resolution ocean modeling at all scales with Oceananigans"

by Gregory L. Wagner, Simone Silvestri, Navid C. Constantinou, Ali Ramadhan, Jean-Michel Campin, Chris Hill, Tomas Chor, Jago Strong-Wright, Xin Kai Lee, Francis Poulin, Andre Souza, Keaton J. Burns, Siddhartha Bishnu, John Marshall, and Raffaele Ferrari

arXiv:2502.14148

bibtex

@article{Oceananigans-overview-paper-2025,
  title = {{High-level, high-resolution ocean modeling at all scales with Oceananigans}},
  author = {G. L. Wagner and S. Silvestri and N. C. Constantinou and A. Ramadhan and J.-M. Campin and C. Hill and T. Chor and J. Strong-Wright and X. K. Lee and F. Poulin and A. Souza and K. J. Burns and S. Bishnu and J. Marshall and R. Ferrari},
  journal = {arXiv preprint},
  year = {2025},
  archivePrefix = {arXiv},
  eprint = {2502.14148},
  doi = {10.48550/arXiv.2502.14148},
  notes = {submitted to the Journal of Advances in Modeling Earth Systems},
}
source

Atmospheres

NumericalEarth.Atmospheres.PrescribedAtmosphereType
PrescribedAtmosphere(grid, times=[zero(grid)];
                     clock = Clock{Float64}(time = 0),
                     surface_layer_height = 10, # meters
                     boundary_layer_height = 512, # meters
                     thermodynamics_parameters = AtmosphereThermodynamicsParameters(eltype(grid)),
                     velocities      = default_atmosphere_velocities(grid, times),
                     tracers         = default_atmosphere_tracers(grid, times),
                     pressure        = default_atmosphere_pressure(grid, times),
                     freshwater_flux = default_freshwater_flux(grid, times))

Return a prescribed, time-evolving atmospheric state with data on grid and at given times.

!!! compat Radiation component The downwelling shortwave / longwave radiation part of the top-level radiation component (see PrescribedRadiation, JRA55PrescribedRadiation).

source
NumericalEarth.Atmospheres.PrescribedPrecipitationFluxType
PrescribedPrecipitationFlux(; rain=nothing, snow=nothing)
PrescribedPrecipitationFlux(rain, snow)

Container for prescribed precipitation fluxes. Either component may be nothing to indicate that the corresponding precipitation type is not represented by the atmosphere (e.g. rain-only datasets). Used as the freshwater_flux of a PrescribedAtmosphere; downstream callers query the snow component via surface_snowfall_flux so that prognostic atmospheres with or without snow can dispatch on this type as well.

source

Bathymetry

NumericalEarth.Bathymetry.ORCAGridFunction
ORCAGrid(arch = CPU(), FT::DataType = Float64;
         dataset,
         halo = (4, 4, 4),
         z = (-6000, 0),
         Nz = 50,
         radius = Oceananigans.defaults.planet_radius,
         with_bathymetry = true,
         active_cells_map = true,
         major_basins = Inf,
         south_rows_to_remove = default_south_rows_to_remove(dataset),
         dir = default_download_directory(dataset))

Construct an OrthogonalSphericalShellGrid with (Periodic, RightFaceFolded, Bounded) topology using coordinate and metric data from a NEMO eORCA mesh_mask file.

The dataset keyword argument specifies which ORCA configuration to use (e.g., ORCA1() or ORCA12()). The mesh mask and bathymetry files are downloaded automatically via the DataWrangling.ORCA metadata interface.

The horizontal grid (including coordinates, scale factors, and areas) is loaded directly from the mesh_mask NetCDF file. If all staggered NEMO fields are present (T, U, V, F points), they are used directly. If only T and F coordinates are available (glamt/gphit/glamf/gphif), staggered coordinates and metrics are reconstructed approximately using Tripolar-style spherical assumptions.

When with_bathymetry = true (the default), the bathymetry is also downloaded and the grid is returned as an ImmersedBoundaryGrid with a GridFittedBottom.

Positional Arguments

  • arch: The architecture (e.g., CPU() or GPU()). Default: CPU().
  • FT: Floating point type. Default: Float64.

Keyword Arguments

  • dataset: The ORCA dataset to use. Default: ORCA1(). ORCA12() is also supported (ORCA1 data from Zenodo; https://doi.org/10.5281/zenodo.4436658).
  • halo: Halo size tuple (Hx, Hy, Hz). Default: (4, 4, 4).
  • z: Vertical coordinate specification. Can be a 2-tuple (z_bottom, z_top), an array of z-interfaces, or, e.g., an ExponentialDiscretization. Default: (-6000, 0).
  • Nz: Number of vertical levels (only used when z is a 2-tuple). Default: 50.
  • radius: Planet radius. Default: Oceananigans.defaults.planet_radius.
  • with_bathymetry: If true, download the bathymetry and return an ImmersedBoundaryGrid with GridFittedBottom. Default: true.
  • active_cells_map: If true and with_bathymetry = true, build an active cells map for efficient kernel execution over wet cells only. Default: true.
  • major_basins: Number of independent connected ocean basins to retain via remove_minor_basins!. Basins are removed from smallest to largest; major_basins = 1 keeps only the largest. Default: Inf (keep all basins).
  • south_rows_to_remove: Number of southern rows to remove from the eORCA grid. The "extended" eORCA grid contains degenerate padding rows near Antarctica that are entirely land. Removing them reduces memory usage and computation.
  • dir: Directory to store and look up ORCA files (mesh_mask and bathymetry). Defaults to the dataset scratch cache via default_download_directory(dataset).
source
NumericalEarth.Bathymetry.regrid_bathymetryMethod
regrid_bathymetry(target_grid, metadata;
                  height_above_water = nothing,
                  minimum_depth = 0,
                  major_basins = 1,
                  interpolation_passes = 1,
                  cache = true)

Return bathymetry that corresponds to metadata onto target_grid.

Arguments

  • target_grid: grid to interpolate the bathymetry onto.

Keyword Arguments

  • height_above_water: limits the maximum height of above-water topography (where $h > 0$) before interpolating. Default: nothing, which implies that the original topography is retained.

  • minimum_depth: minimum depth for the shallow regions, defined as a positive value. h > - minimum_depth is considered land. Default: 0.

  • interpolation_passes: regridding/interpolation passes. The bathymetry is interpolated in interpolation_passes - 1 intermediate steps. The more the interpolation steps the smoother the final bathymetry becomes.

    Example

    Interpolating from a 400 x 200 grid to a 100 x 100 grid in 4 passes involves:

    • 400 x 200 → 325 x 175
    • 325 x 175 → 250 x 150
    • 250 x 150 → 175 x 125
    • 175 x 125 → 100 x 100

    If coarsening the original grid, linear interpolation in passes is equivalent to applying a smoothing filter, with more passes increasing the strength of the filter. If refining the original grid, additional passes do not help and no intermediate steps are performed.

  • major_basins: Number of "independent major basins", or fluid regions fully encompassed by land, that are retained by remove_minor_basins!. Basins are removed by order of size: the smallest basins are removed first. major_basins = 1 retains only the largest basin. If Inf then no basins are removed. Default: 1.

  • cache: If true (default), caches the regridded bathymetry to disk and reuses it on subsequent calls with the same grid and parameters. Set to false to force recomputation.

source

DataWrangling

NumericalEarth.DataWrangling.ColumnType
Column(longitude, latitude; z=nothing, interpolation=Linear())

Create a column region at a single horizontal point (longitude, latitude). When used as a Metadata region, native_grid returns a single-column RectilinearGrid and location reduces horizontal dimensions to Nothing.

Keyword Arguments

  • z: depth range tuple (z_bottom, z_top) for restricting downloads (used by CopernicusMarine/GLORYS). Default: nothing (full depth).
  • interpolation: method for extracting data from the surrounding grid cells. Linear() (default) bilinearly interpolates to the exact point; Nearest() selects the closest grid cell.
source
NumericalEarth.DataWrangling.DatasetRestoringType
DatasetRestoring(metadata::Metadata,
                 arch_or_grid = CPU();
                 rate,
                 mask = 1,
                 time_indices_in_memory = 2, # Not more than this if we want to use GPU!
                 time_indexing = Cyclical(),
                 inpainting = NearestNeighborInpainting(Inf),
                 cache_inpainted_data = true)

Return a forcing term that restores to data from a dataset. The restoring is applied as a forcing on the right hand side of the evolution equations, calculated as:

\[F_ψ = r μ (ψ_{dataset} - ψ)\]

where $μ$ is the mask, $r$ is the restoring rate, $ψ$ is the simulation variable, and $ψ_{dataset}$ is the dataset variable that is linearly interpolated in space and time from the dataset of choice to the simulation grid and time.

Arguments

  • metadata: The medatada for a dataset variable to restore. Choices for variables include:

    • :temperature,
    • :salinity,
    • :u_velocity,
    • :v_velocity,
    • :sea_ice_thickness,
    • :sea_ice_area_fraction,
    • :dissolved_inorganic_carbon,
    • :alkalinity,
    • :nitrate,
    • :phosphate,
    • :dissolved_organic_phosphorus,
    • :particulate_organic_phosphorus,
    • :dissolved_iron,
    • :dissolved_silicate,
    • :dissolved_oxygen.
  • arch_or_grid: Either the architecture of the simulation, or a grid on which the data is pre-interpolated when loaded. If an architecture is provided, such as arch_or_grid = CPU() or arch_or_grid = GPU(), data is interpolated on-the-fly when the forcing tendency is computed. Default: CPU().

Keyword Arguments

  • rate: The restoring rate, i.e., the inverse of the restoring timescale (in s⁻¹).

  • mask: The mask value. Can be a function of (x, y, z, time), an array, or a number.

  • time_indices_in_memory: The number of time indices to keep in memory. The number is chosen based on a trade-off between increased performance (more indices in memory) and reduced memory footprint (fewer indices in memory). Default: 2.

  • time_indexing: The time indexing scheme for the field time series. Default: Cyclical().

  • inpainting: inpainting algorithm, see inpaint_mask!. Default: NearestNeighborInpainting(Inf).

  • cache_inpainted_data: If true, the data is cached to disk after inpainting for later retrieving. Default: true.

source
NumericalEarth.DataWrangling.LinearlyTaperedPolarMaskMethod
LinearlyTaperedPolarMask(; northern = (70,   75),
                           southern = (-75, -70),
                           z = (-20, 0))

Build a mask that is linearly tapered in latitude between the northern and southern edges. The mask is constant in depth between the z and equals zero everywhere else. The mask is limited to lie between (0, 1). The mask has the following functional form:

n = 1 / (northern[2] - northern[1]) * (φ - northern[1])
s = 1 / (southern[1] - southern[2]) * (φ - southern[2])

valid_depth = (z[1] < z < z[2])

mask = valid_depth ? clamp(max(n, s), 0, 1) : 0
source
NumericalEarth.DataWrangling.MetadataMethod
Metadata(variable_name;
         dataset,
         dates = all_dates(dataset, variable_name),
         dir = default_download_directory(dataset),
         region = nothing,
         filename = nothing,
         start_date = nothing,
         end_date = nothing)

Metadata holding a specific dataset information.

Argument

  • variable_name: a symbol representing the name of the variable (for example, :temperature, :salinity, :u_velocity, etc)

Keyword Arguments

  • dataset: Supported datasets are ETOPO2022(), ECCO2Monthly(), ECCO2Daily(), ECCO4Monthly(), EN4Monthly(), GLORYSDaily(), GLORYSMonthly(), RepeatYearJRA55(), and MultiYearJRA55().

  • dates: The dates of the dataset (Dates.AbstractDateTime or CFTime.AbstractCFDateTime). Note that dates can either be a range or a vector of dates, representing a time-series. For a single date, use Metadatum.

  • start_date: If dates = nothing, we can prescribe the first date of metadata as a date (Dates.AbstractDateTime or CFTime.AbstractCFDateTime). If outside the date range of the dataset, the first allowable date is chosen. Default: nothing.

  • end_date: If dates = nothing, we can prescribe the last date of metadata as a date (Dates.AbstractDateTime or CFTime.AbstractCFDateTime). If outside the date range of the dataset, the last allowable date is chosen. Default: nothing.

  • region: Specifies the spatial region of the dataset. Can be a BoundingBox for a rectangular region, a Column for a single horizontal location, or nothing for the full domain.

  • filename: The filename(s) for the dataset. If nothing, the filename is computed from the dataset type. Can be a String (single file for all dates) or a DatewiseFilename (one file per date).

  • dir: The directory where the dataset is stored.

source
NumericalEarth.DataWrangling.MetadataSetMethod
MetadataSet(variable_names::Symbol...;
            dataset,
            dates = all_dates(dataset, first(variable_names)),
            date = nothing,
            dir = default_download_directory(dataset),
            region = nothing,
            filenames = nothing,
            start_date = nothing,
            end_date = nothing)

A bundle of Metadata for many variables that share dataset, dates, region, and dir — differing only in variable name.

Each element of an mset::MetadataSet, e.g., mset[name] (or equivalently mset.name or mset[i]) is itself a Metadata — or a Metadatum when dates is a single date. Iteration walks the variable axis, yielding one Metadata per variable.

Arguments

  • variable_names: one or more Symbols naming the dataset variables to bundle (e.g. :temperature, :salinity). Verbose dataset-internal names — no aliases. It can also be a tuple of Symbols, e.g., (:temperature, :salinity), but not a vector of Symbols.

Keyword Arguments

  • dataset: the shared dataset (e.g. ECCO4Monthly(), ERA5HourlyPressureLevels()).
  • dates: shared date axis. Either a single AbstractDateTime/AbstractCFDateTime or an AbstractVector of dates. Defaults to all_dates(dataset, first(variable_names)).
  • date: convenience scalar form; cannot be used together with dates.
  • region: shared spatial region — BoundingBox, Column, or nothing.
  • dir: shared download directory.
  • filenames: an optional NamedTuple keyed by variable_names overriding the auto-computed per-variable filenames.
  • start_date, end_date: optional date cropping, matching Metadata.

Example

using NumericalEarth, Dates

mset = MetadataSet(:temperature, :salinity;
                   dataset = ECCO4Monthly(),
                   date = DateTime(1995, 1, 1))

mset[2] # Metadata for :salinity

using NumericalEarth, Dates

mset = MetadataSet(:temperature, :salinity;
                   dataset = ECCO4Monthly(),
                   date = DateTime(1995, 1, 1))

mset[2] # Metadata for :salinity

# output
Metadatum{ECCO4Monthly, DateTime}:
├── name: salinity
├── dataset: ECCO4Monthly
├── dates: 1995-01-01 00:00:00
├── filename: SALT_1995_01.nc
└── dir: /.julia/scratchspaces/904d977b-046a-4731-8b86-9235c0d1ef02/ECCO/v4

See also Metadata, Metadatum.

source
NumericalEarth.DataWrangling.MetadatumMethod
Metadatum(variable_name;
          dataset,
          region = nothing,
          date = first_date(dataset, variable_name),
          filename = nothing,
          dir = default_download_directory(dataset))

A specialized constructor for a Metadata object with a single date, representative of a snapshot in time.

source
NumericalEarth.DataWrangling.SurfaceFluxRestoringType
SurfaceFluxRestoring(dataset_restoring::DatasetRestoring)

A thin wrapper around a DatasetRestoring that converts a 3D restoring tendency into a 2D surface flux boundary condition.

When used as a boundary condition (via getbc), the wrapped DatasetRestoring is evaluated at the top cell (k = Nz) and the resulting tendency G is converted to a surface flux as -G * Δz, consistent with the Oceananigans top-flux sign convention (tendency contribution = -J / Δz).

This is intended for use with the additional_surface_fluxes keyword argument of ocean_simulation, allowing a DatasetRestoring to contribute an additional flux at the surface without overwriting the coupled exchange fluxes.

Example

using ..NumericalEarth

restoring = DatasetRestoring(metadata, grid; rate = 1 / 30days)
ocean = ocean_simulation(grid;
    additional_surface_fluxes = (; S = SurfaceFluxRestoring(restoring)))
source
NumericalEarth.DataWrangling.native_gridFunction
native_grid(metadata::Metadata, arch=CPU(); halo = (3, 3, 3))

Return the native grid corresponding to metadata with halo size. Returns a LatitudeLongitudeGrid for global or BoundingBox regions, and a column RectilinearGrid for Column regions.

source
NumericalEarth.DataWrangling.validate_dataset_coverageMethod
validate_dataset_coverage(grid, metadata)

Check that grid lies within the spatial coverage of metadata's dataset. Throws an error if the grid extends outside the dataset's domain. The default implementation does nothing (all grids are accepted). Dataset-specific methods can override this to enforce coverage constraints.

source

DataWrangling.ECCO

NumericalEarth.DataWrangling.ECCO.ECCOPrescribedAtmosphereFunction
ECCOPrescribedAtmosphere([architecture = CPU()];
                          dataset = ECCO4Monthly(),
                          start_date = first_date(dataset, :air_temperature),
                          end_date = last_date(dataset, :air_temperature),
                          dir = default_download_directory(dataset),
                          time_indices_in_memory = 10,
                          time_indexing = Cyclical(),
                          surface_layer_height = 2,  # meters
                          other_kw...)

Return a PrescribedAtmosphere representing ECCO state estimate data. The atmospheric data will be held in FieldTimeSeries objects containing

  • velocities: u, v
  • air temperature and humidity: T, q
  • surface pressure: p
  • freshwater flux: rain

Note: downwelling shortwave / longwave radiation is now part of the top-level radiation component. Use ECCOPrescribedRadiation to load ECCO SW/LW into a PrescribedRadiation.

source
NumericalEarth.DataWrangling.ECCO.ECCOPrescribedRadiationFunction
ECCOPrescribedRadiation([architecture = CPU()];
                        dataset = ECCO4Monthly(),
                        start_date = first_date(dataset, :downwelling_shortwave),
                        end_date = last_date(dataset, :downwelling_shortwave),
                        dir = default_download_directory(dataset),
                        time_indices_in_memory = 10,
                        time_indexing = Cyclical(),
                        ocean_surface = SurfaceRadiationProperties(0.05, 0.97),
                        sea_ice_surface = SurfaceRadiationProperties(0.7, 1.0),
                        stefan_boltzmann_constant = default_stefan_boltzmann_constant,
                        other_kw...)

Return a PrescribedRadiation backed by ECCO downwelling shortwave and longwave fields.

source

DataWrangling.EN4

DataWrangling.ERA5

NumericalEarth.DataWrangling.ERA5.ERA5HourlyPressureLevelsMethod
ERA5HourlyPressureLevels(; pressure_levels = ERA5_all_pressure_levels, z = nothing)

Hourly ERA5 pressure-levels dataset metadata. By default (z = nothing) the native grid's vertical coordinate is a 3-D PressureLevelVerticalDiscretization built from the instantaneous geopotential Φ(λ,φ,p)/g with sub-surface levels clipped at the surface — the right thing over terrain and synoptic Φ-displacement (issue #236).

To pin the z-coordinate to something static — e.g. for a quick test without an extra Φ download — pass a precomputed vector: z = standard_atmosphere_z_interfaces(pressure_levels).

source
NumericalEarth.DataWrangling.ERA5.pressure_fieldFunction
pressure_field(metadata::ERA5PressureMetadatum, arch=CPU(); halo=(3,3,3))

Return a Field{Nothing, Nothing, Center} on the native grid of metadata holding the pressure value (Pa) at each vertical level. Levels are ordered bottom-to-top (k=1 is the highest pressure level). The Nothing horizontal locations make this field broadcast against full 3-D fields without copying.

source

DataWrangling.ETOPO

DataWrangling.GEBCO

NumericalEarth.DataWrangling.GEBCO.GEBCO2024Type
GEBCO2024

General Bathymetric Chart of the Oceans 2024 release. Global bathymetry and topography at 15 arc-second resolution.

The GEBCO_2024 Grid is a global terrain model for ocean and land, providing elevation data on a 15 arc-second interval grid.

Reference: GEBCO Compilation Group (2024) GEBCO 2024 Grid Data source: https://www.gebco.net/dataandproducts/griddedbathymetrydata/

source

DataWrangling.GLORYS

DataWrangling.IBCAO

NumericalEarth.DataWrangling.IBCAO.IBCAOv5Type
IBCAOv5

International Bathymetric Chart of the Arctic Ocean Version 5.1 (2025). 100m resolution bathymetry for the Arctic Ocean (north of 64°N), including Greenland ice sheet surface elevation. Data is provided in Polar Stereographic projection (EPSG:3996) and reprojected to WGS84 geographic coordinates at 0.01° resolution on first use.

Reference: Jakobsson et al. (2024), https://doi.org/10.1038/s41597-024-04278-w Data source: https://www.gebco.net/data-products/gridded-bathymetry-data/arctic-ocean/

source

DataWrangling.IBCSO

NumericalEarth.DataWrangling.IBCSO.IBCSOv2Type
IBCSOv2

International Bathymetric Chart of the Southern Ocean Version 2 (2024 Annual Release). High-resolution (500m) bathymetry for the Southern Ocean (south of 50°S).

Reference: Dorschel et al. (2022), https://doi.org/10.1594/PANGAEA.937574 Data source: https://ibcso.org/ibcso-2024-annual-release/

source

DataWrangling.JRA55

NumericalEarth.DataWrangling.JRA55.JRA55PrescribedAtmosphereFunction
JRA55PrescribedAtmosphere([architecture = CPU()];
                          dataset = RepeatYearJRA55(),
                          start_date = first_date(dataset, :temperature),
                          end_date = last_date(dataset, :temperature),
                          dir = download_JRA55_cache,
                          time_indices_in_memory = 10,
                          time_indexing = Cyclical(),
                          surface_layer_height = 10,  # meters
                          region = nothing,
                          other_kw...)

Return a PrescribedAtmosphere representing JRA55 reanalysis data. Each atmospheric field is constructed via FieldTimeSeries(::JRA55Metadata), which uses a DatasetBackend parameterised by JRA55 metadata so that the JRA55-specific set! (chunked-yearly NetCDF) is dispatched. The region keyword restricts the atmosphere to a sub-domain of the global JRA55 grid.

Note: downwelling shortwave / longwave radiation is now part of the top-level radiation component. Use JRA55PrescribedRadiation to load JRA55 SW/LW into a PrescribedRadiation.

source
NumericalEarth.DataWrangling.JRA55.JRA55PrescribedLandFunction
JRA55PrescribedLand([architecture = CPU()];
                    dataset = RepeatYearJRA55(),
                    start_date = first_date(dataset, :river_freshwater_flux),
                    end_date = last_date(dataset, :river_freshwater_flux),
                    dir = download_JRA55_cache,
                    time_indices_in_memory = 10,
                    time_indexing = Cyclical(),
                    region = nothing,
                    other_kw...)

Return a PrescribedLand representing JRA55 reanalysis land surface data (river runoff and iceberg calving freshwater fluxes).

source
NumericalEarth.DataWrangling.JRA55.JRA55PrescribedRadiationFunction
JRA55PrescribedRadiation([architecture = CPU()];
                         dataset = RepeatYearJRA55(),
                         start_date = first_date(dataset, :downwelling_shortwave_radiation),
                         end_date = last_date(dataset, :downwelling_shortwave_radiation),
                         dir = download_JRA55_cache,
                         time_indices_in_memory = 10,
                         time_indexing = Cyclical(),
                         ocean_surface = SurfaceRadiationProperties(0.05, 0.97),
                         sea_ice_surface = SurfaceRadiationProperties(0.7, 1.0),
                         stefan_boltzmann_constant = default_stefan_boltzmann_constant,
                         region = nothing,
                         other_kw...)

Return a PrescribedRadiation backed by JRA55 downwelling shortwave and longwave FieldTimeSeries. Surface radiative properties (albedo, emissivity) for ocean and sea-ice surfaces default to standard values; pass ocean_surface = nothing (or sea_ice_surface = nothing) to omit a surface.

source

DataWrangling.ORCA

DataWrangling.OSPapa

NumericalEarth.DataWrangling.OSPapa.OSPapaPrescribedAtmosphereFunction
OSPapaPrescribedAtmosphere(architecture = CPU(), FT = Float32;
                            start_date = first_date(OSPapaHourly(), :air_temperature),
                            end_date   = last_date(OSPapaHourly(), :air_temperature),
                            dir = download_OSPapa_cache,
                            surface_layer_height = 2.5,
                            max_gap_hours = 72)

Construct a PrescribedAtmosphere from Ocean Station Papa buoy observations.

Data is automatically downloaded from the NOAA/PMEL AWS S3 bucket if not already cached locally.

Radiation and albedo

The buoy SW and LW variables are downwelling fluxes. When this atmosphere is used with OceanOnlyModel, ClimaOcean applies its own ocean albedo (default α = 0.05) to compute net absorbed shortwave, and computes upwelling longwave from the model SST via Stefan-Boltzmann. This means the resulting net heat flux will differ from the COARE-computed QNET available via os_papa_prescribed_fluxes. If you need the exact observed net fluxes, use os_papa_prescribed_flux_boundary_conditions instead.

Keyword Arguments

  • start_date: start of the time range
  • end_date: end of the time range
  • dir: directory for cached data files
  • surface_layer_height: measurement height in meters (default: 2.5, matching the buoy's temperature/humidity instruments)
  • max_gap_hours: maximum gap size (in hours) to fill by linear interpolation (default: 72)
source
NumericalEarth.DataWrangling.OSPapa.OSPapaPrescribedRadiationFunction
OSPapaPrescribedRadiation(architecture = CPU(), FT = Float32;
                          start_date = first_date(OSPapaHourly(), :shortwave_radiation),
                          end_date   = last_date(OSPapaHourly(), :shortwave_radiation),
                          dir = download_OSPapa_cache,
                          max_gap_hours = 72,
                          ocean_surface = SurfaceRadiationProperties(0.05, 0.97),
                          sea_ice_surface = nothing,
                          stefan_boltzmann_constant = default_stefan_boltzmann_constant)

Construct a PrescribedRadiation from Ocean Station Papa buoy SW/LW observations on a single-column grid.

source
NumericalEarth.DataWrangling.OSPapa.os_papa_prescribed_flux_boundary_conditionsMethod
os_papa_prescribed_flux_boundary_conditions(fluxes;
                                            arch=nothing,
                                            ρ₀=1020.0,
                                            cₚ=3991.0,
                                            u_momentum_flux_correction=no_correction,
                                            v_momentum_flux_correction=no_correction,
                                            temperature_flux_correction=no_correction,
                                            salinity_flux_correction=no_correction)

Create Oceananigans FluxBoundaryConditions for u, v, T, S from prescribed OS Papa flux data. Returns a NamedTuple of FieldBoundaryConditions that can be passed directly to ocean_simulation or HydrostaticFreeSurfaceModel via the boundary_conditions keyword argument.

If arch is provided, the flux FieldTimeSeries are first moved to that architecture.

Arguments

Keyword Arguments

  • arch: target architecture for the flux FieldTimeSeries. If provided, each flux FTS is moved to arch via on_architecture. Defaults to nothing, which keeps them on whatever architecture they were built on.
  • ρ₀: reference ocean density (default: 1020 kg/m³)
  • cₚ: ocean heat capacity (default: 3991 J/(kg·K))
  • u_momentum_flux_correction: discrete-form correction function added to the zonal momentum flux boundary condition (default: no_correction)
  • v_momentum_flux_correction: discrete-form correction function added to the meridional momentum flux boundary condition (default: no_correction)
  • temperature_flux_correction: discrete-form correction function added to the temperature flux boundary condition (default: no_correction)
  • salinity_flux_correction: discrete-form correction function added after converting freshwater flux from EMP into the salinity flux boundary condition (default: no_correction)

Each correction function must have the signature (i, j, grid, clock, model_fields, p) and return a value in the same units as the corresponding boundary condition.

Examples

# Basic usage on GPU:
fluxes = os_papa_prescribed_fluxes(GPU(); start_date, end_date)
bcs = os_papa_prescribed_flux_boundary_conditions(fluxes)
ocean = ocean_simulation(grid; Δt=10minutes, boundary_conditions=bcs)

# With a uniform heat flux correction of +5 W/m² to close the heat budget:
heat_correction = (i, j, grid, clock, model_fields, p) -> 5.0 / (p.ρ₀ * p.cₚ)
bcs = os_papa_prescribed_flux_boundary_conditions(fluxes; temperature_flux_correction=heat_correction)
source
NumericalEarth.DataWrangling.OSPapa.os_papa_prescribed_fluxesFunction
os_papa_prescribed_fluxes(architecture = CPU(), FT = Float64;
                          start_date = first_date(OSPapaFluxHourly(), :net_heat_flux),
                          end_date   = last_date(OSPapaFluxHourly(), :net_heat_flux),
                          dir = download_OSPapa_cache,
                          max_gap_hours = 72)

Download precomputed air-sea fluxes for Ocean Station Papa from the PMEL ERDDAP ocs_papa_flux dataset (computed with COARE 3.0b) and return them as a NamedTuple of FieldTimeSeries on a Flat, Flat, Flat grid at architecture.

The returned tuple contains:

  • Qnet: net heat flux (W/m², positive into ocean)
  • Qlat: latent heat flux (W/m², positive upward = ocean cooling)
  • Qsen: sensible heat flux (W/m², positive upward = ocean cooling)
  • SWnet: net shortwave radiation (W/m², positive downward)
  • LWnet: net longwave radiation (W/m², positive upward = ocean cooling)
  • τx, τy: zonal and meridional wind stress (N/m²)
  • evap, rain, EMP: evaporation, precipitation, E-P (mm/hr)
  • Tsk: skin temperature (°C)

Keyword Arguments

  • start_date: start of the time range
  • end_date: end of the time range
  • dir: directory for cached data files
  • max_gap_hours: maximum gap size (in hours) to fill by linear interpolation (default: 72)
source

DataWrangling.WOA

Diagnostics

NumericalEarth.Diagnostics.MixedLayerDepthFieldMethod
MixedLayerDepthField(bm, grid, tracers;
                     difference_criterion = 2.87e-4,
                     reference_depth = 10.0)

Mixed-layer-depth diagnostic using the de Boyer Montégut et al. (2004) DR003 buoyancy-difference criterion. The default difference_criterion = 2.87e-4 m s⁻² corresponds to Δσ = 0.03 kg m⁻³ at ρ₀ = 1025 kg m⁻³. The default reference_depth = 10 m is the standard dBM reference; buoyancy at that depth is obtained by linear interpolation between adjacent cell centers, so the diagnostic is insensitive to the vertical resolution near the surface.

References

  • de Boyer Montégut, C., G.Madec, A. S.Fischer, A.Lazar, and D.Iudicone (2004), Mixed layer depth over the global ocean: An examination of profile data and a profile-based climatology, J. Geophys. Res., 109, C12003, doi:10.1029/2004JC002378.
source
NumericalEarth.Diagnostics.meridional_heat_transportMethod
meridional_heat_transport(esm::EarthSystemModel;
                          reference_temperature = 0)

Return the meridional heat transport for the coupled esm::EarthSystemModel by computing the meridional heat flux.

The meridional heat transport is computed via:

\[\mathrm{MHT} ≡ ρᵒᶜ cᵒᶜ ∫ v (T - T_{\rm ref}) \, \mathrm{d}x \, \mathrm{d}z\]

Above, $T_{\rm ref}$ is a reference temperature and $ρᵒᶜ$ and $cᵒᶜ$ are the ocean reference density and specific heat capacity respectively.

Only works on LatitudeLongitudeGrid

The meridional_heat_transport diagnostic currently is only supported only on LongitudeLatitudeGrids.

Arguments

  • esm: An EarthSystemModel.

Keyword Arguments

  • reference_temperature: The reference temperature (in ᵒC) used for the calculation; default: 0 ᵒC.

    Reference temperature

    The reference temperature is only relevant when we compute the meridional heat transport over a section where there is a net volume transport. If we are computing the diagnostic globally, i.e., around a whole latitude circle, then by necessity there is no net volume transport and thus the reference temperature value is irrelevant. Section-averaged transport could also be considered as a reference temperature to remove residual barotropic volume fluxes in basin-scale/regional analyses where a net volume transport is present.

Example

using ..NumericalEarth
using Oceananigans

grid = RectilinearGrid(size = (4, 5, 2), extent = (1, 1, 1),
                       topology = (Periodic, Bounded, Bounded))

ocean = ocean_simulation(grid;
                         momentum_advection = nothing,
                         tracer_advection = nothing,
                         closure = nothing,
                         coriolis = nothing)

sea_ice = sea_ice_simulation(grid, ocean)

atmosphere = PrescribedAtmosphere(grid, [0.0])
radiation = PrescribedRadiation(grid)

esm = OceanSeaIceModel(ocean, sea_ice; atmosphere, radiation)

mht = meridional_heat_transport(esm)

# output

Integral of BinaryOperation at (Center, Face, Center) over dims (1, 3)
└── operand: BinaryOperation at (Center, Face, Center)
    └── grid: 4×5×2 RectilinearGrid{Float64, Periodic, Bounded, Bounded} on CPU with 3×3×2 halo
source

EarthSystemModels

NumericalEarth.EarthSystemModels.EarthSystemModelMethod
EarthSystemModel(radiation, atmosphere, land, sea_ice, ocean;
                 clock = Clock{Float64}(time=0),
                 ocean_reference_density = reference_density(ocean),
                 ocean_heat_capacity = heat_capacity(ocean),
                 sea_ice_reference_density = reference_density(sea_ice),
                 sea_ice_heat_capacity = heat_capacity(sea_ice),
                 interfaces = nothing)

Construct a coupled earth system model. Components are passed in struct order (top to bottom): radiation, atmosphere, land, sea_ice, ocean. Pass nothing for components that are absent. For simpler configurations, see OceanOnlyModel, OceanSeaIceModel, and AtmosphereOceanModel.

Arguments

  • radiation: Radiation component, or nothing for a radiatively decoupled surface. Pass a PrescribedRadiation (e.g. JRA55PrescribedRadiation(...)) to enable radiative forcing.
  • atmosphere: A representation of a possibly time-dependent atmospheric state, or nothing. For a prognostic atmosphere, use atmosphere_simulation. For prescribed atmospheric forcing, use JRA55PrescribedAtmosphere or PrescribedAtmosphere.
  • land: Land component, or nothing.
  • sea_ice: A representation of a possibly time-dependent sea ice state, or nothing. For example, the minimalist FreezingLimitedOceanTemperature represents oceanic latent heating during freezing only, but does not evolve sea ice variables. For prognostic sea ice use an Oceananigans.Simulation of ClimaSeaIce.SeaIceModel.
  • ocean: A representation of a possibly time-dependent ocean state. Currently, only Oceananigans.Simulations of Oceananigans.HydrostaticFreeSurfaceModel are tested.

Keyword Arguments

  • clock: Keeps track of time.
  • ocean_reference_density: Reference density for the ocean. Defaults to value from ocean model.
  • ocean_heat_capacity: Heat capacity for the ocean. Defaults to value from ocean model.
  • sea_ice_reference_density: Reference density for sea ice. Defaults to value from sea ice model.
  • sea_ice_heat_capacity: Heat capacity for sea ice. Defaults to value from sea ice model.
  • interfaces: Component interfaces for coupling. Defaults to nothing and will be constructed automatically. To customize the sea ice-ocean heat flux formulation, create interfaces manually using ComponentInterfaces.
source
NumericalEarth.EarthSystemModels.EarthSystemModelMethod
EarthSystemModel(; radiation = nothing,
                   atmosphere = nothing,
                   land = nothing,
                   sea_ice = nothing,
                   ocean = nothing,
                   kw...)

Keyword-only constructor for EarthSystemModel. Equivalent to the positional form, but lets you pass only the components you actually have:

EarthSystemModel(; atmosphere, ocean)                       # ocean + atmosphere
EarthSystemModel(; atmosphere, sea_ice, ocean, radiation)   # full coupled

All keyword arguments accepted by the positional constructor are forwarded.

source
NumericalEarth.EarthSystemModels.OceanOnlyModelMethod
OceanOnlyModel(ocean; atmosphere=nothing, radiation=nothing, land=nothing, kw...)

Construct an ocean-only model without a sea ice component. This is a convenience constructor for EarthSystemModel that sets sea_ice to FreezingLimitedOceanTemperature (a simple freezing limiter that does not evolve sea ice variables).

The atmosphere, radiation, and land keywords can be used to specify prescribed components (e.g., JRA55PrescribedAtmosphere). All other keyword arguments are forwarded to EarthSystemModel.

using NumericalEarth, Oceananigans

grid = LatitudeLongitudeGrid(size = (20, 20, 4),
                             z = (-100, 0),
                             latitude = (-80, 80),
                             longitude = (0, 360),
                             halo = (6, 6, 3))

ocean = ocean_simulation(grid, closure=nothing)
set!(ocean.model, T=20, S=35, u=0.01, v=-0.005)

ocean = OceanOnlyModel(ocean)
# output

EarthSystemModel{CPU}(time = 0 seconds, iteration = 0)
├── radiation: Nothing
├── atmosphere: Nothing
├── land: Nothing
├── sea_ice: FreezingLimitedOceanTemperature{ClimaSeaIce.SeaIceThermodynamics.LinearLiquidus{Float64}}
├── ocean: HydrostaticFreeSurfaceModel{CPU, LatitudeLongitudeGrid}(time = 0 seconds, iteration = 0)
└── interfaces: ComponentInterfaces
source
NumericalEarth.EarthSystemModels.OceanSeaIceModelMethod
OceanSeaIceModel(ocean, sea_ice; atmosphere=nothing, radiation=nothing, land=nothing, kw...)

Construct a coupled ocean–sea ice model.

This is a convenience constructor for EarthSystemModel with explicit ocean and sea ice components and optional prescribed atmosphere, prescribed radiation, and prescribed land.

Positional arguments are ocean then sea_ice.

Example

using NumericalEarth, Oceananigans

grid = LatitudeLongitudeGrid(size = (20, 20, 4),
                             z = (-100, 0),
                             latitude = (-80, 80),
                             longitude = (0, 360),
                             halo = (6, 6, 3))

ocean = ocean_simulation(grid, closure=nothing)
set!(ocean.model, T=20, S=35, u=0.01, v=-0.005)

sea_ice = sea_ice_simulation(grid, ocean)

hi(λ, φ) = φ > 70 || φ < -70
set!(sea_ice.model, h=hi, ℵ=hi)

coupled_model = OceanSeaIceModel(ocean, sea_ice)

# output

EarthSystemModel{CPU}(time = 0 seconds, iteration = 0)
├── radiation: Nothing
├── atmosphere: Nothing
├── land: Nothing
├── sea_ice: SeaIceModel
├── ocean: HydrostaticFreeSurfaceModel{CPU, LatitudeLongitudeGrid}(time = 0 seconds, iteration = 0)
└── interfaces: ComponentInterfaces
source

EarthSystemModels.InterfaceComputations

NumericalEarth.EarthSystemModels.InterfaceComputations.CoefficientBasedFluxesType
CoefficientBasedFluxes(FT = Oceananigans.defaults.FloatType;
                       transfer_coefficients = (1e-3, 1e-3, 1e-3),
                       solver_stop_criteria = nothing,
                       solver_tolerance = 1e-8,
                       solver_maxiter = 20)

Return the structure for computing turbulent fluxes using bulk transfer coefficients. Used in bulk flux calculations to determine the exchange of momentum, heat, and moisture between the ocean/ice surface and the atmosphere.

Arguments

  • FT: (optional) Float type for the coefficients, defaults to Oceananigans.defaults.FloatType

Keyword Arguments

  • transfer_coefficients: Transfer coefficients for momentum, heat, and moisture. Can be a SimilarityScales, a Tuple, or a NamedTuple with constant or callable entries, or an LargeYeagerTransferCoefficients. Defaults to (1e-3, 1e-3, 1e-3).
  • solver_stop_criteria: Criteria for iterative solver convergence. If nothing, creates new criteria using solver_tolerance and solver_maxiter.
  • solver_tolerance: Tolerance for solver convergence when creating new stop criteria, defaults to 1e-8.
  • solver_maxiter: Maximum iterations for solver when creating new stop criteria, defaults to 20

Example

using Oceananigans
using NumericalEarth

grid = RectilinearGrid(size=3, z=(-1, 0), topology=(Flat, Flat, Bounded))
ocean = ocean_simulation(grid; timestepper = :QuasiAdamsBashforth2)

ao_fluxes = CoefficientBasedFluxes(transfer_coefficients = (1e-2, 1e-3, 1e-3))

interfaces = ComponentInterfaces(nothing, ocean; atmosphere_ocean_fluxes=ao_fluxes)

# output
ComponentInterfaces
source
NumericalEarth.EarthSystemModels.InterfaceComputations.ComponentInterfacesType
ComponentInterfaces(atmosphere, ocean, sea_ice=nothing; kwargs...)

Construct component interfaces for atmosphere-ocean-sea ice coupling.

Keyword Arguments

  • sea_ice_ocean_heat_flux: formulation for sea ice-ocean heat flux. Options are:

    • IceBathHeatFlux(): bulk heat flux with interface at freezing point
    • ThreeEquationHeatFlux(): coupled heat/salt/freezing point system (default)
  • radiation: radiation component. Default: Radiation().

  • radiation: radiation component. Default: nothing.

  • freshwater_density: reference density of freshwater. Default: default_freshwater_density.

  • atmosphere_ocean_fluxes: flux formulation for atmosphere-ocean interface. Default: SimilarityTheoryFluxes().

  • atmosphere_sea_ice_fluxes: flux formulation for atmosphere-sea ice interface. Default: SimilarityTheoryFluxes().

  • atmosphere_ocean_interface_temperature: temperature formulation for atmosphere-ocean interface. Default: BulkTemperature().

  • atmosphere_ocean_interface_specific_humidity: specific humidity formulation. Default: default_ao_specific_humidity(ocean).

  • atmosphere_sea_ice_interface_temperature: temperature formulation for atmosphere-sea ice interface. Default: default_ai_temperature(sea_ice).

  • ocean_reference_density: reference density for the ocean. Default: reference_density(ocean).

  • ocean_heat_capacity: heat capacity for the ocean. Default: heat_capacity(ocean).

  • ocean_temperature_units: temperature units for the ocean. Default: DegreesCelsius().

  • sea_ice_temperature_units: temperature units for sea ice. Default: DegreesCelsius().

  • sea_ice_reference_density: reference density for sea ice. Default: reference_density(sea_ice).

  • sea_ice_heat_capacity: heat capacity for sea ice. Default: heat_capacity(sea_ice).

  • gravitational_acceleration: gravitational acceleration. Default: default_gravitational_acceleration.

source
NumericalEarth.EarthSystemModels.InterfaceComputations.IceBathHeatFluxType
IceBathHeatFlux(FT::DataType = Oceananigans.defaults.FloatType;
                heat_transfer_coefficient = 0.006,
                friction_velocity = 0.02)

Construct an IceBathHeatFlux with the specified parameters.

Keyword Arguments

  • heat_transfer_coefficient: turbulent heat exchange coefficient. Default: 0.006.
  • friction_velocity: friction velocity value or formulation. Default: 0.02.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.IceBathHeatFluxType
IceBathHeatFlux{FT, U}

Bulk formulation for sea ice-ocean heat flux.

The interface temperature is fixed at the freezing point of the surface salinity, and the heat flux is computed using bulk transfer:

\[Q = \rho_o c_o \alpha_h u_* (T - T_m)\]

where $\alpha_h$ is the heat transfer coefficient and $u_*$ is the friction velocity.

Fields

  • heat_transfer_coefficient::FT: turbulent heat exchange coefficient $\alpha_h$ (dimensionless)
  • friction_velocity::U: friction velocity value or formulation (constant Number or MomentumBasedFrictionVelocity)

Example

using NumericalEarth.EarthSystemModels: IceBathHeatFlux

flux = IceBathHeatFlux(heat_transfer_coefficient = 0.006, friction_velocity = 0.002)

# output
IceBathHeatFlux{Float64}
├── heat_transfer_coefficient: 0.006
└── friction_velocity: 0.002

References

  • Holland and Jenkins (1999): Holland, D. M., & Jenkins, A. (1999). Modeling thermodynamic ice–ocean interactions at the base of an ice shelf. Journal of Physical Oceanography, 29(8), 1787-1800.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.LargeYeagerTransferCoefficientsType
LargeYeagerTransferCoefficients{FT, D, SF}

Wind-dependent bulk transfer coefficients following the NCAR/Large & Yeager (2004) algorithm as used in OMIP-2 simulations.

Computes all three transfer coefficients (drag, heat, moisture) from a neutral drag polynomial with stability corrections via Monin-Obukhov similarity functions. The empirical heat and moisture coefficients correspond to the Stanton and Dalton numbers of L&Y eq. 6c-6d.

Pass as transfer_coefficients to CoefficientBasedFluxes:

ly = LargeYeagerTransferCoefficients()
fluxes = CoefficientBasedFluxes(transfer_coefficients = ly,
                                solver_stop_criteria = FixedIterations(5))

References:

  • Large, W.G. & Yeager, S.G. (2004): NCAR/TN-460+STR
source
NumericalEarth.EarthSystemModels.InterfaceComputations.MomentumBasedFrictionVelocityType
MomentumBasedFrictionVelocity

A friction velocity formulation that computes the friction velocity from momentum stresses.

The friction velocity is computed as:

\[u_* = \sqrt{\frac{|\boldsymbol{\tau}|}{\rho_o}}\]

where τ is the magnitude of the momentum stress vector and ρᵒᶜ is the ocean reference density.

Example

using NumericalEarth.EarthSystemModels: MomentumBasedFrictionVelocity

fv = MomentumBasedFrictionVelocity()

# output
MomentumBasedFrictionVelocity (computed from momentum stresses)
source
NumericalEarth.EarthSystemModels.InterfaceComputations.MomentumRoughnessLengthType
MomentumRoughnessLength(FT = Float64;
                        gravitational_acceleration = default_gravitational_acceleration,
                        maximum_roughness_length = 1.0,
                        air_kinematic_viscosity = 1.5e-5,
                        wave_formulation = 0.011,
                        smooth_wall_parameter = 0.11)

Construct a MomentumRoughnessLength object that represents the momentum roughness length that regulates the exchange of momentum, heat, and water vapor between the ocean and the atmosphere.

Keyword Arguments

  • gravitational_acceleration: The gravitational acceleration. Default: default_gravitational_acceleration.
  • maximum_roughness_length: The maximum roughness length. Default: 1e-1.
  • air_kinematic_viscosity: The air kinematic viscosity. Default: 1.5e-5.
  • wave_formulation: The gravity wave parameter. Default: 0.011.
  • smooth_wall_parameter: The smoothwallparameter parameter. Default: 0.11.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.ScalarRoughnessLengthType
ScalarRoughnessLength(FT = Float64;
                      air_kinematic_viscosity = temperature_dependent_viscosity,
                      reynolds_number_scaling_function = empirical_scaling_function,
                      maximum_roughness_length = 1.6e-4)

Construct a ScalarRoughnessLength object that represents the scalar roughness length that regulates the exchange of heat and water vapor between the ocean and the atmosphere.

Keyword Arguments

  • air_kinematic_viscosity::Function: The function to compute the air kinematic viscosity.
  • reynolds_number_scaling_function::Function: The function to compute the Reynolds number scaling factor.
  • maximum_roughness_length::Float: The maximum roughness length value. Defaults to 1.6e-4.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.SimilarityTheoryFluxesType
SimilarityTheoryFluxes(FT::DataType = Float64;
                       gravitational_acceleration = 9.81,
                       von_karman_constant = 0.4,
                       turbulent_prandtl_number = 1,
                       gustiness_parameter = 1.2,
                       minimum_gustiness = 0.01,
                       stability_functions = default_stability_functions(FT),
                       roughness_lengths = default_roughness_lengths(FT),
                       similarity_form = LogarithmicSimilarityProfile(),
                       solver_stop_criteria = nothing,
                       solver_tolerance = 1e-8,
                       solver_maxiter = 100)

SimilarityTheoryFluxes contains parameters and settings to calculate air-interface turbulent fluxes using Monin–Obukhov similarity theory.

Keyword Arguments

  • von_karman_constant: The von Karman constant. Default: 0.4.
  • turbulent_prandtl_number: The turbulent Prandtl number. Default: 1.
  • gustiness_parameter: Scaling factor for convective gustiness velocity. Default: 1.2.
  • minimum_gustiness: Minimum gustiness velocity [m/s], used as a floor in stable conditions where convective gustiness is zero. Default: 0.01.
  • stability_functions: The stability functions. Default: default_stability_functions(FT) that follow the formulation of Edson et al. (2013).
  • roughness_lengths: The roughness lengths used to calculate the characteristic scales for momentum, temperature and water vapor. Default: default_roughness_lengths(FT), formulation taken from Edson et al. (2013).
  • similarity_form: The type of similarity profile used to relate the atmospheric state to the interface fluxes / characteristic scales.
  • solver_tolerance: The tolerance for convergence. Default: 1e-8.
  • solver_maxiter: The maximum number of iterations. Default: 100.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.SkinTemperatureType
struct SkinTemperature

A type to represent the interface temperature used in the flux calculation. The interface temperature is calculated from the flux balance at the interface. In particular, the interface temperature $Tₛ$ is the root of:

\[F(Tₛ) - Jᵀ = 0\]

where $Jᵀ$ are the fluxes at the top of the interface (turbulent + radiative), and $F$ is the internal diffusive flux dependent on the interface temperature itself.

Note that all fluxes positive upwards.

source
NumericalEarth.EarthSystemModels.InterfaceComputations.ThreeEquationHeatFluxType
ThreeEquationHeatFlux(FT::DataType = Oceananigans.defaults.FloatType;
                      heat_transfer_coefficient = 0.0095,
                      salt_transfer_coefficient = heat_transfer_coefficient / 35,
                      friction_velocity = 0.002)

Construct a ThreeEquationHeatFlux with the specified parameters.

Default values follow Shi et al. (2021) with $R = \alpha_h / \alpha_s = 35$.

Keyword Arguments

  • heat_transfer_coefficient: turbulent heat exchange coefficient $\alpha_h$. Default: 0.0095.
  • salt_transfer_coefficient: turbulent salt exchange coefficient $\alpha_s$. Default: $\alpha_h / 35 \approx 0.000271$.
  • friction_velocity: friction velocity value or formulation. Default: 0.002.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.ThreeEquationHeatFluxType
ThreeEquationHeatFlux{FT, U}

Three-equation formulation for sea ice-ocean heat flux.

This formulation solves a coupled system for the interface temperature and salinity:

  1. Heat balance: $\rho c_p \gamma_T (T - T_b) = ℰ q$
  2. Salt balance: $\gamma_S (S - S_b) = q (S_b - S_i)$
  3. Freezing point: $T_b = T_m(S_b)$

where $T_b$ and $S_b$ are the interface temperature and salinity, $\gamma_T = \alpha_h u_*$ and $\gamma_S = \alpha_s u_*$ are turbulent exchange velocities, $L$ is the latent heat of fusion, and $q$ is the melt rate (computed, not input).

Fields

  • heat_transfer_coefficient::FT: turbulent heat exchange coefficient $\alpha_h$ (dimensionless)
  • salt_transfer_coefficient::FT: turbulent salt exchange coefficient $\alpha_s$ (dimensionless)
  • internal_heat_flux::FT: diffusive flux inside the sea ice (ConductiveFlux)
  • friction_velocity::U: friction velocity value or formulation (constant Number or MomentumBasedFrictionVelocity)

Example

using NumericalEarth.EarthSystemModels: ThreeEquationHeatFlux

flux = ThreeEquationHeatFlux()

# output
ThreeEquationHeatFlux{Nothing}
├── heat_transfer_coefficient: 0.0095
├── salt_transfer_coefficient: 0.00027142857142857144
└── friction_velocity: 0.002

References

  • Holland and Jenkins (1999): Holland, D. M., & Jenkins, A. (1999). Modeling thermodynamic ice–ocean interactions at the base of an ice shelf. Journal of Physical Oceanography, 29(8), 1787-1800.
  • Shi et al. (2021): Shi, X., Notz, D., Liu, J., Yang, H., & Lohmann, G. (2021). Sensitivity of Northern Hemisphere climate to ice-ocean interface heat flux parameterizations. Geosci. Model Dev., 14, 4891-4908.
source
NumericalEarth.EarthSystemModels.InterfaceComputations.compute_sea_ice_ocean_fluxes!Method
compute_sea_ice_ocean_fluxes!(coupled_model)

Compute heat, salt, and momentum fluxes at the sea ice-ocean interface.

This function computes:

  • Frazil heat flux: heat released when ocean temperature drops below freezing (all formulations)
  • Interface heat flux: heat flux from ocean to ice, computed using the specified formulation
  • Salt flux: salt exchange due to ice growth/melt
  • Momentum stresses: ice-ocean momentum transfer

The interface heat flux formulation is determined by coupled_model.interfaces.sea_ice_ocean_interface.flux_formulation.

source

Grids

NumericalEarth.Grids.PressureLevelVerticalDiscretizationType
PressureLevelVerticalDiscretization{G, Geo}

A vertical discretization for pressure-level reanalysis data on a LatitudeLongitudeGrid. Per-cell heights are geopotential[i, j, k] / gravitational_acceleration; znode(i, j, k, grid) returns that directly.

Two znodes paths:

  • znodes(grid, ℓ...) returns the column-mean z profile as a 1-D Vector{Float64} of length Nz — the representative axis that plot recipes, Lz, and length consumers expect when only the grid is in hand.
  • znodes(field) returns the column-mean Vector{Float64} when the field has no horizontal extent (Flat topology or Nothing horizontal locations, e.g. after mean(...; dims=(1, 2))); otherwise it returns a 3-D Field of per-cell heights at the field's own location.

geopotential (units m²/s²) is a 3-D Field or a TimeSeriesInterpolation over a FieldTimeSeries. The former gives a static z-coordinate; the latter gives a time-evolving one driven by an attached Clock.

The LatitudeLongitudeGrid constructor needs a value for Lz; we compute it as extrema(geopotential) / g inside generate_coordinate.

source
NumericalEarth.Grids.PressureLevelVerticalDiscretizationMethod
PressureLevelVerticalDiscretization(geopotential;
                                    gravitational_acceleration,
                                    surface_geopotential = nothing)

Build a discretization backed by per-column geopotential (m²/s²). znode divides by gravitational_acceleration at read time.

If surface_geopotential is provided (a 2-D Field, m²/s²), columns are clipped so that geopotential[i,j,k] ≥ surface_geopotential[i,j]. Required when the source is ERA5 pressure-level data, because sub-surface levels are filled with non-physical extrapolations that would break the column-monotonicity assumed by _fractional_indices.

source

InitialConditions

Lands

NumericalEarth.Lands.PrescribedLandMethod
PrescribedLand(freshwater_flux; clock=nothing)

Return a PrescribedLand component from a NamedTuple of FieldTimeSeries representing freshwater fluxes (e.g. rivers, icebergs).

If clock is not provided, defaults to a Clock whose time type matches the element type of freshwater_flux.

source

Oceans

NumericalEarth.Oceans.PrescribedOceanType
PrescribedOcean(grid, times=[zero(grid)];
                density = 1025,
                heat_capacity = 4000,
                clock = Clock{FT}(time=0),
                sea_surface_temperature = default_prescribed_sst(grid, times),
                sea_surface_salinity = default_prescribed_sss(grid, times),
                velocities = default_prescribed_velocities(grid, times))

A prescribed ocean component for EarthSystemModel with sea surface temperature, salinity, and velocities prescribed as FieldTimeSeries.

The ocean state does not evolve in response to surface fluxes — it follows the prescribed data. Surface fluxes are still computed so the atmosphere feels the ocean.

Arguments

  • grid: An Oceananigans grid for the ocean domain.
  • times: Time instances for the prescribed data. Default: [zero(grid)] (constant).

Keyword Arguments

  • density: Seawater density in kg/m³. Default: 1025.
  • heat_capacity: Seawater specific heat capacity in J/(kg·K). Default: 4000.
  • clock: Clock for tracking ocean time. Default: Clock{FT}(time=0).
  • sea_surface_temperature: FieldTimeSeries for SST.
  • sea_surface_salinity: FieldTimeSeries for SSS.
  • velocities: NamedTuple of FieldTimeSeries for (u, v).
source
NumericalEarth.Oceans.SlabOceanMethod
SlabOcean(grid;
          FT = eltype(grid),
          depth = 50,
          density = 1025,
          heat_capacity = 4000,
          clock = Clock{FT}(time=0))

A slab ocean component for EarthSystemModel. Represents the ocean as a single well-mixed layer of fixed depth, with temperature evolving in response to net surface heat fluxes:

∂T/∂t = Q / (ρ cₚ H)

where Q is the net downward surface heat flux (W/m²), ρ is the seawater density, cₚ is the heat capacity, and H is the slab depth.

Internally, the EarthSystemModel coupling assembles a temperature flux Jᵀ (in K m/s), so the tendency is computed as ∂T/∂t = -Jᵀ / H.

Use set!(ocean, T=value) to initialize the slab temperature after construction.

Arguments

  • grid: An Oceananigans grid for the slab ocean domain.

Keyword Arguments

  • FT: Floating point type for the slab ocean fields. Default: eltype(grid).
  • depth: Depth of the slab in meters. Default: 50.
  • density: Seawater density in kg/m³. Default: 1025.
  • heat_capacity: Seawater specific heat capacity in J/(kg·K). Default: 4000.
  • clock: Clock for tracking slab ocean time. Default: Clock{FT}(time=0) where FT is the float type of the grid.
source
NumericalEarth.Oceans.ocean_simulationMethod
ocean_simulation(grid; model = :hydrostatic, kwargs...)

Construct and return an ocean simulation tailored to grid. The model keyword selects the underlying Oceananigans model formulation:

  • :hydrostatic (default) — builds a HydrostaticFreeSurfaceModel-based simulation with a free surface, CATKE vertical mixing, quadratic bottom drag, and a TEOS-10 equation of state. See hydrostatic_ocean_simulation for the full kwarg list.

  • :nonhydrostatic — builds a NonhydrostaticModel-based simulation suitable for LES (full 3D pressure, no free surface, no barotropic forcing). See nonhydrostatic_ocean_simulation for the full kwarg list.

Remaining kwargs are forwarded to the per-model builder; an unknown kwarg for the selected model raises the usual MethodError.

Examples

julia> using NumericalEarth, Oceananigans

julia> grid = RectilinearGrid(size = (10, 10, 6), extent = (1, 1, 1), halo=(6, 6, 5));

julia> ocean = ocean_simulation(grid)
Simulation of HydrostaticFreeSurfaceModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── Next time step: 30 minutes
├── run_wall_time: 0 seconds
├── run_wall_time / iteration: NaN days
├── stop_time: Inf days
├── stop_iteration: Inf
├── wall_time_limit: Inf
├── minimum_relative_step: 0.0
├── callbacks: OrderedDict with 4 entries:
│   ├── stop_time_exceeded => Callback of stop_time_exceeded on IterationInterval(1)
│   ├── stop_iteration_exceeded => Callback of stop_iteration_exceeded on IterationInterval(1)
│   ├── wall_time_limit_exceeded => Callback of wall_time_limit_exceeded on IterationInterval(1)
│   └── nan_checker => Callback of NaNChecker for u on IterationInterval(100)
└── output_writers: OrderedDict with no entries

julia> les = ocean_simulation(grid; model = :nonhydrostatic, Δt = 2)
Simulation of NonhydrostaticModel{CPU, RectilinearGrid}(time = 0 seconds, iteration = 0)
├── Next time step: 2 seconds
├── run_wall_time: 0 seconds
├── run_wall_time / iteration: NaN days
├── stop_time: Inf days
├── stop_iteration: Inf
├── wall_time_limit: Inf
├── minimum_relative_step: 0.0
├── callbacks: OrderedDict with 4 entries:
│   ├── stop_time_exceeded => Callback of stop_time_exceeded on IterationInterval(1)
│   ├── stop_iteration_exceeded => Callback of stop_iteration_exceeded on IterationInterval(1)
│   ├── wall_time_limit_exceeded => Callback of wall_time_limit_exceeded on IterationInterval(1)
│   └── nan_checker => Callback of NaNChecker for u on IterationInterval(100)
└── output_writers: OrderedDict with no entries
source

Radiations

NumericalEarth.Radiations.InterfaceRadiationFluxType
InterfaceRadiationFlux{F}

Container for the diagnostic radiative fluxes at an air–surface interface. The same struct type is instantiated per surface (ocean, sea ice, snow, ...).

Fields

  • upwelling_longwave :: F ϵσT⁴
  • downwelling_longwave :: F ϵℐꜜˡʷ (absorbed by the surface)
  • downwelling_shortwave :: F (1−α)ℐꜜˢʷ (transmitted into the surface)
source
NumericalEarth.Radiations.LatitudeDependentAlbedoType
LatitudeDependentAlbedo([FT::DataType=Float64]; diffuse = 0.069, direct = 0.011)

Constructs a LatitudeDependentAlbedo object. The albedo of the ocean surface is assumed to be a function of the latitude, obeying the following formula (Large and Yeager, 2009):

α(φ) = α.diffuse - α.direct * cos(2φ)

where φ is the latitude, α.diffuse is the diffuse albedo, and α.direct is the direct albedo.

Arguments

  • FT::DataType: The data type of the albedo values. Default is Float64.

Keyword Arguments

  • diffuse: The diffuse albedo value. Default is 0.069.
  • direct: The direct albedo value. Default is 0.011.
source
NumericalEarth.Radiations.PrescribedRadiationType
PrescribedRadiation(grid, times = [zero(grid)]; kwargs...)

Construct a PrescribedRadiation with zero downwelling shortwave and longwave fields on grid. Useful for emission-only mode (the surface radiates ϵσT⁴ but no incoming radiation is absorbed). All other keyword arguments are forwarded to the FTS-form constructor.

source
NumericalEarth.Radiations.PrescribedRadiationType
PrescribedRadiation{G, T, FT, SW, LW, S, IF, TI}

Top-level radiation component holding prescribed downwelling shortwave and longwave radiation as FieldTimeSeries, plus per-surface radiative properties (albedo, emissivity) and the Stefan–Boltzmann constant. Diagnostic radiative fluxes (one InterfaceRadiationFlux per surface) are populated by the apply_air_sea_*_radiative_fluxes! kernels at every step; interface_fluxes is nothing until the radiation is paired with an EarthSystemModel (which allocates the per-surface buffers on the exchange grid).

source
NumericalEarth.Radiations.PrescribedRadiationMethod
PrescribedRadiation(downwelling_shortwave, downwelling_longwave;
                    clock = nothing,
                    ocean_surface = SurfaceRadiationProperties(0.05, 0.97),
                    sea_ice_surface = SurfaceRadiationProperties(0.7, 1.0),
                    snow_surface = nothing,
                    land_surface = nothing,
                    stefan_boltzmann_constant = default_stefan_boltzmann_constant)

Construct a PrescribedRadiation component from FieldTimeSeries of downwelling shortwave and longwave radiation. Grid + times are inferred from the shortwave FTS.

Pass *_surface = nothing to omit that surface from surface_properties.

source
NumericalEarth.Radiations.SeaIceAlbedoType
SeaIceAlbedo{FT, HI, HS, TS}

Sea ice albedo parameterization following the CCSM3 scheme (Briegleb et al. 2004).

Computes broadband albedo as a function of ice thickness, snow depth, and surface temperature. The scheme blends between bare ice and snow-covered albedos, with a temperature-dependent reduction near the melting point to implicitly represent melt pond formation.

Algorithm:

  1. Base cold albedos: bare ice (0.53) and snow-covered (0.82)
  2. Temperature reduction within 1C of melting: Δαice = 0.075, Δαsnow = 0.10
  3. Thin-ice transition to ocean albedo below h_amin = 0.5 m
  4. Snow cover interpolation: full snow albedo at hsnow > hsmin = 0.02 m

References:

  • Briegleb, B.P., C.M. Bitz, E.C. Hunke, W.H. Lipscomb, and M.M. Schramm (2004): Scientific description of the sea ice component in CCSM3. NCAR Tech Note.
  • Briegleb, B.P. and B. Light (2007): NCAR/TN-472+STR.
source
NumericalEarth.Radiations.SeaIceAlbedoMethod
SeaIceAlbedo(ice_thickness, snow_thickness, surface_temperature;
                  ice_albedo = 0.54,
                  snow_albedo = 0.83,
                  ice_melt_reduction = 0.075,
                  snow_melt_reduction = 0.10,
                  melting_temperature = 0.0,
                  temperature_range = 1.0,
                  ocean_albedo = 0.06,
                  minimum_ice_thickness = 0.5,
                  minimum_snow_depth = 0.02)

Construct a CCSM3 sea ice albedo parameterization. Requires references to the sea ice model's ice thickness, snow thickness, and surface temperature fields.

Broadband albedos are approximate averages of the visible and near-IR bands weighted by solar spectrum (52% visible, 48% near-IR):

  • ice_albedo ≈ 0.52 x 0.73 + 0.48 x 0.33 ≈ 0.54
  • snow_albedo ≈ 0.52 x 0.96 + 0.48 x 0.68 ≈ 0.83
source
NumericalEarth.Radiations.SurfaceRadiationPropertiesMethod
SurfaceRadiationProperties(albedo, emissivity)

Bundle the radiative properties of a single surface (ocean, sea ice, snow, land) that participate in radiative flux computation: shortwave reflectivity (albedo) and longwave emissivity (emissivity).

albedo may be a Number, a LatitudeDependentAlbedo, a TabulatedAlbedo, or any other object for which stateindex is defined. emissivity may be a Number or any other stateindex-able object.

source
NumericalEarth.Radiations.TabulatedAlbedoType
TabulatedAlbedo(arch = CPU(), FT = Float64;
                S₀ = convert(FT, 1365),
                α_table  = α_payne,
                φ_values = (0:2:90) ./ 180 * π,
                𝓉_values = 0:0.05:1)

Constructs a TabulatedAlbedo object that interpolated the albedo from a value table α_table that is function of latitude φ and atmospheric transimissivity 𝓉.

Note: TabulatedAlbedo assumes that the latitude and the transissivity in the table are uniformly spaced.

The transmissivity of the atmosphere is calculated as the ratio of the downwelling solar radiation to the maximum possible downwelling solar radiation for a transparent atmosphere, function of hour of the day, latitude, and day in the year.

Arguments

  • arch: The architecture to use. Default: CPU().
  • FT: The floating-point type to use. Default: Float64.

Keyword Arguments

  • S₀: The solar constant. Default: convert(FT, 1365).
  • α_table: The table of albedo values. Default: α_payne.
  • φ_values: The latitude values for the table. Default: (0:2:90) ./ 180 * π.
  • 𝓉_values: The transmissivity values for the table. Default: 0:0.05:1.
source

SeaIces

NumericalEarth.SeaIces.FreezingLimitedOceanTemperatureType
FreezingLimitedOceanTemperature(FT=Float64; liquidus=LinearLiquidus(FT))

The minimal possible sea ice representation, clipping the temperature below to the freezing point. Not really a "model" per se, however, it is the most simple way to make sure that temperature does not dip below freezing.

The melting temperature is a function of salinity and is controlled by the liquidus.

source
NumericalEarth.SeaIces.sea_ice_simulationFunction
sea_ice_simulation(grid, ocean=nothing;
                   Δt = 5minutes,
                   ice_salinity = 4, # psu
                   advection = nothing,
                   tracers = (),
                   ice_heat_capacity = 2100, # J kg⁻¹ K⁻¹
                   ice_consolidation_thickness = 0.05, # m
                   sea_ice_density = 900, # kg m⁻³
                   snow_density = 330, # kg m⁻³
                   dynamics = sea_ice_dynamics(grid, ocean),
                   bottom_heat_boundary_condition = nothing,
                   top_heat_boundary_condition = nothing,
                   timestepper = :SplitRungeKutta3,
                   phase_transitions = PhaseTransitions(eltype(grid);
                                                        heat_capacity=ice_heat_capacity,
                                                        density=sea_ice_density),
                   conductivity = 2, # W m⁻¹ K⁻¹
                   internal_heat_flux = ConductiveFlux(; conductivity),
                   snow_thermodynamics = default_snow_thermodynamics(grid))

Construct a sea ice simulation with the given grid and optional ocean simulation. The sea ice model is configured with a slab thermodynamics, Elasto-Visco-Plastic rheology, and a SplitExplicit Runge-Kutta 3rd order time stepper by default. The thermodynamics include conductive internal heat flux, and the option to specify top and bottom heat boundary conditions. The dynamics include a semi-implicit ocean stress formulation, with the option to specify a free drift velocity.

Arguments

  • grid: the grid on which to build the sea ice model
  • ocean: optional ocean simulation to provide surface velocities and salinity for the sea ice

Keyword Arguments

  • Δt: time step for the sea ice simulation
  • ice_salinity: salinity of the sea ice (psu)
  • advection: optional advection scheme for the sea ice model; if nothing (default), no advection is applied and only thermodynamics evolve the sea ice state
  • tracers: optional tracers to include in the sea ice model
  • ice_heat_capacity: heat capacity of the sea ice (J kg⁻¹ K⁻¹)
  • ice_consolidation_thickness: thickness threshold for sea ice consolidation (m)
  • sea_ice_density: density of the sea ice (kg m⁻³)
  • snow_density: density of the snow (kg m⁻³)
  • dynamics: sea ice dynamics model to use (default is sea_ice_dynamics(grid, ocean))
  • bottom_heat_boundary_condition: heat boundary condition at the ice-ocean interface (default is IceWaterThermalEquilibrium with ocean surface salinity)
  • top_heat_boundary_condition: heat boundary condition at the ice-atmosphere interface (default is a prescribed temperature calculated in the flux computation)
  • timestepper: time stepper to use for the sea ice model (default is :SplitRungeKutta3)
  • phase_transitions: phase transition properties for the sea ice (default is a PhaseTransitions with specified heat capacity and density)
  • conductivity: thermal conductivity for the internal heat flux (W m⁻¹ K⁻¹)
  • internal_heat_flux: internal heat flux formulation for the sea ice (default is a ConductiveFlux with specified conductivity)
  • snow_thermodynamics: thermodynamics for the snow layer (default is a slab thermodynamics with specified conductivity and prescribed temperature)
source

NumericalEarthBreezeExt

NumericalEarthSpeedyWeatherExt

NumericalEarthVerosExt