OceanGreensFunctionMethods

TracerInverseGaussian(Γ,Δ)

using LinearAlgebra: NumberArray

The tracer inverse Gaussian distribution with mean Γ and width Δ has probability density function

\[G(𝐱, \tau) = \sqrt{\frac{\Gamma^3 }{4 \pi \Delta^2 \tau^3 }} \exp \left( - \frac{\Gamma (\tau - \Gamma)^2}{4 \Delta ^2 \tau}\right) \]

TracerInverseGaussian()              # Tracer Inverse Gaussian distribution with unit mean and unit width, i.e. TracerInverseGaussian(1, 1)
TracerInverseGaussian(Γ, Δ)          # Tracer Inverse Gaussian distribution with mean Γ and width Δ

params(d)           # Get the parameters, i.e. (Γ, Δ)
mean(d)             # Get the mean parameter, i.e. Γ
shape(d)            # Get the shape parameter, i.e. Δ

External links

• Compare to Inverse Gaussian distribution on Wikipedia

abyssal_overturning(Ψ,model_dims)

Set volume flux, Ψ, in an abyssal overturning loop that satisfies the conservation of volume. Return a structure of Fluxes.

advective_diffusive_flux(C, Fv; ρ)

Advective-diffusive flux of tracer C given volume fluxes Fv and optional density ρ.

Arguments

• C::VectorDimArray: tracer distribution
• Fv::VectorDimArray: volume fluxes
• ρ::Number=1035kg/m^3: uniform density

Returns

• Fc::VectorDimArray: tracer flux

advective_diffusive_flux(C, Fv; ρ)

Advective-diffusive flux of tracer C given volume fluxes Fv and optional density ρ.

Arguments

• C::VectorDimArray: tracer distribution
• Fv::Fluxes: volume fluxes
• ρ::Number=1035kg/m^3: uniform density

Returns

• Fc::Fluxes: tracer flux

boundary_dimensions()

Define labels for the boundary's physical dimensions, as well as labels for the box names, consistently with the model dimensions. Use the format of DimensionalData.jl. Permits numerical quantities to be bundled with their meta-data. Dimensions are Unordered to avoid issues related to the alphabetical order.

boundary_flux(f::VectorDimArray, C::VectorDimArray, Fb::VectorDimArray)

Convergence or net effect of boundary fluxes.

Arguments

• f::VectorDimArray: Dirichlet boundary condition
• C::VectorDimArray: tracer distribution
• Fb::VectorDimArray: boundary exchange volume flux

Returns

• Jb::Fluxes: boundary tracer flux

boundary_propagator(τ, A, B; alg=:forward)

Forward and adjoint boundary propagators.

Forward boundary propagator

boundary_propagator(τ, A, B, alg=:forward)

The (forward) boundary propagator is the box model surface-to-interior transit time distribution (TTD) over transit times τ = t - t′, as given by equation 88 of Haine et al. (2024):

\[{\bf G}' (\tau) = {\bf G} (\tau) ~ {\bf B}\]

The N × Nₛ 𝐆′(τ) matrix quantifies transfer from the Nₛ components of the surface forcing to the N boxes with transit time τ.

Adjoint boundary propagator

boundary_propagator(τ, A, B, alg=:adjoint)

The box model adjoint boundary propagator (interior-to-surface TTD over transit time τ† = t″ - t, where t″ ≥ t is the time of the adjoint source; equation 93 of Haine et al., 2024) is

\[{\bf G}'^{\dagger} (\tau^{\dagger} ) = {\bf B}^{T}~ {\bf G} (\tau^{\dagger}).\]

This Nₛ × N 𝐆′†(τ†) matrix quantifies transfer from the N interior boxes to the Nₛ surface boxes with transit time τ†.

boundary_propagator_adjoint(t,A,B)

boundary_propagator_forward(t,A,B)

convergence(J)

Convergence of fluxes J of type Fluxes. This is a computational method that depends on proper slices and broadcasting and thus currently requires using parent on the left hand side below.

evolve_concentration(C₀, A, B, tlist, source_history; halflife = nothing)

Integrate forcing vector over time to compute the concentration history. Find propagator by analytical expression using eigen-methods.

Arguments

• C₀: initial tracer concentration
• A: tracer transport information used in matrix calculations
• B: boundary condition information used in matrix calculations
• tlist: list of times to save tracer concentration
• source_history::Function: returns Dirichlet boundary condition at a given time
• halflife=nothing: radioactive half life (optional)

forcing_integrand(t, tf, μ, V, B, source_history)

Integrand for boundary condition term in equation 10 (Haine et al., 2024).

Arguments

• t: time
• tf: final time
• μ: eigenvalue diagonal matrix
• V: eigenvector matrix
• B: boundary condition matrix
• source_history::Function: returns Dirichlet boundary condition at a given time

global_ttd(t, A, B; alg=:forward)

Forward and adjoint global transit time distributions (TTDs).

Forward Global TTD

The forward global (total) TTD is the sum of surface-to-interior TTDs (equation 90 of Haine et al., 2024):

\[{\cal G} (t) = {\bf G} (t) ~ {\bf B} ~ {\bf 1}_{N_S},\]

where the product with the Ns × 1 column vector of ones (i.e., last matrix in previous equation) computes the sum over surface boxes. This expression yields an N × 1 column vector that is normalized for each box.

Adjoint Global TTD

The adjoint global (total) TTD is the sum of interior-to-surface TTDs.

global_ttd_adjoint(t, A, B)

global_ttd_forward(t, A, B)

greens_function(τ,A)

Green's function for a box model (for steady transport given by the matrix 𝐀 for response at time t to a source at time t′ where τ = t - t′): the matrix exponential function of the elapsed time between the source time and field time:

\[{\bf G}(\tau) = e^{ {\bf A} \tau}\]

where 𝐆(t) is a N × N matrix with the spatial locations of field points (boxes) down its N rows and source points (boxes) along its N columns. Thus, the element 𝐆{i,j}(τ) quantifies transfer from a source at time t′ in box j to receiver at time t in box i.

ideal_age(A, B; alg= :forward)

ideal_age_adjoint(A, B)

ideal_age_forward(A, B)

integrate_forcing(t0, tf, μ, V, B, source_history)

Integrate boundary condition term in equation 10 (Haine et al., 2024).

Arguments

• t0: initial time
• tf: final time
• μ: eigenvalue diagonal matrix
• V: eigenvector matrix
• B: boundary condition matrix
• source_history::Function: returns Dirichlet boundary condition at a given time

intermediate_overturning(Ψ,model_dims)

Set the volume flux, Ψ, in an intermediate overturning loop that satisfies the conservation of volume. Return a structure of Fluxes.

linear_probe(funk, x, args...)

Probe a function to determine its linear response in matrix form. Assumes units are needed and available. A simpler function to handle cases without units would be nice.

Arguments

• funk: function to be probed
• x: input (independent) variable
• halflife::Number: radioactive half life
• args: the arguments that follow x in funk

Returns

• A::MatrixDimArray: labeled transport information used in matrix operations

location_transient_tracer_histories()

URL of iodine-129 source history file.

location_transient_tracer_histories()

URL of tracer source history file.

mass(V; ρ)

Seawater mass derived from the volume V and an optional input of density ρ.

mass_convergence(Fv)

Convergence of volume derived from a field of volume fluxes Fv, translated into a mass flux convergence with the assumption of uniform density.

maximum_timescale(μ)

Return Tmax for the eigenvalues μ. The matrix exponential of 𝐀τ has asymptotic properties because G(t) must eventually decay exponentially with timescale

\[T_{max} = -1/\mu_{min}, \]

where μmin is the eigenvalue with smallest real part. Thus, the Green's function has a maximum timescale of Tmax which is larger than all other transport timescales.

mean_age(μ, V, B; alg=:forward)

Mean age of the forward TTDs, adjoint TTDs, and residence-time distributions.

Arguments

• μ: eigenvalues vector
• V: eigenvector matrix
• B: boundary matrix
• alg=:forward: algorithm (optional)

Forward mean age

mean_age(μ, V, B, alg=:forward)

The mean transit time 𝚪 (mean age) is (equation 92 of Haine et al., 2004),

\[{\bf \Gamma} = {\bf V} ~ \mu^{-2} ~ {\bf V}^{-1} ~ {\bf B} ~ {\bf 1}_{N_S},\]

which is an N × 1 vector for each box (and which also equals the ideal age).

Adjoint mean age

mean_age(μ, V, B, alg=:adjoint)

Residence-time mean age

mean_age(μ, V, B, alg=:residence)

mean_age_adjoint(μ, V, B)

mean_age_forward(μ, V, B)

mean_age_residence(μ, V, B)

model_dimensions()

Define labels for the model's physical dimensions, as well as labels for the box names. Use the format of DimensionalData.jl. Permits numerical quantities to be bundled with their meta-data. Dimensions are Unordered to avoid issues related to the alphabetical order.

normalized_exponential_decay(t,Tmax)

path_density(μ, V, B, t, mbox, vbox)

Arguments

• μ: eigenvalue diagonal matrix
• V: eigenvector matrix
• B: boundary matrix
• t: time
• mbox: name of meridional box of interest
• vbox: name of vertical box of interest

Returns

• E: path density

The path density 𝐄_i(τ) for i ∈ 1 ... N is (equation 96 of Haine et al., 2024):

\[{\bf E}_i (\tau) = \frac{1}{N} \int_{t - \tau}^{t} {\bf G}'^{\dagger} (t^* + \tau - t) ~ {\bf D}_i ~ {\bf G}' (t - t^*) ~ d t ^* , \]

where 𝐃i is the N × N matrix unit of zeros with a single one at the i-th row and i-th column. Therefore,

\[{\bf E}_i (\tau) = \frac{1}{N} \int_{0}^{\tau} {\bf G}'^{\dagger} (t') ~ {\bf D}_i ~ {\bf G}' (\tau - t') ~ d t '\]

and

\[{\bf E}_i (\tau) = \frac{1}{N} {\bf B}^{T} \int_{0}^{\tau} ~ e^{{\bf A} t'} ~ {\bf D}_i e^{{\bf A} ( au - t')} ~ d t ' {\bf B}\]

and

\[{\bf E}_i (\tau) = \frac{1}{N}{\bf B}^{T} ~ {\bf V} \left( \overline{\bf D}_i \circ \Phi (t) \right) {\bf V}^{-1} ~ {\bf B}\]

where ϕ is defined in equation (100) of Haine et al. (2024). For a particular interior box i, 𝐄_i(τ) is the density of pathways between all combinations of surface entry and surface exit boxes over total residence time τ.

phi_function(t, μ)

radioactive_decay(C, halflife)

Radioactive decay rate of tracer C with half life of halflife.

read_iodine129_history()

read_transient_tracer_histories()

Read transient tracer source histories and save as a DimArray.

residence_time(t, A, B)

The surface-to-surface residence-time distribution (RTD) is (equations 94 and 95 of Haine et al., 2024):

\[{\bf R} (\tau) = \frac{1}{N} \int_{t - \tau}^{t} {\bf G}'^{\dagger} (t^* + \tau - t) ~ {\bf G}' (t - t^*) ~ d t ^*\]

or

\[{\bf R} (\tau) = \frac{\tau}{N} {\bf B}^{T} ~ {\bf G} (\tau) ~ {\bf B},\]

where N is the number of boxes, G(τ) is the forward Green's function and B is the boundary matrix. The Ns × Ns R(τ) matrix quantifies transfer from the Ns surface boxes back to the Ns surface boxes with residence time τ (element R{i,j}(τ) quantifies transfer from entry box j to exit box i).

Note: not normalized by number of boxes in this code: consistent with manuscript?

steady_tracer_timeseries(tracername, A, B, halflife, tlist, mbox1, vbox1)

Simulate non-transient tracers and return tracer timeseries from one box.

Arguments

• tracername: name of tracer
• A: tracer transport matrix
• B: boundary condition matrix
• halflife: radioactive half life
• BD: Dirichlet boundary condition
• tlist: list of times to save tracer concentration
• mbox: name of meridional box of interest
• vbox: name of vertical box of interest

timestep_initial_condition(C, μ, V, ti, tf)

Arguments

• C::VectorDimArray: tracer distribution at ti
• μ: eigenvalue diagonal matrix
• V: eigenvector matrix
• ti: initial time
• tf: final time

Returns

• Cf::VectorDimArray: tracer distribution at tf

tracer_point_source_history(tracername, BD)

Return a function that yields transient tracer source history (such as CFCs) at given time.

Arguments

• tracername: name of tracer in source history file
• BD::DimArray: Dirichlet boundary condition compendium for many tracers

tracer_source_history(t, tracername, box2_box1_ratio, BD = nothing)

Return source history values for all boundary points.

Arguments

• t: time
• tracername: name of tracer in source history file
• box2_box1_ratio: ratio of boundary condition value in Mid-latitudes to High Latitudes
• BD::DimArray=nothing: Dirichlet boundary condition compendium (optional)

tracer_tendency(f, C, Fv, Fb, V)

Tracer tendency ∂C/∂t for a boundary flux f, for use with finding B boundary matrix.

Arguments

• f::VectorDimArray: Dirichlet boundary condition
• C::VectorDimArray: tracer distribution
• Fv::Fluxes: volume fluxes
• Fb::Fluxes: volume fluxes
• V::VectorDimArray: box volume

Returns

• dCdt::VectorDimArray: tracer tendency

tracer_tendency(C)

Tracer tendency ∂C/∂t for the radioactive decay of a tracer C with half life halflife, for use with finding the radioactive contribution to a tracer transport matrix.

Arguments

• C::VectorDimArray: tracer distribution
• halflife::Number: radioactive half life

Returns

• dCdt::VectorDimArray: tracer tendency

tracer_tendency(C, f, Fv, Fb, V)

Tracer tendency ∂C/∂t for a tracer C, especially useful for finding a tracer transport matrix.

Arguments

• C::VectorDimArray: tracer distribution
• f::VectorDimArray: Dirichlet boundary condition
• Fv::Fluxes: volume fluxes
• Fb::VectorDimArray: boundary flux convergence
• V::VectorDimArray: box volume

Returns

• dCdt::VectorDimArray: tracer tendency

tracer_timeseries(tracername, A, B, tlist, mbox1, vbox1; BD=nothing, halflife=nothing)

Simulate tracers and return tracer timeseries from one box.

Arguments

• tracername: name of tracer
• A: tracer transport matrix
• B: boundary condition matrix
• tlist: list of times to save tracer concentration
• mbox: name of meridional box of interest
• vbox: name of vertical box of interest
• BD=nothing: Dirichlet boundary condition
• halflife=nothing: radioactive half life

transient_tracer_timeseries(tracername, A, B, BD, tlist, mbox1, vbox1; halflife = nothing)

Simulate transient tracers and return tracer timeseries from one box.

Arguments

• tracername: name of tracer
• A: tracer transport matrix
• B: boundary condition matrix
• BD: Dirichlet boundary condition
• tlist: list of times to save tracer concentration
• mbox: name of meridional box of interest
• vbox: name of vertical box of interest
• halflife=nothing: radioactive half life

ttd_width(μ, V, B; alg=:forward)

Width of the forward TTDs, adjoint TTDs, and residence-time distributions.

Arguments

• μ: eigenvalue diagonal matrix
• V: eigenvector matrix
• B: boundary matrix
• alg=:forward: algorithm (optional)

Returns

• Δ: TTD width

Width of forward TTD

ttd_width(μ, V, B, alg=:forward)

The TTD width is given by (equation 92 of Haine et al., 2024),

\[2 {\bf \Delta}^2 = -2 ~ {\bf V} ~ \mu^{-3} ~ {\bf V}^{-1} ~ {\bf B} ~ {\bf 1}_{N_S} - {\bf \Gamma}^2,\]

which is a N × 1 vector for each box.

Adjoint mean age

mean_age(μ, V, B, alg=:adjoint)

Residence-time mean age

mean_age(μ, V, B, alg=:residence)

ttd_width_adjoint(μ, V, B)

ttd_width_forward(μ, V, B)

ttd_width_residence(μ, V, B)

vertical_diffusion(Fv_exchange,model_dims)

Set vertical diffusive-like exchange flux Fv_exchange. Return a structure of Fluxes.

watermass_fraction(μ, V, B; alg=:forward)

Forward, adjoint, and residence-time water-mass fractions.

Forward water-mass fraction

watermass_fraction(μ, V, B, alg=:forward)

The water mass fractions are (equation 89 of Haine et al., 2024)

\[{\bf a} = \int_0^{\infty} {\bf G} (\tau) ~ {\bf B} ~ d \tau\]

or

\[{\bf a} = -{\bf V} ~ \mu^{-1} ~ {\bf V}^{-1} ~ {\bf B} , \]

which is an N × Ns matrix with the interior boxes down the rows and the surface sources across the rows.

Adjoint water-mass fraction

watermass_fraction(μ, V, B, alg=:adjoint)

Fraction of water that will return to the surface in a particular box.

Residence-time water-mass fraction

watermass_fraction(μ, V, B, alg=:residence)

Fraction of water that leaves a particular box and returns in another box.

watermass_fraction_adjoint(μ, V, B)

watermass_fraction_forward(μ, V, B)

watermass_fraction_residence(μ, V, B)