Library
The public user interface.
Architectures
PlanktonIndividuals.Architectures.Architecture
— TypeArchitecture
Abstract type for architectures supported by PlanktonIndividuals.
PlanktonIndividuals.Architectures.CPU
— TypeCPU <: Architecture
Run PlanktonIndividuals on one CPU node.
PlanktonIndividuals.Architectures.GPU
— TypeGPU <: Architecture
Run PlanktonIndividuals on one CUDA GPU node.
Grids
PlanktonIndividuals.Grids.AbstractGrid
— TypeAbstractGrid{FT, TX, TY, TZ}
Abstract type for grids with elements of type FT
and topology {TX, TY, TZ}
.
PlanktonIndividuals.Grids.Bounded
— TypeBounded
Grid topology for bounded dimensions.
PlanktonIndividuals.Grids.Periodic
— TypePeriodic
Grid topology for periodic dimensions.
PlanktonIndividuals.Grids.RectilinearGrid
— MethodRectilinearGrid(;size, x, y, z,
FT = Float32,
topology = (Periodic, Periodic, Bounded),
landmask = nothing,
halo = (2, 2, 2))
Creats a RectilinearGrid
struct with size = (Nx, Ny, Nz)
grid points. x
and y
directions must be regular spaced, z
direction can be vertically stretched or regular spaced.
Keyword Arguments (Required)
size
: A tuple prescribing the number of grid points.size
is a 3-tuple no matter for 3D, 2D, or 1D model.x
andy
: A 2-tuple that specify the start and end points of the domain.z
: is either a (1) 1D array that specifies the locations of cell faces in z direction, or (2) 2-tuples that specify the start and end points of the domain. Vertical indexing starts from surface and use negative numbers for depth.
Keyword Arguments (Optional)
FT
: Floating point data type. Default:Float32
.topology
: A 3-tuple specifying the topology of the domain. The topology can be either Periodic or Bounded in each direction.landmask
: a 3-dimentional array to indicate where the land is.halo
: A tuple of integers that specifies the size of the halo region of cells surrounding the physical interior for each direction.halo
is a 3-tuple no matter for 3D, 2D, or 1D model. At least 2 halo points are needed for DST3FL advection scheme.
PlanktonIndividuals.Grids.LatLonGrid
— MethodLatLonGrid(;size, lat, lon, z,
FT = Float32,
radius = 6370.0e3,
landmask = nothing,
halo = (2, 2, 2))
Creats a LatLonGrid
struct with size = (Nx, Ny, Nz)
grid points.
Keyword Arguments (Required)
size
: A tuple prescribing the number of grid points.size
is a 3-tuple no matter for 3D, 2D, or 1D model.lat
: A 2-tuple specifying the startind and ending points in latitudinal direction. Possible values are from -80 (80S) to 80 (80N).lon
: A 2-tuple specifying the startind and ending points in longitudinal direction. Possible values are from -180 (180W) to 180 (180E).z
: is either a (1) 1D array that specifies the locations of cell faces in z direction, or (2) 2-tuples that specify the start and end points of the domain. Vertical indexing starts from surface and use negative numbers for depth.
Keyword Arguments (Optional)
FT
: Floating point data type. Default:Float32
.radius
: Specify the radius of the Earth used in the model, 6370.0e3 meters by default.landmask
: a 3-dimentional array to indicate where the land is.halo
: A tuple of integers that specifies the size of the halo region of cells surrounding the physical interior for each direction.halo
is a 3-tuple no matter for 3D, 2D, or 1D model. At least 2 halo points are needed for DST3FL advection scheme.
PlanktonIndividuals.Grids.LoadLatLonGrid
— MethodLoadLatLonGrid(;grid_info, size, lat, lon,
FT = Float32,
landmask = nothing,
halo=(2,2,2))
Creats a LatLonGrid
struct with size = (Nx, Ny, Nz)
grid points.
Keyword Arguments (Required)
grid_info
: A NamedTuple contains external grid information (e.g. from MITgcm), please refer to documentation for the required format.size
: A tuple prescribing the number of grid points.size
is a 3-tuple no matter for 3D, 2D, or 1D model.lat
: A 2-tuple specifying the startind and ending points in latitudinal direction. Possible values are from -80 (80S) to 80 (80N).lon
: A 2-tuple specifying the startind and ending points in longitudinal direction. Possible values are from -180 (180W) to 180 (180E).
Keyword Arguments (Optional)
FT
: Floating point data type. Default:Float32
.landmask
: a 3-dimentional array to indicate where the land is.halo
: A tuple of integers that specifies the size of the halo region of cells surrounding the physical interior for each direction.halo
is a 3-tuple no matter for 3D, 2D, or 1D model. At least 2 halo points are needed for DST3FL advection scheme.
Biogeochemistry
PlanktonIndividuals.Biogeochemistry.Field
— MethodField(arch::Architecture, grid::AbstractGrid, FT::DataType; bcs = default_bcs())
Construct a Field
on grid
with data and boundary conditions on architecture arch
with DataType FT
.
PlanktonIndividuals.Biogeochemistry.default_tracer_init
— Methoddefault_tracer_init()
Generate defalut bgc tracer initial conditions.
PlanktonIndividuals.Biogeochemistry.generate_tracers
— Methodgenerate_tracers(arch, grid, source, FT)
Set up initial bgc tracer fields according to grid
.
Arguments
arch
:CPU()
orGPU()
. The computer architecture used to time-stepmodel
.grid
: The resolution and discrete geometry on which nutrient fields are solved.source
: ANamedTuple
containing 10 numbers each of which is the uniform initial condition of one tracer, or aDict
containing the file paths pointing to the files of nutrient initial conditions.FT
: Floating point data type. Default:Float32
.
Parameters
PlanktonIndividuals.Parameters.abiotic_params_default
— Methodabiotic_params_default(N::Int64)
Generate default abiotic particle parameter values based on species number N
.
PlanktonIndividuals.Parameters.bgc_params_default
— Methodbgc_params_default(FT::DataType)
Generate default biogeochemical parameter values
PlanktonIndividuals.Parameters.default_PARF
— Methoddefault_PARF(grid, ΔT, iterations)
Generate default hourly surface PAR.
PlanktonIndividuals.Parameters.default_temperature
— Methoddefault_temperature(grid, ΔT, iterations)
Generate default hourly temperature.
PlanktonIndividuals.Parameters.phyt_params_default
— Methodphyt_params_default(N::Int64, mode::AbstractMode)
Generate default phytoplankton parameter values based on AbstractMode
and species number N
.
PlanktonIndividuals.Parameters.phyt_params_default
— Methodphyt_params_default(N::Int64, mode::AbstractMode)
Generate default phytoplankton parameter values based on AbstractMode
and species number N
.
PlanktonIndividuals.Parameters.phyt_params_default
— Methodphyt_params_default(N::Int64, mode::AbstractMode)
Generate default phytoplankton parameter values based on AbstractMode
and species number N
.
PlanktonIndividuals.Parameters.update_abiotic_params
— Methodupdate_abiotic_params(tmp::Dict, FT::DataType; N::Int64)
Update parameter values based on a Dict
provided by user Keyword Arguments =================
tmp
is aDict
containing the parameters needed to be upadatedFT
: Floating point data type. Default:Float32
.N
is aInt64
indicating the number of species
PlanktonIndividuals.Parameters.update_bgc_params
— Methodupdate_bgc_params(tmp::Ditc, FT::DataType)
Update parameter values based on a Dict
provided by user
Keyword Arguments
tmp
: aDict
containing the parameters needed to be upadatedFT
: Floating point data type. Default:Float32
.
PlanktonIndividuals.Parameters.update_phyt_params
— Methodupdate_phyt_params(tmp::Dict, FT::DataType; N::Int64, mode::AbstractMode)
Update parameter values based on a Dict
provided by user Keyword Arguments =================
tmp
is aDict
containing the parameters needed to be upadatedFT
: Floating point data type. Default:Float32
.N
is aInt64
indicating the number of speciesmode
is the mode of phytoplankton physiology resolved in the model
Diagnostics
PlanktonIndividuals.Diagnostics.PlanktonDiagnostics
— MethodPlanktonDiagnostics(model; tracer=(:PAR, :NH4, :NO3, :DOC),
phytoplankton = (:num, :graz, :mort, :dvid, :ptc),
abiotic_particle = (:num),
time_interval = 1)
Generate a PlanktonDiagnostics
structure.
Keyword Arguments (Optional)
tracer
: aTuple
containing the names of nutrient fields to be diagnosed.phytoplankton
: aTuple
containing the names of physiological processes of phytoplankton individuals to be diagnosed.abiotic_particle
: aTuple
containing the names of state variables of abiotic particles to be diagnosed.iteration_interval
: The number of timesteps that diagnostics is averaged, 1 iteration by default.
Model
PlanktonIndividuals.CarbonMode
— TypeCarbonMode <: AbstractMode
Type for the phytoplankton physiology mode which only resolves carbon quota.
PlanktonIndividuals.IronEnergyMode
— TypeIronEnergyMode <: AbstractMode
Type for the phytoplankton physiology mode which resolves carbon, nitrogen, phosphorus, and iron quotas. This mode also resolves energy.
PlanktonIndividuals.MacroMolecularMode
— TypeMacroMolecularMode <: AbstractMode
Type for the phytoplankton physiology mode which resolves marco-molecules.
PlanktonIndividuals.QuotaMode
— TypeQuotaMode <: AbstractMode
Type for the phytoplankton physiology mode which resolves carbon, nitrogen, and phosphorus quotas.
PlanktonIndividuals.Model.PlanktonModel
— MethodPlanktonModel(arch::Architecture, grid::AbstractGrid;
FT = Float32,
mode = QuotaMode(),
max_individuals::Int = 8*1024,
bgc_params = nothing,
tracer_initial = default_tracer_init(),
phyto = nothing,
abiotic = nothing,
t::AbstractFloat = 0.0f0,
)
Generate a PlanktonModel
data structure.
Keyword Arguments (Required)
arch
:CPU()
orGPU()
. Computer architecture being used to run the model.grid
: aAbstractGrid
structure. Discrete grid for the model (resolution and geometry).
Keyword Arguments (Optional)
FT
: Floating point data type. Default:Float32
.mode
: Phytoplankton physiology mode, choose among CarbonMode(), QuotaMode(), or MacroMolecularMode().max_individuals
: Maximum number of individuals for each species the model can hold, usually take the maximum of all the species and apply a factor to account for the growth of individuals during one simulation.bgc_params
: Parameter set for biogeochemical processes modeled in the model, use default ifnothing
, useDict
to update parameters, the format and names of parameters can be found by runningbgc_params_default()
.tracer_initial
: The source of initial conditions of tracer fields, should be either aNamedTuple
or aDict
containing the file paths pointing to the files of nutrient initial conditions.phyto
: nothing or aphyto_setup
. Whether to use default setup of phytoplankton in the model. If yes, it should be a NamedTuple like thisphyto = phyto_setup(params = nothing, N = [2^10, 2^10], Nsp = 2)
.abiotic
: nothing or aabiotic_setup
. Whether to include abiotic particles in the model. If yes, it should be a NamedTuple like thisabiotic = abiotic_setup(params = nothing, N = [2^10, 2^10], Nsa = 2, palat = [(:sp1, :sa1)])
.t
: Model time, start from 0 by default, in second.
Simulation
PlanktonIndividuals.Simulation.PlanktonSimulation
— MethodPlanktonSimulation(model; ΔT, iterations,
PARF = default_PARF(model.grid),
temp = default_temperature(model.grid),
diags = nothing,
vels = (;),
ΔT_vel = ΔT,
ΔT_PAR::AbstractFloat = 3600.0f0,
ΔT_temp::AbstractFloat = 3600.0f0,
output_writer = nothing,
)
Generate a PlanktonSimulation
data structure.
Keyword Arguments (Required)
ΔT
: time step in second.iterations
: run the simulation for this number of iterations.
Keyword Arguments (Optional)
PARF
: External forcings of surface PAR. Hourly PAR of a single day is provided by default.temp
: External forcings of temperature. Hourly data of a single day is provided by default.diags
: Diagnostics of the simulation generated byPlanktonDiagnostics
.vels
: The velocity fields for tracer fields and individuals.nothing
means no velocities will be applied in the simulation. Otherwise,vels
mush be aNamedTuple
containing allu
,v
, andw
. Each ofu
,v
, andw
must be an 4D-Array
of(Nx, Ny, Nz, nΔT)
elements, excluding halo points.N+1
is required for bounded direction.ΔT_vel
: time step of velocities provided externally (in seconds).ΔT_PAR
: time step of surface PAR provided externally (in seconds).ΔT_temp
: time step of temperature provided externally (in seconds).output_writer
: Output writer of the simulation generated byPlanktonOutputWriter
.
PlanktonIndividuals.Simulation.update!
— Methodupdate!(sim::PlanktonSimulation; time_offset = (vels = true, PAFR = true, temp = true))
Update the PlanktonSimulation
for sim.iterations
time steps. time_offset
is used when velocities (or PARF or temperature) starts from timestep 1, but model.t is not. It is usually used when velocity fields are too large and need to be broken down into several parts. Only one part of the whole velocity fields can be constructed into a PlanktonSimulation
, so in this PlanktonSimulation model.iteration
might no be 1, but the velocity fields need to start from 1 (same for PARF or temperature fields).
PlanktonIndividuals.Simulation.vel_copy!
— Methodvel_copy!(vel::NamedTuple, u, v, w, g::AbstractGrid{FT, TX, TY, TZ}) where {FT, TX, TY, TZ}
Copy external velocities into PlanktonModel
Output
PlanktonIndividuals.Output.PlanktonOutputWriter
— MethodPlanktonOutputWriter(;dir = "./results",
diags_prefix = "diags",
phytoplankton_prefix = "phytoplankton",
abiotic_particle_prefix = "abiotic_particle",
write_log = false,
save_diags = false,
save_phytoplankton = false,
save_abiotic_particle = false,
phytoplankton_include = (:x, :y, :z),
abiotic_particle_include = (:x, :y, :z),
phytoplankton_iteration_interval = 1,
abiotic_particle_iteration_interval = 1,
max_filesize = Inf,
)
Generate a PlanktonOutputWriter
structure which includes settings for model outputs
Keyword Arguments (Optional)
dir
: The directory to store model outputs, "./results" by defaultdiags_prefix
: Descriptive filename prefixed to diagnostic output files.phytoplankton_prefix
: Descriptive filename prefixed to phytoplankton output files.write_log
: write model logs which contain global averages of simulated phytoplankton, default:false
.save_diags
: write diagnostics to disk, default:false
.save_phytoplankton
: write phytoplankton to disk, default:false
.phytoplankton_include
: list of phytoplankton properties to save, default:(:x, :y, :z, :Sz)
.phytoplankton_iteration_interval
: The time interval that phytoplankton are saved, 1 timestep by default.save_abiotic_particle
: write abiotic_particle to disk, default:false
.abiotic_particle_include
: list of abiotic_particle properties to save, default:(:x, :y, :z, :Sz)
.abiotic_particle_iteration_interval
: The time interval that abiotic_particle are saved, 1 timestep by default.max_filesize
: The writer will stop writing to the output file once the file size exceedsmax_filesize
, and write to a new one with a consistent naming scheme ending inpart1
,part2
, etc. default:Inf
.