Public Documentation

<!– Title –> <h1 align="center"> ClimaOcean.jl </h1>

<!– description –> <p align="center"> <strong>🌎 A framework for realistic ocean-only and coupled ocean + sea-ice simulations driven by prescribed atmospheres and based on <a href=https://github.com/CliMA/Oceananigans.jl>Oceananigans</a> and <a href=https://github.com/CliMA/ClimaSeaIce.jl>ClimaSeaIce</a></strong>. </p>

Installation

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

julia> using Pkg

julia> Pkg.add("ClimaOcean")

julia> Pkg.instantiate()

Use Pkg.add("url=https://github.com/CliMA/ClimaOcean.jl.git", rev="main") to install the latest version of ClimaOcean. For more information, see the documentation for Pkg.jl.

Why? What's the difference between ClimaOcean and Oceananigans?

Oceananigans is a general-purpose library for ocean-flavored fluid dynamics. ClimaOcean implements a framework for driving realistic Oceananigans simulations with prescribed atmospheres, and coupling them to prognostic sea ice simulations.

A core abstraction: ClimaOcean.OceanSeaIceModel

Our system for realistic modeling is anchored by ClimaOcean.OceanSeaIceModel, which encapsulates the ocean simulation, sea ice simulation, prescribed atmospheric state, and specifies how the three communicate. To illustrate how OceanSeaIceModel works we set up a simulation on a grid with 10 vertical levels and 1/4-degree horizontal resolution:

using Oceananigans
using Oceananigans.Units
using Dates, CFTime
import ClimaOcean

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

bathymetry = ClimaOcean.regrid_bathymetry(grid) # builds gridded bathymetry based on ETOPO1
grid = ImmersedBoundaryGrid(grid, GridFittedBottom(bathymetry))

# Build an ocean simulation initialized to the ECCO state estimate on Jan 1, 1993
ocean = ClimaOcean.ocean_simulation(grid)
date  = DateTimeProlepticGregorian(1993, 1, 1)
set!(ocean.model, T = ClimaOcean.ECCOMetadata(:temperature; date),
                  S = ClimaOcean.ECCOMetadata(:salinity; date))

# Build and run an OceanSeaIceModel (with no sea ice component) forced by JRA55 reanalysis
atmosphere = ClimaOcean.JRA55_prescribed_atmosphere(arch)
coupled_model = ClimaOcean.OceanSeaIceModel(ocean; atmosphere)
simulation = Simulation(coupled_model, Δt=5minutes, stop_time=30days)
run!(simulation)

The simulation above 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 in our scripts. For example, to plot the surface speed at the end of the simulation we write

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

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

which produces

Additional features: a utility for ocean_simulations and data wrangling

A second core abstraction in ClimaOcean is ocean_simulation. ocean_simulation configures an Oceananigans model for realistic simulations including temperature and salinity, the TEOS-10 equation of state, boundary conditions to store computed air-sea fluxes, the automatically-calibrated turbulence closure CATKEVerticalDiffusivity, and the WENOVectorInvariant advection scheme for mesoscale-turbulence-resolving simulations.

ClimaOcean also provides convenience features for wrangling datasets of bathymetry, ocean temperature, salinity, ocean velocity fields, and prescribed atmospheric states.

ClimaOcean is built on top of Oceananigans and ClimaSeaIce, so it's important that ClimaOcean users become proficient with Oceananigans. Note that though ClimaOcean is currently focused on hydrostatic modeling with Oceananigans.HydrostaticFreeSurfaceModel, realistic nonhydrostatic modeling is also within the scope of this package.

@distribute for i in iterable
    ...
end

Distribute a for loop among different ranks

@handshake exs...

perform exs on all ranks, but only one rank at a time, where ranks r2 > r1 wait for rank r1 to finish before executing exs

@onrank rank, exs...

Perform exp only on rank rank Other ranks will wait for the root rank to finish before continuing. The expression is run anyways if MPI in not initialized

@root exs...

Perform exs only on rank 0, otherwise know as "root" rank. Other ranks will wait for the root rank to finish before continuing

MixedLayerDepthField(grid, tracers, buoyancy_model; Δb = 3e-4, field_kw...)

Return a reduced Field{Center, Center, Nothing} that represents mixed layer depth for model, based on a buoyancy differential criterion. The mixed layer depth is defined as the depth $h$ for which

\[b(z=0) - b(z=-h) = Δb\]

This criterion is solved by integrating downwards and linearly interpolating to find h, assuming that $b$ decreases with depth.

Keyword arguments

• Δb: Buoyancy differential used to calculate mixed layer depth
• field_kw: Keyword arguments passed to Field.

Example

h = MixedLayerDepth(model)
compute!(h) # compute mixed layer depth

struct ECCOMetadata{D, V}

Metadata information for an ECCO dataset:

• name: The name of the dataset.
• dates: The dates of the dataset, in an AbstractCFDateTime format.
• version: The version of the dataset, could be ECCO2Monthly, ECCO2Daily, or ECCO4Monthly.
• dir: The directory where the dataset is stored.

ECCOMetadata(name::Symbol; 
             dates = DateTimeProlepticGregorian(1993, 1, 1),
             version = ECCO4Monthly(),
             dir = download_ECCO_cache)

Construct an ECCOMetadata object with the specified parameters.

Arguments

• name::Symbol: The name of the metadata.

Keyword Arguments

• dates: The date(s) of the metadata. Note this can either be a single date, representing a snapshot, or a range of dates, representing a time-series. Default: DateTimeProlepticGregorian(1993, 1, 1).

• version: The data version. Supported versions are ECCO2Monthly(), ECCO2Daily(), or ECCO4Monthly().

• dir: The directory of the data file. Default: download_ECCO_cache.

ECCORestoring([arch=CPU(),]
              variable_name::Symbol;
              version=ECCO4Monthly(),
              dates=all_ECCO_dates(version),
              dates = all_ECCO_dates(version), 
              time_indices_in_memory = 2, 
              time_indexing = Cyclical(),
              mask = 1,
              rate = 1,
              grid = nothing,
              inpainting = NearestNeighborInpainting(prod(size(metadata))))

Create a forcing term that restores to values stored in an ECCO field time series. The restoring is applied as a forcing on the right hand side of the evolution equations calculated as

F = mask ⋅ rate ⋅ (ECCO_variable - simulation_variable[i, j, k])

where ECCO_variable is linearly interpolated in space and time from the ECCO dataset of choice to the simulation grid and time.

Arguments

• arch: The architecture. Typically CPU() or GPU(). Default: CPU().

• variable_name: The name of the variable to restore. Choices include:

  • :temperature,
  • :salinity,
  • :u_velocity,
  • :v_velocity,
  • :sea_ice_thickness,
  • :sea_ice_area_fraction.

Keyword Arguments

• version: The version of the ECCO dataset. Default: ECCO4Monthly().

• dates: The dates to use for the ECCO dataset. Default: all_ECCO_dates(version).

• time_indices_in_memory: The number of time indices to keep in memory; trade-off between performance and memory footprint.

• time_indexing: The time indexing scheme for the field time series≥

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

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

• time_indices_in_memory: how many time instances are loaded in memory; the remaining are loaded lazily.

• grid: if not a nothing, the ECCO data is directly interpolated on the grid.

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

• grid: If isnothing(grid), ECCO data is interpolated on-the-fly to the simulation grid. If !isnothing(grid), ECCO data is pre-interpolated to grid. Default: nothing.

It is possible to also pass an ECCOMetadata type as the first argument without the need for the variable_name argument and the version and dates keyword arguments.

LinearlyTaperedPolarMask(; northern = (70,   75),
                           southern = (-75, -70),
                           z = (-20, 0))

Build a mask that is linearly tapered in latitude inbetween the northern and southern edges. The mask is constant in depth between the z and is equal to zero everywhere else. The mask has the following functional form:

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

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

mask = within_depth ? max(n, s, 0) : 0

ECCO_field(metadata::ECCOMetadata;
           architecture = CPU(),
           horizontal_halo = (3, 3))

Retrieve the ecco field corresponding to metadata. The data is loaded from filename on architecture with horizontal_halo in the x and y direction. The halo in the z-direction is one.

ECCO_immersed_grid(metadata, architecture = CPU())

Compute the ImmersedBoundaryGrid for metadata with a bottom height field that is defined by the first non-missing value from the bottom up.

ECCO_mask(architecture = CPU(); minimum_value = Float32(-1e5))

A boolean field where true represents a missing value in the ECCO dataset.

regrid_bathymetry(target_grid;
                  url = "https://www.ngdc.noaa.gov/thredds/fileServer/global/ETOPO2022/60s/60s_surface_elev_netcdf", 
                  height_above_water = <none>,
                  minimum_depth = 0,
                  dir = download_cache,
                  filename = "ETOPO_2022_v1_60s_N90W180_surface.nc")

Regrid bathymetry associated with the NetCDF file at path = joinpath(dir, filename) to target_grid. If path does not exist, then a download is attempted from joinpath(url, filename).

Arguments

• target_grid: grid to interpolate onto

Keyword Arguments

• height_above_water: limits the maximum height of above-water topography (where h > 0). If nothing the original topography is retained

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

• dir: directory of the bathymetry-containing file

• filename: file containing bathymetric data. Must be netcdf with fields: (1) lat vector of latitude nodes (2) lon vector of longitude nodes (3) z matrix of depth values

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

  Example: interpolating from a 400x200 grid to a 100x100 grid in 4 passes involves:

  • 400x200 → 325x175
  • 325x175 → 250x150
  • 250x150 → 175x125
  • 175x125 → 100x100

  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 will not help and no intermediate steps will be 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 will retain only the largest basin. Default: Inf, which does not remove any basins.

retrieve_bathymetry(grid, filename; kw...)

Retrieve the bathymetry data from a file or generate it using a grid and save it to a file.

Arguments

• grid: The grid used to generate the bathymetry data.
• filename: The name of the file to read or save the bathymetry data.
• kw...: Additional keyword arguments.

Returns

• bottom_height: The retrieved or generated bathymetry data.

If the specified file exists, the function reads the bathymetry data from the file. Otherwise, it generates the bathymetry data using the provided grid and saves it to the file before returning it.