Fixed error in SLatissima and Particles setup

src/Light/2band.jl

@@ -13,6 +13,8 @@ using Oceananigans.Utils: launch!
 
     ∫chlᵉʳ = -(P[i, j, grid.Nz]*params.Rd_chl/params.r_pig)^params.e_r*z
     ∫chlᵉᵇ = -(P[i, j, grid.Nz]*params.Rd_chl/params.r_pig)^params.e_b*z
+
+    # TODO: Change this to an interpolation to the midpoint between the cells/make this actually valid for all grid types
     PAR[i, j, grid.Nz] =  sp*(exp(params.k_r0 * z - params.Χ_rp * ∫chlᵉʳ).+ exp(params.k_b0 * z - params.Χ_bp * ∫chlᵉᵇ))/2
     for k=grid.Nz-1:-1:1
         z = grid.zᵃᵃᶜ[k]

src/Models/Individuals/SLatissima.jl

@@ -105,8 +105,8 @@ function equations(x::AbstractFloat, y::AbstractFloat, z::AbstractFloat, t::Abst
         r = _r(T, μ, j_NO₃ + j_NH₄, params)
 
         dA = (μ - ν) * A / (60*60*24)
-        dN = ((j_NO₃ + j_NH₄ - p*e*14/(12*6.56)) / params.K_A - μ * (N + params.N_struct)) / (60*60*24)
-        dC = ((p* (1 - e) - r) / params.K_A - μ * (C + params.C_struct)) / (60*60*24)
+        dN = ((j_NO₃ + j_NH₄ - p * e * 14/(12*6.56)) / params.K_A - μ * (N + params.N_struct)) / (60*60*24)
+        dC = ((p * (1 - e) - r) / params.K_A - μ * (C + params.C_struct)) / (60*60*24)
 
         A_new = A+dA*Δt 
         N_new = N+dN*Δt 
@@ -117,7 +117,7 @@ function equations(x::AbstractFloat, y::AbstractFloat, z::AbstractFloat, t::Abst
             C_new = params.C_min
         end
 
-        pp = (p-r)*A / (60*60*24*12*0.001) #gC/dm^2/hr to mmol C/s
+        pp = (p - r)*A / (60*60*24*12*0.001) #gC/dm^2/hr to mmol C/s
         e *= p*A / (60*60*24*12*0.001)#mmol C/s
         νⁿ = ν*params.K_A*A*(params.N_struct + N) / (60*60*24*14*0.001)#1/hr to mmol N/s
         νᶜ = ν*params.K_A*A*(params.C_struct + C) / (60*60*24*12*0.001)#1/hr to mmol C/s
@@ -249,6 +249,8 @@ struct Particle
     NO₃  :: AbstractFloat
     NH₄  :: AbstractFloat
     PAR :: AbstractFloat
+    T :: AbstractFloat
+    S :: AbstractFloat
 end
 
 default_tracers = (NO₃ = :NO₃, PAR=:PAR) # tracer = property
@@ -269,7 +271,7 @@ function defineparticles(initials, n)
             throw(ArgumentError("Invalid initial values given for $var, must be a single number or vector of length n"))
         end
     end
-    return StructArray{Particle}((x̄₀[1], x̄₀[2], x̄₀[3], zeros(n), zeros(n), zeros(n), x̄₀[4], x̄₀[5], x̄₀[6], [zeros(n) for i in 1:10]...))
+    return StructArray{Particle}((x̄₀[1], x̄₀[2], x̄₀[3], zeros(n), zeros(n), zeros(n), x̄₀[4], x̄₀[5], x̄₀[6], [zeros(n) for i in 1:12]...))
 end
 
 @inline no_dynamics(args...) = nothing
@@ -339,7 +341,7 @@ function setup(n, x₀, y₀, z₀, A₀, N₀, C₀, latitude;
         if tracer in optional_tracer_dependencies
             tracers = merge(tracers, NamedTuple{(tracer, )}((getproperty(optional_tracer_dependencies, tracer), )))
         end
-        if tracer in optional_tracer_coupling
+        if tracer in keys(optional_tracer_coupling)
             coupled = merge(coupled, NamedTuple{(tracer, )}((getproperty(optional_tracer_coupling, tracer), )))
         end
     end

src/Particles/Particles.jl

@@ -132,7 +132,7 @@ function ActiveLagrangianParticles(particles::StructArray;
                                                         scalefactor::AbstractFloat = 1.0, 
                                                         custom_dynamics=no_dynamics)
 
-    required_fields = (equation_arguments..., prognostic..., diagnostic..., keys(tracked_fields)..., keys(coupled_fields)...)
+    required_fields = (equation_arguments..., prognostic..., diagnostic..., keys(tracked_fields)..., values(coupled_fields)...)
     for property in required_fields
         if !(property in propertynames(particles))
             throw(ArgumentError("$property is a required particle property but $(eltype(particles)) has no property $property."))

test/runtests.jl

@@ -3,3 +3,5 @@ include("test_light.jl")
 include("test_simulations.jl")
 include("test_auxiliary_field_forcing.jl")
 include("test_LOBSTER.jl")
+include("test_gasexchange.jl")
+include("test_slatissima.jl")

test/test_slatissima.jl

@@ -0,0 +1,39 @@
+using Test
+using OceanBioME: SLatissima
+
+@testset "Sugar Kelp Particle Setup" begin
+    # Initial properties
+
+    kelp_particles = SLatissima.setup(1, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0)
+
+    @test length(kelp_particles) == 1
+    @test all(kelp_particles.properties.x .== 1.0)
+
+    kelp_particles = SLatissima.setup(2, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0)
+
+    @test length(kelp_particles) == 2
+    @test all(kelp_particles.properties.x .== 1.0)
+
+    kelp_particles = SLatissima.setup(2, [1.0, 2.0], 1.0, 1.0, 1.0, 1.0, 1.0, 1.0)
+
+    @test length(kelp_particles) == 2
+    @test all(kelp_particles.properties.x == [1.0, 2.0])
+    @test all(kelp_particles.properties.y .== 1.0)
+
+    # Temperature and salinity function passing
+
+    @test !(:T in keys(kelp_particles.parameters.equation_parameters)) && !(:S in keys(kelp_particles.parameters.equation_parameters))
+
+    kelp_particles = SLatissima.setup(1, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0; T = (x, y, z, t) -> 1.0, S = (x, y, z, t) -> 35.0, urel = 0.2)
+    @test (:T in keys(kelp_particles.parameters.equation_parameters)) && (:S in keys(kelp_particles.parameters.equation_parameters))
+
+    # Optional tracers
+
+    @test kelp_particles.parameters.coupled_fields == (NO₃ = :j_NO₃,)
+    @test kelp_particles.parameters.tracked_fields == (NO₃ = :NO₃, PAR = :PAR)
+
+    kelp_particles = SLatissima.setup(1, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0; optional_tracers = (:NH₄, :DIC, :DD, :DDᶜ, :OXY, :DOM))
+
+    @test kelp_particles.parameters.tracked_fields == (NO₃ = :NO₃, PAR = :PAR, NH₄ = :NH₄)
+    @test kelp_particles.parameters.coupled_fields == (NO₃ = :j_NO₃, NH₄ = :j_NH₄, DIC = :j_DIC, DD = :νⁿ, DDᶜ = :νᶜ, OXY = :j_OXY, DOM = :e)
+end