From 67536a9c8dd1248f7ff042c3623becb37baf3e81 Mon Sep 17 00:00:00 2001 From: Jago Stong-Wright Date: Sun, 24 Sep 2023 19:16:31 +0100 Subject: [PATCH] trying to not modify the sediment model --- src/Boundaries/Sediments/simple_multi_G.jl | 130 ++++++++++++--------- 1 file changed, 74 insertions(+), 56 deletions(-) diff --git a/src/Boundaries/Sediments/simple_multi_G.jl b/src/Boundaries/Sediments/simple_multi_G.jl index 9d629222a..77e804910 100644 --- a/src/Boundaries/Sediments/simple_multi_G.jl +++ b/src/Boundaries/Sediments/simple_multi_G.jl @@ -29,17 +29,18 @@ end """ SimpleMultiG(; grid - fast_decay_rate = 2/day, - slow_decay_rate = 0.2/day, - fast_redfield = 0.1509, - slow_redfield = 0.13, - fast_fraction = 0.74, - slow_fraction = 0.26, - refactory_fraction = 0.1, - nitrate_oxidation_params = (- 1.9785, 0.2261, -0.0615, -0.0289, - 0.36109, - 0.0232), - denitrification_params = (- 3.0790, 1.7509, 0.0593, - 0.1923, 0.0604, 0.0662), - anoxic_params = (- 3.9476, 2.6269, - 0.2426, -1.3349, 0.1826, - 0.0143), - solid_dep_params = (0.233, 0.336, 982.0, - 1.548)) + fast_decay_rate::FT = 2/day, + slow_decay_rate::FT = 0.2/day, + fast_redfield::FT = 0.1509, + slow_redfield::FT = 0.13, + fast_fraction::FT = 0.74, + slow_fraction::FT = 0.26, + refactory_fraction::FT = 0.1, + nitrate_oxidation_params::P1 = (A = - 1.9785, B = 0.2261, C = -0.0615, D = -0.0289, E = - 0.36109, F = - 0.0232), + denitrification_params::P2 = (A = - 3.0790, B = 1.7509, C = 0.0593, D = - 0.1923, E = 0.0604, F = 0.0662), + anoxic_params::P3 = (A = - 3.9476, B = 2.6269, C = - 0.2426, D = -1.3349, E = 0.1826, F = - 0.0143), + depth = abs(znodes(grid, Face())[1]), + solid_dep_params::P4 = (A = 0.233, B = 0.336, C = 982, D = - 1.548, depth = depth)) Return a single-layer "multi G" sediment model (`SimpleMultiG`) on `grid`, where parameters can be optionally specified. @@ -70,17 +71,37 @@ Single-layer multi-G sediment model (Float64) ``` """ function SimpleMultiG(; grid, - fast_decay_rate = 2/day, - slow_decay_rate = 0.2/day, - fast_redfield = 0.1509, - slow_redfield = 0.13, - fast_fraction = 0.74, - slow_fraction = 0.26, - refactory_fraction = 0.1, - nitrate_oxidation_params = (- 1.9785, 0.2261, -0.0615, -0.0289, - 0.36109, - 0.0232), - denitrification_params = (- 3.0790, 1.7509, 0.0593, - 0.1923, 0.0604, 0.0662), - anoxic_params = (- 3.9476, 2.6269, - 0.2426, -1.3349, 0.1826, - 0.0143), - solid_dep_params = (0.233, 0.336, 982.0, - 1.548)) + fast_decay_rate::FT = 2/day, + slow_decay_rate::FT = 0.2/day, + fast_redfield::FT = 0.1509, + slow_redfield::FT = 0.13, + fast_fraction::FT = 0.74, + slow_fraction::FT = 0.26, + refactory_fraction::FT = 0.1, + nitrate_oxidation_params::P1 = (A = - 1.9785, + B = 0.2261, + C = -0.0615, + D = -0.0289, + E = - 0.36109, + F = - 0.0232), + denitrification_params::P2 = (A = - 3.0790, + B = 1.7509, + C = 0.0593, + D = - 0.1923, + E = 0.0604, + F = 0.0662), + anoxic_params::P3 = (A = - 3.9476, + B = 2.6269, + C = - 0.2426, + D = -1.3349, + E = 0.1826, + F = - 0.0143), + depth = abs(znodes(grid, Face())[1]), + solid_dep_params::P4 = (A = 0.233, + B = 0.336, + C = 982, + D = - 1.548, + depth = depth)) where {FT, P1, P2, P3, P4} @warn "Sediment models are an experimental feature and have not yet been validated." @info "This sediment model is currently only compatible with models providing NH₄, NO₃, O₂, and DIC." @@ -91,7 +112,7 @@ function SimpleMultiG(; grid, fields = NamedTuple{tracer_names}(Tuple(CenterField(grid) for tracer in tracer_names)) tendencies = (Gⁿ = NamedTuple{tracer_names}(Tuple(CenterField(grid) for tracer in tracer_names)), G⁻ = NamedTuple{tracer_names}(Tuple(CenterField(grid) for tracer in tracer_names))) - + bottom_indices = calculate_bottom_indices(grid) F = typeof(fields) @@ -111,20 +132,19 @@ function SimpleMultiG(; grid, end adapt_structure(to, sediment::SimpleMultiG) = - SimpleMultiG(adapt(to, sediment.fast_decay_rate), - adapt(to, sediment.slow_decay_rate), - adapt(to, sediment.fast_redfield), - adapt(to, sediment.slow_redfield), - adapt(to, sediment.fast_fraction), - adapt(to, sediment.slow_fraction), - adapt(to, sediment.refactory_fraction), - adapt(to, sediment.nitrate_oxidation_params), - adapt(to, sediment.denitrification_params), - adapt(to, sediment.anoxic_params), - adapt(to, sediment.solid_dep_params), + SimpleMultiG(sediment.fast_decay_rate, + sediment.slow_decay_rate, + sediment.fast_redfield, + sediment.slow_redfield, + sediment.fast_fraction, + sediment.slow_fraction, + sediment.refactory_fraction, + sediment.nitrate_oxidation_params, + sediment.denitrification_params, + sediment.anoxic_params, + sediment.solid_dep_params, adapt(to, sediment.fields), - nothing, - adapt(to, sediment.bottom_indices)) + adapt(to, sediment.tendencies)) sediment_tracers(::SimpleMultiG) = (:C_slow, :C_fast, :C_ref, :N_slow, :N_fast, :N_ref) sediment_fields(model::SimpleMultiG) = (C_slow = model.fields.C_slow, @@ -136,7 +156,7 @@ sediment_fields(model::SimpleMultiG) = (C_slow = model.fields.C_slow, @inline bottom_index_array(sediment::SimpleMultiG) = sediment.bottom_indices -@kernel function _calculate_tendencies!(sediment::SimpleMultiG, bgc, grid, advection, tracers, tendencies, sediment_tendencies) +@kernel function _calculate_tendencies!(sediment::SimpleMultiG, bgc, grid, advection, tracers, timestepper) i, j = @index(Global, NTuple) k = bottom_index(i, j, sediment) @@ -159,30 +179,28 @@ sediment_fields(model::SimpleMultiG) = (C_slow = model.fields.C_slow, Cᵐⁱⁿ = C_min_slow + C_min_fast Nᵐⁱⁿ = N_min_slow + N_min_fast - reactivity = Cᵐⁱⁿ * day / (sediment.fields.C_slow[i, j, 1] + sediment.fields.C_fast[i, j, 1]) + k = Cᵐⁱⁿ * day / (sediment.fields.C_slow[i, j, 1] + sediment.fields.C_fast[i, j, 1]) # sediment evolution - sediment_tendencies.C_slow[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.slow_fraction * carbon_deposition - C_min_slow - sediment_tendencies.C_fast[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.fast_fraction * carbon_deposition - C_min_fast - sediment_tendencies.C_ref[i, j, 1] = sediment.refactory_fraction * carbon_deposition + sediment.tendencies.Gⁿ.C_slow[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.slow_fraction * carbon_deposition - C_min_slow + sediment.tendencies.Gⁿ.C_fast[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.fast_fraction * carbon_deposition - C_min_fast + sediment.tendencies.Gⁿ.C_ref[i, j, 1] = sediment.refactory_fraction * carbon_deposition - sediment_tendencies.N_slow[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.slow_fraction * nitrogen_deposition - N_min_slow - sediment_tendencies.N_fast[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.fast_fraction * nitrogen_deposition - N_min_fast - sediment_tendencies.N_ref[i, j, 1] = sediment.refactory_fraction * nitrogen_deposition + sediment.tendencies.Gⁿ.N_slow[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.slow_fraction * nitrogen_deposition - N_min_slow + sediment.tendencies.Gⁿ.N_fast[i, j, 1] = (1 - sediment.refactory_fraction) * sediment.fast_fraction * nitrogen_deposition - N_min_fast + sediment.tendencies.Gⁿ.N_ref[i, j, 1] = sediment.refactory_fraction * nitrogen_deposition # efflux/influx - O₂ = tracers.O₂[i, j, k] - NO₃ = tracers.NO₃[i, j, k] - NH₄ = tracers.NH₄[i, j, k] - - A, B, C, D, E, F = sediment.nitrate_oxidation_params - - pₙᵢₜ = exp(A + - B * log(Cᵐⁱⁿ * day) * log(O₂) + - C * log(Cᵐⁱⁿ * day) ^ 2 + - D * log(reactivity) * log(NH₄) + - E * log(Cᵐⁱⁿ * day) + - F * log(Cᵐⁱⁿ * day) * log(NH₄)) / (Nᵐⁱⁿ * day) + O₂ = tracers.O₂[i, j, 1] + NO₃ = tracers.NO₃[i, j, 1] + NH₄ = tracers.NH₄[i, j, 1] + + pₙᵢₜ = exp(sediment.nitrate_oxidation_params.A + + sediment.nitrate_oxidation_params.B * log(Cᵐⁱⁿ * day) * log(O₂) + + sediment.nitrate_oxidation_params.C * log(Cᵐⁱⁿ * day) ^ 2 + + sediment.nitrate_oxidation_params.D * log(k) * log(NH₄) + + sediment.nitrate_oxidation_params.E * log(Cᵐⁱⁿ * day) + + sediment.nitrate_oxidation_params.F * log(Cᵐⁱⁿ * day) * log(NH₄)) / (Nᵐⁱⁿ * day) #= pᵈᵉⁿⁱᵗ = exp(sediment.denitrification_params.A +