From 279bf2887679d616d33a7bbfbc4814f5c4088fdb Mon Sep 17 00:00:00 2001 From: jbisits Date: Wed, 20 Nov 2024 14:32:45 +1100 Subject: [PATCH 01/12] Move background related things to its own scrtip --- src/StaircaseShenanigans.jl | 1 + src/set_staircase_initial_conditions.jl | 4 +- src/staircase_background.jl | 137 +++++++++++++++++++++++ src/staircase_restoring.jl | 141 +----------------------- 4 files changed, 142 insertions(+), 141 deletions(-) create mode 100644 src/staircase_background.jl diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index 89aaf02..adb315a 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -51,6 +51,7 @@ include("staircase_noise.jl") include("staircase_diagnostics.jl") include("staircase_model.jl") include("set_staircase_initial_conditions.jl") +include("staircase_background.jl") include("staircase_restoring.jl") include("alt_lineareos.jl") include("makie_plotting_functions.jl") diff --git a/src/set_staircase_initial_conditions.jl b/src/set_staircase_initial_conditions.jl index d1b5cf0..3299a85 100644 --- a/src/set_staircase_initial_conditions.jl +++ b/src/set_staircase_initial_conditions.jl @@ -51,8 +51,8 @@ function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs) return nothing end -"Here the `BackgroundField` behaves as the `initial_condition` and noise is added to the -tracer fields to create an instability." +"Currently only random noise can be added as the initial condition in the tracer field +when there is a `BackgroundField`." set_staircase_initial_conditions!(model, ics::PeriodoicSingleInterfaceICs) = nothing "Fallback --- don't set any noise." diff --git a/src/staircase_background.jl b/src/staircase_background.jl new file mode 100644 index 0000000..5ae4244 --- /dev/null +++ b/src/staircase_background.jl @@ -0,0 +1,137 @@ +Oceananigans.BackgroundField(bf::BackgroundFunction) = + BackgroundField(bf.func, parameters = bf.parameters) +""" + mutable struct BackgroundTanh{F, T, P} +Container for a tanh background field. +""" +mutable struct BackgroundTanh{F, T} <: BackgroundFunction + "tanh function" + func :: F + "Scale the steepness of the `tanh` change" + scale :: T + "Parameters for the tanh background field" + parameters :: NamedTuple +end +BackgroundTanh() = BackgroundTanh(tanh_background, 100, NamedTuple()) +BackgroundTanh(scale) = BackgroundTanh(tanh_background, scale, NamedTuple()) +""" + mutable struct BackgroundLinear{F, P} +Container for a linear background field. +""" +mutable struct BackgroundLinear{F} <: BackgroundFunction + "Linear function" + func :: F + "Parameters for the tanh background field" + parameters :: NamedTuple +end +BackgroundLinear() = BackgroundLinear(linear_background, NamedTuple()) + +function Base.show(io::IO, bf::BackgroundFunction) + if bf isa BackgroundTanh + println(io, "BackgroundTanh") + println(io, "┣━━━ function: $(bf.func)") + println(io, "┣━━━━━━ scale: $(bf.scale)") + print(io, "┗━ parameters: $(bf.parameters)") + elseif bf isa BackgroundLinear + println(io, "BackgroundLinear") + println(io, "┣━━━ function: $(bf.func)") + print(io, "┗━ parameters: $(bf.parameters)") + end +end +Base.summary(bt::BackgroundTanh) = "$(bt.func)" +Base.summary(bl::BackgroundLinear) = "$(bl.func)" + +""" + function S_and_T_background_fields(initial_conditions) +Set background fields for the `S` and `T` tracer fields where the domain is triply periodic. +""" +function S_and_T_background_fields(ics::PeriodicSTSingleInterfaceInitialConditions, Lz) + + get_parameters!(ics, :salinity_values, Lz) + S_background = BackgroundField(ics.background_state) + get_parameters!(ics, :temperature_values, Lz) + T_background = BackgroundField(ics.background_state) + + return (S = S_background, T = T_background) +end +"Sets a background state that is hyperbolic tangent. There is also a method to save an +`Array` of this backgorund state to output." +@inline tanh_background(x, y, z, t, p) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.Lz)) +tanh_background(z, ΔC, Cₗ, Lz, z_interface, D) = Cₗ - 0.5 * ΔC * (1 + tanh(D * (z - z_interface) / Lz)) +@inline linear_background(x, y, z, t, p) = p.Cᵤ - p.ΔC * z / p.Lz +linear_background(z, ΔC, Cᵤ, Lz) = Cᵤ - ΔC * z / Lz + +function get_parameters!(ics::PeriodicSTSingleInterfaceInitialConditions, tracer::Symbol, Lz) + + z_interface = ics.depth_of_interface + C = Array(getproperty(ics, tracer)) + ΔC = diff(C)[1] + Cᵤ, Cₗ = C + update_parameters!(ics.background_state, ΔC, Cᵤ, Cₗ, abs(Lz), z_interface) + + return nothing +end + +function update_parameters!(backgound_state::BackgroundTanh, ΔC, Cᵤ, Cₗ, Lz, z_interface) + + backgound_state.parameters = (; ΔC, Cᵤ, Cₗ, Lz, z_interface, D = backgound_state.scale) + return nothing +end +function update_parameters!(backgound_state::BackgroundLinear, ΔC, Cᵤ, Cₗ, Lz, z_interface) + + backgound_state.parameters = (; ΔC, Cᵤ, Lz) + return nothing +end + +""" + function save_background_state!(simulation, sdns) +Where there is `BackgroundField` (currently this assumes that is for periodic simulations) +save the background state for the tracers and density so this can be used later +""" +save_background_state!(simulation, sdns) = save_background_state!(simulation, sdns.model, sdns.initial_conditions) +save_background_state!(simulation, model, initial_conditions) = nothing +function save_background_state!(simulation, model, initial_conditions::PeriodoicSingleInterfaceICs) + + S_background = Field(model.background_fields.tracers.S) + compute!(S_background) + S_background_array = Array(interior(S_background, :, :, :)) + T_background = Field(model.background_fields.tracers.T) + compute!(T_background) + T_background_array = Array(interior(T_background, :, :, :)) + σ_background = Field(seawater_density(model, temperature = T_background, salinity = S_background, + geopotential_height = 0)) + compute!(σ_background) + σ_background_array = Array(interior(σ_background, :, :, :)) + + if simulation.output_writers[:tracers] isa NetCDFOutputWriter + + NCDataset(simulation.output_writers[:tracers].filepath, "a") do ds + defVar(ds, "S_background", S_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for salinity", + "units" => "gkg⁻¹")) + defVar(ds, "T_background", T_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for temperature", + "units" => "°C")) + end + + NCDataset(simulation.output_writers[:computed_output].filepath, "a") do ds + defVar(ds, "σ_background", σ_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for potential density (0dbar) computed from the `S` and `T` background fields", + "units" => "kgm⁻³")) + end + + elseif simulation.output_writers[:tracers] isa JLD2OutputWriter + + jldopen(simulation.output_writers[:tracers].filepath, "a+") do f + f["S_background"] = S_background_array + f["T_background"] = T_background_array + end + + jldopen(simulation.output_writers[:computed_output].filepath, "a+") do f + f["σ_background"] = σ_background_array + end + + end + + return nothing +end diff --git a/src/staircase_restoring.jl b/src/staircase_restoring.jl index c61f07b..e1a7517 100644 --- a/src/staircase_restoring.jl +++ b/src/staircase_restoring.jl @@ -1,144 +1,7 @@ -Oceananigans.BackgroundField(bf::BackgroundFunction) = - BackgroundField(bf.func, parameters = bf.parameters) -""" - mutable struct BackgroundTanh{F, T, P} -Container for a tanh background field. -""" -mutable struct BackgroundTanh{F, T} <: BackgroundFunction - "tanh function" - func :: F - "Scale the steepness of the `tanh` change" - scale :: T - "Parameters for the tanh background field" - parameters :: NamedTuple -end -BackgroundTanh() = BackgroundTanh(tanh_background, 100, NamedTuple()) -BackgroundTanh(scale) = BackgroundTanh(tanh_background, scale, NamedTuple()) -""" - mutable struct BackgroundLinear{F, P} -Container for a linear background field. -""" -mutable struct BackgroundLinear{F} <: BackgroundFunction - "Linear function" - func :: F - "Parameters for the tanh background field" - parameters :: NamedTuple -end -BackgroundLinear() = BackgroundLinear(linear_background, NamedTuple()) - -function Base.show(io::IO, bf::BackgroundFunction) - if bf isa BackgroundTanh - println(io, "BackgroundTanh") - println(io, "┣━━━ function: $(bf.func)") - println(io, "┣━━━━━━ scale: $(bf.scale)") - print(io, "┗━ parameters: $(bf.parameters)") - elseif bf isa BackgroundLinear - println(io, "BackgroundLinear") - println(io, "┣━━━ function: $(bf.func)") - print(io, "┗━ parameters: $(bf.parameters)") - end -end -Base.summary(bt::BackgroundTanh) = "$(bt.func)" -Base.summary(bl::BackgroundLinear) = "$(bl.func)" - -""" - function S_and_T_background_fields(initial_conditions) -Set background fields for the `S` and `T` tracer fields where the domain is triply periodic. -""" -function S_and_T_background_fields(ics::PeriodicSTSingleInterfaceInitialConditions, Lz) - - get_parameters!(ics, :salinity_values, Lz) - S_background = BackgroundField(ics.background_state) - get_parameters!(ics, :temperature_values, Lz) - T_background = BackgroundField(ics.background_state) - - return (S = S_background, T = T_background) -end -"Sets a background state that is hyperbolic tangent. There is also a method to save an -`Array` of this backgorund state to output." -@inline tanh_background(x, y, z, t, p) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.Lz)) -tanh_background(z, ΔC, Cₗ, Lz, z_interface, D) = Cₗ - 0.5 * ΔC * (1 + tanh(D * (z - z_interface) / Lz)) -@inline linear_background(x, y, z, t, p) = p.Cᵤ - p.ΔC * z / p.Lz -linear_background(z, ΔC, Cᵤ, Lz) = Cᵤ - ΔC * z / Lz - -function get_parameters!(ics::PeriodicSTSingleInterfaceInitialConditions, tracer::Symbol, Lz) - - z_interface = ics.depth_of_interface - C = Array(getproperty(ics, tracer)) - ΔC = diff(C)[1] - Cᵤ, Cₗ = C - update_parameters!(ics.background_state, ΔC, Cᵤ, Cₗ, abs(Lz), z_interface) - - return nothing -end - -function update_parameters!(backgound_state::BackgroundTanh, ΔC, Cᵤ, Cₗ, Lz, z_interface) - - backgound_state.parameters = (; ΔC, Cᵤ, Cₗ, Lz, z_interface, D = backgound_state.scale) - return nothing -end -function update_parameters!(backgound_state::BackgroundLinear, ΔC, Cᵤ, Cₗ, Lz, z_interface) - - backgound_state.parameters = (; ΔC, Cᵤ, Lz) - return nothing -end - -""" - function save_background_state!(simulation, sdns) -Where there is `BackgroundField` (currently this assumes that is for periodic simulations) -save the background state for the tracers and density so this can be used later -""" -save_background_state!(simulation, sdns) = save_background_state!(simulation, sdns.model, sdns.initial_conditions) -save_background_state!(simulation, model, initial_conditions) = nothing -function save_background_state!(simulation, model, initial_conditions::PeriodoicSingleInterfaceICs) - - S_background = Field(model.background_fields.tracers.S) - compute!(S_background) - S_background_array = Array(interior(S_background, :, :, :)) - T_background = Field(model.background_fields.tracers.T) - compute!(T_background) - T_background_array = Array(interior(T_background, :, :, :)) - σ_background = Field(seawater_density(model, temperature = T_background, salinity = S_background, - geopotential_height = 0)) - compute!(σ_background) - σ_background_array = Array(interior(σ_background, :, :, :)) - - if simulation.output_writers[:tracers] isa NetCDFOutputWriter - - NCDataset(simulation.output_writers[:tracers].filepath, "a") do ds - defVar(ds, "S_background", S_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for salinity", - "units" => "gkg⁻¹")) - defVar(ds, "T_background", T_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for temperature", - "units" => "°C")) - end - - NCDataset(simulation.output_writers[:computed_output].filepath, "a") do ds - defVar(ds, "σ_background", σ_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for potential density (0dbar) computed from the `S` and `T` background fields", - "units" => "kgm⁻³")) - end - - elseif simulation.output_writers[:tracers] isa JLD2OutputWriter - - jldopen(simulation.output_writers[:tracers].filepath, "a+") do f - f["S_background"] = S_background_array - f["T_background"] = T_background_array - end - - jldopen(simulation.output_writers[:computed_output].filepath, "a+") do f - f["σ_background"] = σ_background_array - end - - end - - return nothing -end - # The following functions are to be used as `BoundaryConditions` so that tracers can # re-enter the domain with the initial gradient added effectively allowing the gradient to -# be maintained. Another option is to add background tracer fields and only evolve the anomaly +# be maintained. This is effectively a restoring (or more maintaining) which is why it +# is in this script. """ T_reentry(i, j, grid, clock, model_fields, ΔT) From 30afbce63199bacc0aca1c2ac7bbdf7f3b8dacd2 Mon Sep 17 00:00:00 2001 From: jbisits Date: Wed, 20 Nov 2024 16:14:50 +1100 Subject: [PATCH 02/12] Rename to `AbstractBackgroundFunction` --- src/staircase_background.jl | 8 ++++---- src/staircase_initial_conditions.jl | 2 +- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/src/staircase_background.jl b/src/staircase_background.jl index 5ae4244..3ee0189 100644 --- a/src/staircase_background.jl +++ b/src/staircase_background.jl @@ -1,10 +1,10 @@ -Oceananigans.BackgroundField(bf::BackgroundFunction) = +Oceananigans.BackgroundField(bf::AbstractBackgroundFunction) = BackgroundField(bf.func, parameters = bf.parameters) """ mutable struct BackgroundTanh{F, T, P} Container for a tanh background field. """ -mutable struct BackgroundTanh{F, T} <: BackgroundFunction +mutable struct BackgroundTanh{F, T} <: AbstractBackgroundFunction "tanh function" func :: F "Scale the steepness of the `tanh` change" @@ -18,7 +18,7 @@ BackgroundTanh(scale) = BackgroundTanh(tanh_background, scale, NamedTuple()) mutable struct BackgroundLinear{F, P} Container for a linear background field. """ -mutable struct BackgroundLinear{F} <: BackgroundFunction +mutable struct BackgroundLinear{F} <: AbstractBackgroundFunction "Linear function" func :: F "Parameters for the tanh background field" @@ -26,7 +26,7 @@ mutable struct BackgroundLinear{F} <: BackgroundFunction end BackgroundLinear() = BackgroundLinear(linear_background, NamedTuple()) -function Base.show(io::IO, bf::BackgroundFunction) +function Base.show(io::IO, bf::AbstractBackgroundFunction) if bf isa BackgroundTanh println(io, "BackgroundTanh") println(io, "┣━━━ function: $(bf.func)") diff --git a/src/staircase_initial_conditions.jl b/src/staircase_initial_conditions.jl index 4a84bd0..0942832 100644 --- a/src/staircase_initial_conditions.jl +++ b/src/staircase_initial_conditions.jl @@ -57,7 +57,7 @@ struct PeriodicSTSingleInterfaceInitialConditions{T, A, F} <: AbstractStaircaseI temperature_values :: A "Initial R_ρ at the interface" R_ρ :: T - "Function used to define the background state about which an anomaly is evolved." + "`BackgroundField` about which anomaly is advected. Should be an `AbstractBackgroundFunction`." background_state :: F end function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature, background_state) From 84e0f071034947705eb26da3d7216f584a7489cd Mon Sep 17 00:00:00 2001 From: jbisits Date: Wed, 20 Nov 2024 21:19:37 +1100 Subject: [PATCH 03/12] Break up sinlgle interface and staircase ics --- src/StaircaseShenanigans.jl | 4 +- src/single_interfaces_initial_conditions.jl | 102 ++++++++++++++++ src/staircase_initial_conditions.jl | 127 +++----------------- 3 files changed, 122 insertions(+), 111 deletions(-) create mode 100644 src/single_interfaces_initial_conditions.jl diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index adb315a..f507d30 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -22,7 +22,8 @@ import SeawaterPolynomials.SecondOrderSeawaterPolynomials: RoquetSeawaterPolynom abstract type AbstractStaircaseModel end abstract type AbstractInitialConditions end abstract type AbstractNoise end -abstract type BackgroundFunction end +abstract type AbstractBackgroundFunction end +abstract type AbstractSmoothingFunction end export StaircaseDNS, PeriodicStaircaseDNS, DNSModel, SDNS_simulation_setup @@ -47,6 +48,7 @@ export animate_tracers, animate_density, visualise_initial_conditions, visualise animate_tracers_anomaly, animate_density_anomaly include("staircase_initial_conditions.jl") +include("single_interfaces_initial_conditions.jl") include("staircase_noise.jl") include("staircase_diagnostics.jl") include("staircase_model.jl") diff --git a/src/single_interfaces_initial_conditions.jl b/src/single_interfaces_initial_conditions.jl new file mode 100644 index 0000000..0c4fb55 --- /dev/null +++ b/src/single_interfaces_initial_conditions.jl @@ -0,0 +1,102 @@ +abstract type AbstractSingleInterfaceInitialConditions <: AbstractInitialConditions end + +""" + struct STSingleInterfaceInitialConditions +Initial conditions for a single interface (i.e. two layers with uniform `S` and `T` seperated +by a step change). The property `maintain_interface` is a `Boolean` which if set to `true` will +set [reentrant_boundary_conditions](@ref) so that the interface will be maintained (by not +letting the system run down). +""" +struct STSingleInterfaceInitialConditions{T, A, B} <: AbstractSingleInterfaceInitialConditions + "The depth of the interface" + depth_of_interface :: T + "Salinity values in each layer" + salinity_values :: A + "Temperature values in each layer" + temperature_values :: A + "Initial R_ρ at the interface" + R_ρ :: T + "Boolean whether or not to set reentrant boundary condtions to approximately maintain the initial + interface gradients" + maintain_interface :: B +end +function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature; + maintain_interface = false) + + eos = model.buoyancy.model.equation_of_state + + R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) + + return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) + +end +function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature; + maintain_interface = false) + + R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) + + return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) + +end +const SingleInterfaceICs = STSingleInterfaceInitialConditions # alias +Base.summary(ics::SingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface)" + +""" + struct PeriodicSTSingleInterfaceInitialConditions +Sets a `BackgroundField` according to `background_State` and uses a triply periodic domain +to evolve salinity and temperature anomalies about the background state. +""" +struct PeriodicSTSingleInterfaceInitialConditions{T, A, F} <: AbstractSingleInterfaceInitialConditions + "The depth of the interface" + depth_of_interface :: T + "Salinity values in each layer" + salinity_values :: A + "Temperature values in each layer" + temperature_values :: A + "Initial R_ρ at the interface" + R_ρ :: T + "`BackgroundField` about which anomaly is advected. Should be an `AbstractBackgroundFunction`." + background_state :: F +end +function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature, background_state) + + R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) + + return PeriodicSTSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, background_state) + +end +const PeriodoicSingleInterfaceICs = PeriodicSTSingleInterfaceInitialConditions # alias +Base.summary(ics::PeriodoicSingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface) on triply periodic domain with $(summary(ics.background_state)) state" + +function compute_R_ρ(salinity, temperature, depth_of_interfaces::Number, eos) + + S_u = S_g = salinity[1] + S_l = S_f = salinity[2] + T_u = T_f = temperature[1] + T_l = T_g = temperature[2] + + ρ_u = ρ(T_u, S_u, depth_of_interfaces, eos) + ρ_l = ρ(T_l, S_l, depth_of_interfaces, eos) + ρ_f = ρ(T_f, S_f, depth_of_interfaces, eos) + ρ_g = ρ(T_g, S_g, depth_of_interfaces, eos) + + return (0.5 * (ρ_f - ρ_u) + 0.5 * (ρ_l - ρ_g)) / (0.5 * (ρ_f - ρ_l) + 0.5 * (ρ_u - ρ_g)) +end + +function Base.show(io::IO, sics::AbstractSingleInterfaceInitialConditions) + if sics isa PeriodicSTSingleInterfaceInitialConditions + println(io, "STSingleInterfaceInitialConditions") + println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") + println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") + println(io, "┣━ temperature_values: $(sics.temperature_values)") + println(io, "┣━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") + println(io, "┗━━━ background_state: $(summary(sics.background_state))") + elseif sics isa STSingleInterfaceInitialConditions + println(io, "STSingleInterfaceInitialConditions") + println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") + println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") + println(io, "┣━ temperature_values: $(sics.temperature_values)") + println(io, "┣━ maintain_interface: $(round.(sics.maintain_interface))") + print(io, "┗━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") + end +end diff --git a/src/staircase_initial_conditions.jl b/src/staircase_initial_conditions.jl index 0942832..9d88d68 100644 --- a/src/staircase_initial_conditions.jl +++ b/src/staircase_initial_conditions.jl @@ -3,73 +3,6 @@ abstract type AbstractStaircaseInitialConditions <: AbstractInitialConditions en Base.iterate(sics::AbstractStaircaseInitialConditions, state = 1) = state > length(fieldnames(sics)) ? nothing : (getfield(sics, state), state + 1) -""" - struct STSingleInterfaceInitialConditions -Initial conditions for a single interface (i.e. two layers with uniform `S` and `T` seperated -by a step change). The property `maintain_interface` is a `Boolean` which if set to `true` will -set [reentrant_boundary_conditions](@ref) so that the interface will be maintained (by not -letting the system run down). -""" -struct STSingleInterfaceInitialConditions{T, A, B} <: AbstractStaircaseInitialConditions - "The depth of the interface" - depth_of_interface :: T - "Salinity values in each layer" - salinity_values :: A - "Temperature values in each layer" - temperature_values :: A - "Initial R_ρ at the interface" - R_ρ :: T - "Boolean whether or not to set reentrant boundary condtions to approximately maintain the initial - interface gradients" - maintain_interface :: B -end -function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature; - maintain_interface = false) - - eos = model.buoyancy.model.equation_of_state - - R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - - return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) - -end -function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature; - maintain_interface = false) - - R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - - return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) - -end -const SingleInterfaceICs = STSingleInterfaceInitialConditions # alias -Base.summary(ics::SingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface)" -""" - struct PeriodicSTSingleInterfaceInitialConditions -Sets a `BackgroundField` according to `background_State` and uses a triply periodic domain -to evolve salinity and temperature anomalies about the background state. -""" -struct PeriodicSTSingleInterfaceInitialConditions{T, A, F} <: AbstractStaircaseInitialConditions - "The depth of the interface" - depth_of_interface :: T - "Salinity values in each layer" - salinity_values :: A - "Temperature values in each layer" - temperature_values :: A - "Initial R_ρ at the interface" - R_ρ :: T - "`BackgroundField` about which anomaly is advected. Should be an `AbstractBackgroundFunction`." - background_state :: F -end -function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature, background_state) - - R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - - return PeriodicSTSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, background_state) - -end -const PeriodoicSingleInterfaceICs = PeriodicSTSingleInterfaceInitialConditions # alias -Base.summary(ics::PeriodoicSingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface) on triply periodic domain with $(summary(ics.background_state)) state" - "Container for initial conditions that have well mixed layers seperated by sharp step interfaces." struct STStaircaseInitialConditions{T} <: AbstractStaircaseInitialConditions "Number of interfaces in the initial state" @@ -95,6 +28,22 @@ end const StaircaseICs = STStaircaseInitialConditions # alias Base.summary(ics::StaircaseICs) = "Multiple S-T interfaces at z = $(ics.depth_of_interfaces)" + +# TODO: implement this so it is an option. +"Container for initial conditions that have well mixed layers seperated by smoothed step interfaces." +struct SmoothSTStaircaseInitialConditions{T} <: AbstractStaircaseInitialConditions + "Number of interfaces in the initial state" + number_of_interfaces :: Int + "The depth of the interfaces, **note:** length (depth_of_interfaces) == number_of_interfaces" + depth_of_interfaces :: T + "Salinity values in each layer" + salinity_values :: T + "Temperature values in each layer" + temperature_values :: T + "Function to smooth step transition" + smoothing_funciton :: Function +end + """ function compute_R_ρ(salinity, temperature, depth_of_interfaces, eos) Compute the density ratio, ``R_{\rho}``, at a diffusive interface with a non-linear equation of state @@ -121,51 +70,9 @@ function compute_R_ρ(salinity, temperature, depth_of_interfaces::Array, eos) return R_ρ end -function compute_R_ρ(salinity, temperature, depth_of_interfaces::Number, eos) - - S_u = S_g = salinity[1] - S_l = S_f = salinity[2] - T_u = T_f = temperature[1] - T_l = T_g = temperature[2] - - ρ_u = ρ(T_u, S_u, depth_of_interfaces, eos) - ρ_l = ρ(T_l, S_l, depth_of_interfaces, eos) - ρ_f = ρ(T_f, S_f, depth_of_interfaces, eos) - ρ_g = ρ(T_g, S_g, depth_of_interfaces, eos) - - return (0.5 * (ρ_f - ρ_u) + 0.5 * (ρ_l - ρ_g)) / (0.5 * (ρ_f - ρ_l) + 0.5 * (ρ_u - ρ_g)) -end -# TODO: implement this so it is an option. -"Container for initial conditions that have well mixed layers seperated by smoothed step interfaces." -struct SmoothSTStaircaseInitialConditions{T} <: AbstractStaircaseInitialConditions - "Number of interfaces in the initial state" - number_of_interfaces :: Int - "The depth of the interfaces, **note:** length (depth_of_interfaces) == number_of_interfaces" - depth_of_interfaces :: T - "Salinity values in each layer" - salinity_values :: T - "Temperature values in each layer" - temperature_values :: T - "Function to smooth step transition" - smoothing_funciton :: Function -end function Base.show(io::IO, sics::AbstractStaircaseInitialConditions) - if sics isa PeriodicSTSingleInterfaceInitialConditions - println(io, "STSingleInterfaceInitialConditions") - println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") - println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") - println(io, "┣━ temperature_values: $(sics.temperature_values)") - println(io, "┣━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") - println(io, "┗━━━ background_state: $(summary(sics.background_state))") - elseif sics isa STSingleInterfaceInitialConditions - println(io, "STSingleInterfaceInitialConditions") - println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") - println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") - println(io, "┣━ temperature_values: $(sics.temperature_values)") - println(io, "┣━ maintain_interface: $(round.(sics.maintain_interface))") - print(io, "┗━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") - elseif sics isa STStaircaseInitialConditions + if sics isa STStaircaseInitialConditions println(io, "STStaircaseInitialConditions") println(io, "┣━ number_of_interfaces: $(sics.number_of_interfaces)") println(io, "┣━━ depth_of_interfaces: $(sics.depth_of_interfaces)") From e1d07f041aeb91033f9d017682510b3630a78c51 Mon Sep 17 00:00:00 2001 From: jbisits Date: Sat, 23 Nov 2024 21:28:57 +1100 Subject: [PATCH 04/12] Add option for interface smoothing function Plus remove remanants of boundary restoring. --- src/StaircaseShenanigans.jl | 3 +- src/interface_smoothing.jl | 20 +++++++ src/single_interfaces_initial_conditions.jl | 61 ++++++++++++--------- src/staircase_initial_conditions.jl | 4 +- src/staircase_model.jl | 3 +- src/staircase_restoring.jl | 58 -------------------- 6 files changed, 59 insertions(+), 90 deletions(-) create mode 100644 src/interface_smoothing.jl diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index f507d30..d94cf74 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -23,7 +23,7 @@ abstract type AbstractStaircaseModel end abstract type AbstractInitialConditions end abstract type AbstractNoise end abstract type AbstractBackgroundFunction end -abstract type AbstractSmoothingFunction end +abstract type AbstractInterfaceSmoothing end export StaircaseDNS, PeriodicStaircaseDNS, DNSModel, SDNS_simulation_setup @@ -47,6 +47,7 @@ export compute_R_ρ! export animate_tracers, animate_density, visualise_initial_conditions, visualise_initial_density, animate_tracers_anomaly, animate_density_anomaly +include("interface_smoothing.jl") include("staircase_initial_conditions.jl") include("single_interfaces_initial_conditions.jl") include("staircase_noise.jl") diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl new file mode 100644 index 0000000..bf169bb --- /dev/null +++ b/src/interface_smoothing.jl @@ -0,0 +1,20 @@ + +""" + struct TanhInterfaceSmoothing +Container to set a hyperbolic tangent over a single interface as an initial condition. +""" +struct TanhInterfaceSmoothing{T, A} <: AbstractInterfaceSmoothing + "Tracer in the deeper of the two layers seperated by an interface" + Cₗ :: T + "Change in tracer content over an interface" + ΔC :: T + "Steepness of the tanh change" + D :: T + "Location of the interface" + z_interface :: T + "Total z range of the two layers seperated by a given interface" + z_range :: A +end +@inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.Lz)) + +Base.summary(p::TanhInterfaceSmoothing) = "tanh smoothing." diff --git a/src/single_interfaces_initial_conditions.jl b/src/single_interfaces_initial_conditions.jl index 0c4fb55..7679654 100644 --- a/src/single_interfaces_initial_conditions.jl +++ b/src/single_interfaces_initial_conditions.jl @@ -1,5 +1,7 @@ abstract type AbstractSingleInterfaceInitialConditions <: AbstractInitialConditions end +Base.iterate(sics::AbstractSingleInterfaceInitialConditions, state = 1) = + state > length(fieldnames(sics)) ? nothing : (getfield(sics, state), state + 1) """ struct STSingleInterfaceInitialConditions Initial conditions for a single interface (i.e. two layers with uniform `S` and `T` seperated @@ -7,35 +9,35 @@ by a step change). The property `maintain_interface` is a `Boolean` which if set set [reentrant_boundary_conditions](@ref) so that the interface will be maintained (by not letting the system run down). """ -struct STSingleInterfaceInitialConditions{T, A, B} <: AbstractSingleInterfaceInitialConditions +struct STSingleInterfaceInitialConditions{T, A, IS} <: AbstractSingleInterfaceInitialConditions "The depth of the interface" - depth_of_interface :: T + depth_of_interface :: T "Salinity values in each layer" - salinity_values :: A + salinity_values :: A "Temperature values in each layer" - temperature_values :: A + temperature_values :: A + "Initial smoothing function over interface" + interface_smoothing :: IS "Initial R_ρ at the interface" - R_ρ :: T - "Boolean whether or not to set reentrant boundary condtions to approximately maintain the initial - interface gradients" - maintain_interface :: B + R_ρ :: T end -function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature; - maintain_interface = false) +function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature) eos = model.buoyancy.model.equation_of_state R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) + return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, + interface_smoothing, R_ρ) end -function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature; - maintain_interface = false) +function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, + salinity, temperature) R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, maintain_interface) + return STSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, + interface_smoothing, R_ρ) end const SingleInterfaceICs = STSingleInterfaceInitialConditions # alias @@ -46,23 +48,27 @@ Base.summary(ics::SingleInterfaceICs) = "Single S-T interface at z = $(ics.depth Sets a `BackgroundField` according to `background_State` and uses a triply periodic domain to evolve salinity and temperature anomalies about the background state. """ -struct PeriodicSTSingleInterfaceInitialConditions{T, A, F} <: AbstractSingleInterfaceInitialConditions +struct PeriodicSTSingleInterfaceInitialConditions{T, A, IS, BF} <: AbstractSingleInterfaceInitialConditions "The depth of the interface" - depth_of_interface :: T + depth_of_interface :: T "Salinity values in each layer" - salinity_values :: A + salinity_values :: A "Temperature values in each layer" - temperature_values :: A + temperature_values :: A + "Initial smoothing function over interface" + interface_smoothing :: IS "Initial R_ρ at the interface" - R_ρ :: T + R_ρ :: T "`BackgroundField` about which anomaly is advected. Should be an `AbstractBackgroundFunction`." - background_state :: F + background_state :: BF end -function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature, background_state) +function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, + salinity, temperature, interface_smoothing, background_state) R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) - return PeriodicSTSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, R_ρ, background_state) + return PeriodicSTSingleInterfaceInitialConditions(depth_of_interface, salinity, temperature, + interface_smoothing, R_ρ, background_state) end const PeriodoicSingleInterfaceICs = PeriodicSTSingleInterfaceInitialConditions # alias @@ -89,14 +95,15 @@ function Base.show(io::IO, sics::AbstractSingleInterfaceInitialConditions) println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") println(io, "┣━ temperature_values: $(sics.temperature_values)") + println(io, "┣━ interface_smoothing: $(sumary(interface_smoothing))") println(io, "┣━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") println(io, "┗━━━ background_state: $(summary(sics.background_state))") elseif sics isa STSingleInterfaceInitialConditions println(io, "STSingleInterfaceInitialConditions") - println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") - println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") - println(io, "┣━ temperature_values: $(sics.temperature_values)") - println(io, "┣━ maintain_interface: $(round.(sics.maintain_interface))") - print(io, "┗━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") + println(io, "┣━━ depth_of_interface: $(sics.depth_of_interface)") + println(io, "┣━━━━━ salinity_values: $(sics.salinity_values)") + println(io, "┣━━ temperature_values: $(sics.temperature_values)") + println(io, "┣━ interface_smoothing: $(sumary(interface_smoothing))") + print(io, "┗━━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") end end diff --git a/src/staircase_initial_conditions.jl b/src/staircase_initial_conditions.jl index 9d88d68..ddfd001 100644 --- a/src/staircase_initial_conditions.jl +++ b/src/staircase_initial_conditions.jl @@ -31,7 +31,7 @@ Base.summary(ics::StaircaseICs) = "Multiple S-T interfaces at z = $(ics.depth_of # TODO: implement this so it is an option. "Container for initial conditions that have well mixed layers seperated by smoothed step interfaces." -struct SmoothSTStaircaseInitialConditions{T} <: AbstractStaircaseInitialConditions +struct SmoothSTStaircaseInitialConditions{T, F} <: AbstractStaircaseInitialConditions "Number of interfaces in the initial state" number_of_interfaces :: Int "The depth of the interfaces, **note:** length (depth_of_interfaces) == number_of_interfaces" @@ -41,7 +41,7 @@ struct SmoothSTStaircaseInitialConditions{T} <: AbstractStaircaseInitialConditio "Temperature values in each layer" temperature_values :: T "Function to smooth step transition" - smoothing_funciton :: Function + smoothing_funciton :: F end """ diff --git a/src/staircase_model.jl b/src/staircase_model.jl index e01d7d6..571e208 100644 --- a/src/staircase_model.jl +++ b/src/staircase_model.jl @@ -34,9 +34,8 @@ The initial conditions are set after building the `model`. """ function StaircaseDNS(model_setup::NamedTuple, initial_conditions, initial_noise) - boundary_conditions = reentrant_boundary_conditions(initial_conditions) architecture, diffusivities, domain_extent, resolution, eos = model_setup - model = DNSModel(architecture, domain_extent, resolution, diffusivities, eos; boundary_conditions) + model = DNSModel(architecture, domain_extent, resolution, diffusivities, eos) sdns = StaircaseDNS(model, initial_conditions, initial_noise) set_staircase_initial_conditions!(sdns) diff --git a/src/staircase_restoring.jl b/src/staircase_restoring.jl index e1a7517..4f0b85c 100644 --- a/src/staircase_restoring.jl +++ b/src/staircase_restoring.jl @@ -1,61 +1,3 @@ -# The following functions are to be used as `BoundaryConditions` so that tracers can -# re-enter the domain with the initial gradient added effectively allowing the gradient to -# be maintained. This is effectively a restoring (or more maintaining) which is why it -# is in this script. - -""" - T_reentry(i, j, grid, clock, model_fields, ΔT) -Discrete boundary condition to set the temperature on a vertical boundary to tracer value at -the **top** of the domain and adds ΔT. That is a reentrant T condition to be used on the -**bottom** of the domain. -""" -@inline T_reentry(i, j, grid, clock, model_fields, ΔT) = - @inbounds ℑzᵃᵃᶠ(i, j, grid.Nz+1, model_fields.T) + ΔT -""" - S_reentry(i, j, grid, clock, model_fields, ΔS) -As for [T_reentry](@ref) but using salinity tracer instead. -""" -@inline S_reentry(i, j, grid, clock, model_fields, ΔS) = - @inbounds ℑzᵃᵃᶠ(i, j, grid.Nz+1, model_fields.S) + ΔS - -"Access velocity field at the _bottom_ of the domain for use as `BoundaryCondition`." -@inline _w_bottom(i, j, grid, clock, model_fields) = @inbounds model_fields.w[i, j, 1] -"Access velocity field at the _top_ of the domain for use as `BoundaryCondition`." -@inline _w_top(i, j, grid, clock, model_fields) = @inbounds model_fields.w[i, j, grid.Nz+1] - -"Interpolated velocity between top and bottom (vertical) faces of domain." -@inline w_top_bottom_interpolate(i, j, grid, clock, model_fields) = - @inbounds 0.5 * (model_fields.w[i, j, 1] + model_fields.w[i, j, grid.Nz+1]) - -"""" - function reentrant_boundary_conditions(ics::SingleInterfaceICs) -Setup boundary conditions to maintain a gradient (ΔC) between the two layers of a -`SingleInterface` model. The boundary conditions for the tracers are re-entrant from top to -bottom with ΔC added to maintain a ΔC difference between the top and bottom layer. The -vertical velocitiy field also has modified vertical boundary conditions where the velocity -on the bottom face is set to the velocity on the top face and vice versa. -""" -function reentrant_boundary_conditions(ics::SingleInterfaceICs) - - bcs = if ics.maintain_interface - ΔT = diff(ics.temperature_values)[1] - ΔS = diff(ics.salinity_values)[1] - T_bottom_reentry = ValueBoundaryCondition(T_reentry, discrete_form=true, parameters = ΔT) - S_bottom_reentry = ValueBoundaryCondition(S_reentry, discrete_form=true, parameters = ΔS) - T_bcs = FieldBoundaryConditions(bottom = T_bottom_reentry) - S_bcs = FieldBoundaryConditions(bottom = S_bottom_reentry) - - w_top = OpenBoundaryCondition(_w_bottom, discrete_form=true) - w_bottom = OpenBoundaryCondition(_w_top, discrete_form=true) - w_bcs = FieldBoundaryConditions(top = w_top, bottom = w_bottom) - - (T=T_bcs, S=S_bcs, w=w_bcs) - else - NamedTuple() - end - - return bcs -end """ struct OuterStairMask{T} Container for location of first and last stair for `mask` in `Relaxation`. This need not be From 3885bacab002dc3d0f39b5f0f0795d7a5d3a6bcc Mon Sep 17 00:00:00 2001 From: jbisits Date: Sat, 23 Nov 2024 21:51:05 +1100 Subject: [PATCH 05/12] Done for `nothing` in SIICS --- examples/single_interface.jl | 2 +- src/StaircaseShenanigans.jl | 2 ++ src/interface_smoothing.jl | 1 + src/set_staircase_initial_conditions.jl | 19 ++++++++++++++++++- src/single_interfaces_initial_conditions.jl | 9 +++++---- src/staircase_model.jl | 4 ++-- 6 files changed, 29 insertions(+), 8 deletions(-) diff --git a/examples/single_interface.jl b/examples/single_interface.jl index 09eca30..bb1642f 100644 --- a/examples/single_interface.jl +++ b/examples/single_interface.jl @@ -11,7 +11,7 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.58, 34.70] temperature = [-1.5, 0.5] -interface_ics = SingleInterfaceICs(eos, depth_of_interface, salinity, temperature, maintain_interface = true) +interface_ics = SingleInterfaceICs(eos, depth_of_interface, salinity, temperature) velocity_noise = VelocityNoise(0.0, 0.0, 1e-7) ## setup model diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index d94cf74..f83d4d1 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -32,6 +32,8 @@ export STStaircaseInitialConditions, StaircaseICs, SmoothSTStaircaseInitialCondi PeriodicSTSingleInterfaceInitialConditions, PeriodoicSingleInterfaceICs, set_staircase_initial_conditions! +export TanhInterfaceSmoothing + export BackgroundTanh, BackgroundLinear, tanh_background, linear_background export AbstractNoise, VelocityNoise, TracerNoise diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl index bf169bb..093d1cf 100644 --- a/src/interface_smoothing.jl +++ b/src/interface_smoothing.jl @@ -18,3 +18,4 @@ end @inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.Lz)) Base.summary(p::TanhInterfaceSmoothing) = "tanh smoothing." +Base.summary(p::Nothing) = "No smoothing." diff --git a/src/set_staircase_initial_conditions.jl b/src/set_staircase_initial_conditions.jl index 3299a85..168ecd5 100644 --- a/src/set_staircase_initial_conditions.jl +++ b/src/set_staircase_initial_conditions.jl @@ -33,7 +33,9 @@ function set_staircase_initial_conditions!(model, ics::SmoothSTStaircaseInitialC return nothing end -function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs) +set_staircase_initial_conditions!(model, ics::SingleInterfaceICs) = + set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, ics.interface_smoothing) +function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, interface_smoothing::Nothing) depth_of_interface = ics.depth_of_interface z = znodes(model.grid, Center()) @@ -51,6 +53,21 @@ function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs) return nothing end +function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, + interface_smoothing::TanhInterfaceSmoothing) + + depth_of_interface = ics.depth_of_interface + Lz = model.grid.Lz + S = ics.salinity_values + T = ics.temperature_values + + S₀ = TanhInterfaceSmoothing(S[2], diff(S)[1], 100, depth_of_interface, abs(Lz)) + T₀ = TanhInterfaceSmoothing(T[2], diff(S)[1], 100, depth_of_interface, abs(Lz)) + + set!(model, S = S₀, T = T₀) + + return nothing +end "Currently only random noise can be added as the initial condition in the tracer field when there is a `BackgroundField`." set_staircase_initial_conditions!(model, ics::PeriodoicSingleInterfaceICs) = nothing diff --git a/src/single_interfaces_initial_conditions.jl b/src/single_interfaces_initial_conditions.jl index 7679654..164fb8d 100644 --- a/src/single_interfaces_initial_conditions.jl +++ b/src/single_interfaces_initial_conditions.jl @@ -21,7 +21,8 @@ struct STSingleInterfaceInitialConditions{T, A, IS} <: AbstractSingleInterfaceIn "Initial R_ρ at the interface" R_ρ :: T end -function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature) +function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, temperature; + interface_smoothing = nothing) eos = model.buoyancy.model.equation_of_state @@ -32,7 +33,7 @@ function STSingleInterfaceInitialConditions(model, depth_of_interface, salinity, end function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, - salinity, temperature) + salinity, temperature; interface_smoothing = nothing) R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) @@ -95,7 +96,7 @@ function Base.show(io::IO, sics::AbstractSingleInterfaceInitialConditions) println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") println(io, "┣━ temperature_values: $(sics.temperature_values)") - println(io, "┣━ interface_smoothing: $(sumary(interface_smoothing))") + println(io, "┣━ interface_smoothing: $(summary(sics.interface_smoothing))") println(io, "┣━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") println(io, "┗━━━ background_state: $(summary(sics.background_state))") elseif sics isa STSingleInterfaceInitialConditions @@ -103,7 +104,7 @@ function Base.show(io::IO, sics::AbstractSingleInterfaceInitialConditions) println(io, "┣━━ depth_of_interface: $(sics.depth_of_interface)") println(io, "┣━━━━━ salinity_values: $(sics.salinity_values)") println(io, "┣━━ temperature_values: $(sics.temperature_values)") - println(io, "┣━ interface_smoothing: $(sumary(interface_smoothing))") + println(io, "┣━ interface_smoothing: $(summary(sics.interface_smoothing))") print(io, "┗━━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") end end diff --git a/src/staircase_model.jl b/src/staircase_model.jl index 571e208..8eb7a55 100644 --- a/src/staircase_model.jl +++ b/src/staircase_model.jl @@ -1,5 +1,5 @@ struct StaircaseDNS{NHM <: NonhydrostaticModel, - SIC <: AbstractStaircaseInitialConditions, + SIC <: AbstractInitialConditions, AN <: Union{AbstractNoise, Nothing}} <: AbstractStaircaseModel "An [Oceananigans.jl `NonhydrostaticModel`](https://clima.github.io/OceananigansDocumentation/dev/appendix/library/#Oceananigans.Models.NonhydrostaticModels.NonhydrostaticModel-Tuple{})" model :: NHM @@ -32,7 +32,7 @@ Build the model from `model_setup` then return a `StaircaseDNS`, mainly used to the `model` with [reentrant_boundary_conditions](@ref) `boundary_conditions` from `initial_conditions`. The initial conditions are set after building the `model`. """ -function StaircaseDNS(model_setup::NamedTuple, initial_conditions, initial_noise) +function StaircaseDNS(model_setup::NamedTuple, initial_conditions::SingleInterfaceICs, initial_noise) architecture, diffusivities, domain_extent, resolution, eos = model_setup model = DNSModel(architecture, domain_extent, resolution, diffusivities, eos) From 3247f5f8c099d4d521e0f63b632e0362a47c997a Mon Sep 17 00:00:00 2001 From: jbisits Date: Sat, 23 Nov 2024 22:34:33 +1100 Subject: [PATCH 06/12] Implemented for non-perioidic single_interface --- examples/single_interface.jl | 4 +--- src/StaircaseShenanigans.jl | 2 +- src/interface_smoothing.jl | 9 +++++---- src/set_staircase_initial_conditions.jl | 6 +++--- 4 files changed, 10 insertions(+), 11 deletions(-) diff --git a/examples/single_interface.jl b/examples/single_interface.jl index bb1642f..5727abc 100644 --- a/examples/single_interface.jl +++ b/examples/single_interface.jl @@ -11,7 +11,7 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.58, 34.70] temperature = [-1.5, 0.5] -interface_ics = SingleInterfaceICs(eos, depth_of_interface, salinity, temperature) +interface_ics = SingleInterfaceICs(eos, depth_of_interface, salinity, temperature, interface_smoothing = Tanh) velocity_noise = VelocityNoise(0.0, 0.0, 1e-7) ## setup model @@ -28,7 +28,5 @@ simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, output_path, max_Δt = 5) ## Run run!(simulation) - -## Compute density ratio compute_R_ρ!(simulation.output_writers[:computed_output].filepath, simulation.output_writers[:tracers].filepath, eos) diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index f83d4d1..6abbb5b 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -32,7 +32,7 @@ export STStaircaseInitialConditions, StaircaseICs, SmoothSTStaircaseInitialCondi PeriodicSTSingleInterfaceInitialConditions, PeriodoicSingleInterfaceICs, set_staircase_initial_conditions! -export TanhInterfaceSmoothing +export TanhInterfaceSmoothing, Tanh export BackgroundTanh, BackgroundLinear, tanh_background, linear_background diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl index 093d1cf..291725f 100644 --- a/src/interface_smoothing.jl +++ b/src/interface_smoothing.jl @@ -3,7 +3,7 @@ struct TanhInterfaceSmoothing Container to set a hyperbolic tangent over a single interface as an initial condition. """ -struct TanhInterfaceSmoothing{T, A} <: AbstractInterfaceSmoothing +struct TanhInterfaceSmoothing{T} <: AbstractInterfaceSmoothing "Tracer in the deeper of the two layers seperated by an interface" Cₗ :: T "Change in tracer content over an interface" @@ -13,9 +13,10 @@ struct TanhInterfaceSmoothing{T, A} <: AbstractInterfaceSmoothing "Location of the interface" z_interface :: T "Total z range of the two layers seperated by a given interface" - z_range :: A + z_range :: T end -@inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.Lz)) +const Tanh = TanhInterfaceSmoothing{T} where {T} +@inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.z_range)) -Base.summary(p::TanhInterfaceSmoothing) = "tanh smoothing." +Base.summary(p::Type{<:TanhInterfaceSmoothing}) = "tanh smoothing." Base.summary(p::Nothing) = "No smoothing." diff --git a/src/set_staircase_initial_conditions.jl b/src/set_staircase_initial_conditions.jl index 168ecd5..df9963a 100644 --- a/src/set_staircase_initial_conditions.jl +++ b/src/set_staircase_initial_conditions.jl @@ -54,15 +54,15 @@ function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, inter return nothing end function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, - interface_smoothing::TanhInterfaceSmoothing) + interface_smoothing::Type{<:Tanh}) depth_of_interface = ics.depth_of_interface Lz = model.grid.Lz S = ics.salinity_values T = ics.temperature_values - S₀ = TanhInterfaceSmoothing(S[2], diff(S)[1], 100, depth_of_interface, abs(Lz)) - T₀ = TanhInterfaceSmoothing(T[2], diff(S)[1], 100, depth_of_interface, abs(Lz)) + S₀(x, y, z) = Tanh(S[2], diff(S)[1], 100.0, depth_of_interface, abs(Lz))(x, y, z) + T₀(x, y, z) = Tanh(T[2], diff(T)[1], 100.0, depth_of_interface, abs(Lz))(x, y, z) set!(model, S = S₀, T = T₀) From 55f242860f2670349eeb62d8bc95031a9c78c7ce Mon Sep 17 00:00:00 2001 From: jbisits Date: Sat, 23 Nov 2024 23:10:59 +1100 Subject: [PATCH 07/12] Changes for periodic initial conditions Looks to all be setup correctly --- ...ngle_interface_periodic_tanh_background.jl | 5 +- src/interface_smoothing.jl | 2 +- src/set_staircase_initial_conditions.jl | 16 ++-- src/single_interfaces_initial_conditions.jl | 16 ++-- src/staircase_background.jl | 77 ++++++++++--------- src/staircase_model.jl | 4 +- 6 files changed, 62 insertions(+), 58 deletions(-) diff --git a/examples/single_interface_periodic_tanh_background.jl b/examples/single_interface_periodic_tanh_background.jl index 7af175e..4d195a7 100644 --- a/examples/single_interface_periodic_tanh_background.jl +++ b/examples/single_interface_periodic_tanh_background.jl @@ -11,8 +11,7 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.56, 34.70] temperature = [-1.5, 0.5] -interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, - BackgroundTanh()) +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, interface_smoothing = Tanh) tracer_noise = TracerNoise(1e-6, 1e-6) ## setup model @@ -30,6 +29,6 @@ simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, ## Run run!(simulation) -## Compute density ratio +# Compute density ratio compute_R_ρ!(simulation.output_writers[:computed_output].filepath, simulation.output_writers[:tracers].filepath, eos) diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl index 291725f..aa42b14 100644 --- a/src/interface_smoothing.jl +++ b/src/interface_smoothing.jl @@ -19,4 +19,4 @@ const Tanh = TanhInterfaceSmoothing{T} where {T} @inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.z_range)) Base.summary(p::Type{<:TanhInterfaceSmoothing}) = "tanh smoothing." -Base.summary(p::Nothing) = "No smoothing." +Base.summary(p::Nothing) = "$p" diff --git a/src/set_staircase_initial_conditions.jl b/src/set_staircase_initial_conditions.jl index df9963a..c192096 100644 --- a/src/set_staircase_initial_conditions.jl +++ b/src/set_staircase_initial_conditions.jl @@ -33,9 +33,9 @@ function set_staircase_initial_conditions!(model, ics::SmoothSTStaircaseInitialC return nothing end -set_staircase_initial_conditions!(model, ics::SingleInterfaceICs) = - set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, ics.interface_smoothing) -function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, interface_smoothing::Nothing) +set_staircase_initial_conditions!(model, ics::Union{SingleInterfaceICs, PeriodoicSingleInterfaceICs}) = + set_staircase_initial_conditions!(model, ics, ics.interface_smoothing) +function set_staircase_initial_conditions!(model, ics, interface_smoothing::Nothing) depth_of_interface = ics.depth_of_interface z = znodes(model.grid, Center()) @@ -53,8 +53,7 @@ function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, inter return nothing end -function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, - interface_smoothing::Type{<:Tanh}) +function set_staircase_initial_conditions!(model, ics, interface_smoothing::Type{<:Tanh}) depth_of_interface = ics.depth_of_interface Lz = model.grid.Lz @@ -68,9 +67,6 @@ function set_staircase_initial_conditions!(model, ics::SingleInterfaceICs, return nothing end -"Currently only random noise can be added as the initial condition in the tracer field -when there is a `BackgroundField`." -set_staircase_initial_conditions!(model, ics::PeriodoicSingleInterfaceICs) = nothing "Fallback --- don't set any noise." set_noise!(model, noise::Nothing) = nothing @@ -97,8 +93,10 @@ function set_noise!(model, noise::TracerNoise) S, T = model.tracers S_noise = noise.S_magnitude * randn(size(S)) T_noise = noise.T_magnitude * randn(size(T)) + S₀ = S + S_noise + T₀ = T + T_noise - set!(model, S = S_noise, T = T_noise) + set!(model, S = S₀, T = T₀) end """ diff --git a/src/single_interfaces_initial_conditions.jl b/src/single_interfaces_initial_conditions.jl index 164fb8d..cc87bb7 100644 --- a/src/single_interfaces_initial_conditions.jl +++ b/src/single_interfaces_initial_conditions.jl @@ -64,7 +64,9 @@ struct PeriodicSTSingleInterfaceInitialConditions{T, A, IS, BF} <: AbstractSingl background_state :: BF end function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, - salinity, temperature, interface_smoothing, background_state) + salinity, temperature; + interface_smoothing = nothing, + background_state = nothing) R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) @@ -92,13 +94,13 @@ end function Base.show(io::IO, sics::AbstractSingleInterfaceInitialConditions) if sics isa PeriodicSTSingleInterfaceInitialConditions - println(io, "STSingleInterfaceInitialConditions") - println(io, "┣━ depth_of_interface: $(sics.depth_of_interface)") - println(io, "┣━━━━ salinity_values: $(sics.salinity_values)") - println(io, "┣━ temperature_values: $(sics.temperature_values)") + println(io, "PeriodicSTSingleInterfaceInitialConditions") + println(io, "┣━━ depth_of_interface: $(sics.depth_of_interface)") + println(io, "┣━━━━━ salinity_values: $(sics.salinity_values)") + println(io, "┣━━ temperature_values: $(sics.temperature_values)") println(io, "┣━ interface_smoothing: $(summary(sics.interface_smoothing))") - println(io, "┣━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") - println(io, "┗━━━ background_state: $(summary(sics.background_state))") + println(io, "┣━━━━ background_state: $(summary(sics.background_state))") + print(io, "┗━━━━━━━━━━━━━━━━━ R_ρ: $(round.(sics.R_ρ; digits = 2))") elseif sics isa STSingleInterfaceInitialConditions println(io, "STSingleInterfaceInitialConditions") println(io, "┣━━ depth_of_interface: $(sics.depth_of_interface)") diff --git a/src/staircase_background.jl b/src/staircase_background.jl index 3ee0189..8305358 100644 --- a/src/staircase_background.jl +++ b/src/staircase_background.jl @@ -40,6 +40,7 @@ function Base.show(io::IO, bf::AbstractBackgroundFunction) end Base.summary(bt::BackgroundTanh) = "$(bt.func)" Base.summary(bl::BackgroundLinear) = "$(bl.func)" +Base.summary(bn::Nothing) = "$(bn)" """ function S_and_T_background_fields(initial_conditions) @@ -92,46 +93,48 @@ save_background_state!(simulation, sdns) = save_background_state!(simulation, sd save_background_state!(simulation, model, initial_conditions) = nothing function save_background_state!(simulation, model, initial_conditions::PeriodoicSingleInterfaceICs) - S_background = Field(model.background_fields.tracers.S) - compute!(S_background) - S_background_array = Array(interior(S_background, :, :, :)) - T_background = Field(model.background_fields.tracers.T) - compute!(T_background) - T_background_array = Array(interior(T_background, :, :, :)) - σ_background = Field(seawater_density(model, temperature = T_background, salinity = S_background, - geopotential_height = 0)) - compute!(σ_background) - σ_background_array = Array(interior(σ_background, :, :, :)) - - if simulation.output_writers[:tracers] isa NetCDFOutputWriter - - NCDataset(simulation.output_writers[:tracers].filepath, "a") do ds - defVar(ds, "S_background", S_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for salinity", - "units" => "gkg⁻¹")) - defVar(ds, "T_background", T_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for temperature", - "units" => "°C")) - end - - NCDataset(simulation.output_writers[:computed_output].filepath, "a") do ds - defVar(ds, "σ_background", σ_background_array, ("xC", "yC", "zC"), - attrib = Dict("longname" => "Background field for potential density (0dbar) computed from the `S` and `T` background fields", - "units" => "kgm⁻³")) - end - - elseif simulation.output_writers[:tracers] isa JLD2OutputWriter + if !isnothing(initial_conditions.background_state) + S_background = Field(model.background_fields.tracers.S) + compute!(S_background) + S_background_array = Array(interior(S_background, :, :, :)) + T_background = Field(model.background_fields.tracers.T) + compute!(T_background) + T_background_array = Array(interior(T_background, :, :, :)) + σ_background = Field(seawater_density(model, temperature = T_background, salinity = S_background, + geopotential_height = 0)) + compute!(σ_background) + σ_background_array = Array(interior(σ_background, :, :, :)) + + if simulation.output_writers[:tracers] isa NetCDFOutputWriter + + NCDataset(simulation.output_writers[:tracers].filepath, "a") do ds + defVar(ds, "S_background", S_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for salinity", + "units" => "gkg⁻¹")) + defVar(ds, "T_background", T_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for temperature", + "units" => "°C")) + end + + NCDataset(simulation.output_writers[:computed_output].filepath, "a") do ds + defVar(ds, "σ_background", σ_background_array, ("xC", "yC", "zC"), + attrib = Dict("longname" => "Background field for potential density (0dbar) computed from the `S` and `T` background fields", + "units" => "kgm⁻³")) + end + + elseif simulation.output_writers[:tracers] isa JLD2OutputWriter + + jldopen(simulation.output_writers[:tracers].filepath, "a+") do f + f["S_background"] = S_background_array + f["T_background"] = T_background_array + end + + jldopen(simulation.output_writers[:computed_output].filepath, "a+") do f + f["σ_background"] = σ_background_array + end - jldopen(simulation.output_writers[:tracers].filepath, "a+") do f - f["S_background"] = S_background_array - f["T_background"] = T_background_array - end - - jldopen(simulation.output_writers[:computed_output].filepath, "a+") do f - f["σ_background"] = σ_background_array end end - return nothing end diff --git a/src/staircase_model.jl b/src/staircase_model.jl index 8eb7a55..f8b846f 100644 --- a/src/staircase_model.jl +++ b/src/staircase_model.jl @@ -52,7 +52,9 @@ function StaircaseDNS(model_setup::NamedTuple, initial_conditions::PeriodoicSing architecture, diffusivities, domain_extent, resolution, eos = model_setup z_topology = Periodic - background_fields = S_and_T_background_fields(initial_conditions, domain_extent.Lz) + background_fields = isnothing(initial_conditions.background_state) ? NamedTuple() : + S_and_T_background_fields(initial_conditions, domain_extent.Lz) + model = DNSModel(architecture, domain_extent, resolution, diffusivities, eos; z_topology, background_fields) From 297e192a8506ff5881d572816a9b7b7b8850bf0a Mon Sep 17 00:00:00 2001 From: jbisits Date: Sat, 23 Nov 2024 23:25:44 +1100 Subject: [PATCH 08/12] Remove undefined function --- src/StaircaseShenanigans.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index 6abbb5b..a143440 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -25,7 +25,7 @@ abstract type AbstractNoise end abstract type AbstractBackgroundFunction end abstract type AbstractInterfaceSmoothing end -export StaircaseDNS, PeriodicStaircaseDNS, DNSModel, SDNS_simulation_setup +export StaircaseDNS, DNSModel, SDNS_simulation_setup export STStaircaseInitialConditions, StaircaseICs, SmoothSTStaircaseInitialConditions, STSingleInterfaceInitialConditions, SingleInterfaceICs, From 5c2af4d157771177e7fef91ed5af4dc37808299a Mon Sep 17 00:00:00 2001 From: jbisits Date: Sun, 24 Nov 2024 13:06:42 +1100 Subject: [PATCH 09/12] small tidy ups to examples before creating some validation type models to check everything looks ok --- examples/four_stair_staircase.jl | 2 +- examples/single_interface.jl | 5 ++--- examples/single_interface_periodic_tanh_background.jl | 7 +++---- 3 files changed, 6 insertions(+), 8 deletions(-) diff --git a/examples/four_stair_staircase.jl b/examples/four_stair_staircase.jl index 4284613..3f155b1 100644 --- a/examples/four_stair_staircase.jl +++ b/examples/four_stair_staircase.jl @@ -26,7 +26,7 @@ set_staircase_initial_conditions!(sdns) Δt = 1e-1 stop_time = 10 * 60 # seconds save_schedule = 10 # seconds -output_path = joinpath(@__DIR__, "output") +output_path = joinpath(@__DIR__, "output_staircase") simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!; output_path, Δt, save_schedule,) diff --git a/examples/single_interface.jl b/examples/single_interface.jl index 5727abc..6560808 100644 --- a/examples/single_interface.jl +++ b/examples/single_interface.jl @@ -21,9 +21,8 @@ sdns = StaircaseDNS(model_setup, interface_ics, velocity_noise) Δt = 1e-1 stop_time = 210 * 60 # seconds save_schedule = 10 # seconds -output_path = joinpath(@__DIR__, "output") -simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, - save_vertical_velocities!; +output_path = joinpath(@__DIR__, "output_single_interface") +simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, output_path, max_Δt = 5) ## Run diff --git a/examples/single_interface_periodic_tanh_background.jl b/examples/single_interface_periodic_tanh_background.jl index 4d195a7..162e30d 100644 --- a/examples/single_interface_periodic_tanh_background.jl +++ b/examples/single_interface_periodic_tanh_background.jl @@ -11,7 +11,7 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.56, 34.70] temperature = [-1.5, 0.5] -interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, interface_smoothing = Tanh) +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity) tracer_noise = TracerNoise(1e-6, 1e-6) ## setup model @@ -21,9 +21,8 @@ sdns = StaircaseDNS(model_setup, interface_ics, tracer_noise) Δt = 1e-1 stop_time = 4 * 60 * 60 # seconds save_schedule = 30 # seconds -output_path = joinpath(@__DIR__, "output_tanh_background") -simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, - save_vertical_velocities!; +output_path = joinpath(@__DIR__, "output_periodic") +simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, output_path, max_Δt = 5) ## Run From 7e456adb15039884bd552a0d411919572a843eff Mon Sep 17 00:00:00 2001 From: jbisits Date: Sun, 24 Nov 2024 14:39:15 +1100 Subject: [PATCH 10/12] Rearrange for validations + slightly nicer info --- ...ngle_interface_periodic_tanh_background.jl | 2 +- src/StaircaseShenanigans.jl | 4 +- src/interface_smoothing.jl | 8 ++- src/single_interfaces_initial_conditions.jl | 8 +-- src/staircase_background.jl | 5 +- ...le_interface_periodic_linear_background.jl | 14 ++--- ...ngle_interface_periodic_tanh_background.jl | 58 +++++++++++++++++++ 7 files changed, 81 insertions(+), 18 deletions(-) rename {examples => validation}/single_interface_periodic_linear_background.jl (73%) create mode 100644 validation/single_interface_periodic_tanh_background.jl diff --git a/examples/single_interface_periodic_tanh_background.jl b/examples/single_interface_periodic_tanh_background.jl index 162e30d..2cb5940 100644 --- a/examples/single_interface_periodic_tanh_background.jl +++ b/examples/single_interface_periodic_tanh_background.jl @@ -11,7 +11,7 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.56, 34.70] temperature = [-1.5, 0.5] -interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity) +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, background_state = BackgroundTanh()) tracer_noise = TracerNoise(1e-6, 1e-6) ## setup model diff --git a/src/StaircaseShenanigans.jl b/src/StaircaseShenanigans.jl index a143440..8950a30 100644 --- a/src/StaircaseShenanigans.jl +++ b/src/StaircaseShenanigans.jl @@ -32,9 +32,9 @@ export STStaircaseInitialConditions, StaircaseICs, SmoothSTStaircaseInitialCondi PeriodicSTSingleInterfaceInitialConditions, PeriodoicSingleInterfaceICs, set_staircase_initial_conditions! -export TanhInterfaceSmoothing, Tanh +export TanhInterfaceSmoothing, Tanh, NoSmoothing -export BackgroundTanh, BackgroundLinear, tanh_background, linear_background +export BackgroundTanh, BackgroundLinear, tanh_background, linear_background, NoBackground export AbstractNoise, VelocityNoise, TracerNoise diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl index aa42b14..44fb1d8 100644 --- a/src/interface_smoothing.jl +++ b/src/interface_smoothing.jl @@ -1,4 +1,3 @@ - """ struct TanhInterfaceSmoothing Container to set a hyperbolic tangent over a single interface as an initial condition. @@ -18,5 +17,8 @@ end const Tanh = TanhInterfaceSmoothing{T} where {T} @inline (p::TanhInterfaceSmoothing)(x, y, z) = p.Cₗ - 0.5 * p.ΔC * (1 + tanh(p.D * (z - p.z_interface) / p.z_range)) -Base.summary(p::Type{<:TanhInterfaceSmoothing}) = "tanh smoothing." -Base.summary(p::Nothing) = "$p" +Base.summary(p::Type{<:TanhInterfaceSmoothing}) = "tanh smoothing" + +struct NoSmoothing <: AbstractInterfaceSmoothing end +NoSmoothing() = nothing +Base.summary(ns::NoSmoothing) = "no interface smoothing" diff --git a/src/single_interfaces_initial_conditions.jl b/src/single_interfaces_initial_conditions.jl index cc87bb7..1ac47bb 100644 --- a/src/single_interfaces_initial_conditions.jl +++ b/src/single_interfaces_initial_conditions.jl @@ -42,7 +42,7 @@ function STSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, dept end const SingleInterfaceICs = STSingleInterfaceInitialConditions # alias -Base.summary(ics::SingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface)" +Base.summary(ics::SingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface) on bounded z domain" """ struct PeriodicSTSingleInterfaceInitialConditions @@ -65,8 +65,8 @@ struct PeriodicSTSingleInterfaceInitialConditions{T, A, IS, BF} <: AbstractSingl end function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfState, depth_of_interface, salinity, temperature; - interface_smoothing = nothing, - background_state = nothing) + interface_smoothing = NoSmoothing(), + background_state = NoBackground()) R_ρ = compute_R_ρ(salinity, temperature, depth_of_interface, eos) @@ -75,7 +75,7 @@ function PeriodicSTSingleInterfaceInitialConditions(eos::BoussinesqEquationOfSta end const PeriodoicSingleInterfaceICs = PeriodicSTSingleInterfaceInitialConditions # alias -Base.summary(ics::PeriodoicSingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface) on triply periodic domain with $(summary(ics.background_state)) state" +Base.summary(ics::PeriodoicSingleInterfaceICs) = "Single S-T interface at z = $(ics.depth_of_interface) on triply periodic domain" function compute_R_ρ(salinity, temperature, depth_of_interfaces::Number, eos) diff --git a/src/staircase_background.jl b/src/staircase_background.jl index 8305358..a1959eb 100644 --- a/src/staircase_background.jl +++ b/src/staircase_background.jl @@ -38,9 +38,12 @@ function Base.show(io::IO, bf::AbstractBackgroundFunction) print(io, "┗━ parameters: $(bf.parameters)") end end + +struct NoBackgroundFunction <: AbstractBackgroundFunction end +const NoBackground = NoBackgroundFunction() = nothing Base.summary(bt::BackgroundTanh) = "$(bt.func)" Base.summary(bl::BackgroundLinear) = "$(bl.func)" -Base.summary(bn::Nothing) = "$(bn)" +Base.summary(bn::NoBackgroundFunction) = "no background" """ function S_and_T_background_fields(initial_conditions) diff --git a/examples/single_interface_periodic_linear_background.jl b/validation/single_interface_periodic_linear_background.jl similarity index 73% rename from examples/single_interface_periodic_linear_background.jl rename to validation/single_interface_periodic_linear_background.jl index cab3fa6..e4aa597 100644 --- a/examples/single_interface_periodic_linear_background.jl +++ b/validation/single_interface_periodic_linear_background.jl @@ -11,25 +11,25 @@ model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) depth_of_interface = -0.5 salinity = [34.56, 34.70] temperature = [-1.5, 0.5] -interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, - BackgroundLinear()) -tracer_noise = TracerNoise(1e-6, 1e-6) + +#### Background linear state with tanh smoothing +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, + background_state = BackgroundLinear(), + interface_smoothing = Tanh) ## setup model -sdns = StaircaseDNS(model_setup, interface_ics, tracer_noise) +sdns = StaircaseDNS(model_setup, interface_ics, nothing) ## Build simulation Δt = 1e-1 stop_time = 4 * 60 * 60 # seconds save_schedule = 30 # seconds -output_path = joinpath(@__DIR__, "output_linear_background") +output_path = joinpath(@__DIR__, "output_periodic_tanh_interface_linear_background") simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, output_path, max_Δt = 5) ## Run run!(simulation) - -## Compute density ratio compute_R_ρ!(simulation.output_writers[:computed_output].filepath, simulation.output_writers[:tracers].filepath, eos) diff --git a/validation/single_interface_periodic_tanh_background.jl b/validation/single_interface_periodic_tanh_background.jl new file mode 100644 index 0000000..9ff5fa5 --- /dev/null +++ b/validation/single_interface_periodic_tanh_background.jl @@ -0,0 +1,58 @@ +using StaircaseShenanigans + +architecture = CPU() +diffusivities = (ν = 1e-5, κ = (S = 1e-7, T = 1e-5)) +domain_extent = (Lx = 0.1, Ly = 0.1, Lz = -1.0) +resolution = (Nx = 5, Ny = 5, Nz = 50) +eos = CustomLinearEquationOfState(-0.5, 34.6) +model_setup = (;architecture, diffusivities, domain_extent, resolution, eos) + +## Initial conditions +depth_of_interface = -0.5 +salinity = [34.56, 34.70] +temperature = [-1.5, 0.5] + +#### Interface smoothing no background state +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, + interface_smoothing = Tanh) + +## setup model +sdns = StaircaseDNS(model_setup, interface_ics, nothing) + +## Build simulation +Δt = 1e-1 +stop_time = 4 * 60 * 60 # seconds +save_schedule = 30 # seconds +output_path = joinpath(@__DIR__, "output_periodic_tanh_interface") +simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; + Δt, save_schedule, + output_path, max_Δt = 5) +## Run +run!(simulation) + +# Compute density ratio +compute_R_ρ!(simulation.output_writers[:computed_output].filepath, + simulation.output_writers[:tracers].filepath, eos) + +#### Background state with tanh interface smoothing +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, + background_state = BackgroundTanh(), + interface_smoothing = Tanh) + +## setup model +sdns = StaircaseDNS(model_setup, interface_ics, nothing) + +## Build simulation +Δt = 1e-1 +stop_time = 4 * 60 * 60 # seconds +save_schedule = 30 # seconds +output_path = joinpath(@__DIR__, "output_periodic_tanh_interface_tanh_background") +simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; + Δt, save_schedule, + output_path, max_Δt = 5) +## Run +run!(simulation) + +# Compute density ratio +compute_R_ρ!(simulation.output_writers[:computed_output].filepath, + simulation.output_writers[:tracers].filepath, eos) From 2b0cd9572d997542e5e399487b3a4bd85108e66a Mon Sep 17 00:00:00 2001 From: jbisits Date: Sun, 24 Nov 2024 14:40:10 +1100 Subject: [PATCH 11/12] Remove output in path name --- validation/single_interface_periodic_linear_background.jl | 2 +- validation/single_interface_periodic_tanh_background.jl | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/validation/single_interface_periodic_linear_background.jl b/validation/single_interface_periodic_linear_background.jl index e4aa597..81f5ffd 100644 --- a/validation/single_interface_periodic_linear_background.jl +++ b/validation/single_interface_periodic_linear_background.jl @@ -24,7 +24,7 @@ sdns = StaircaseDNS(model_setup, interface_ics, nothing) Δt = 1e-1 stop_time = 4 * 60 * 60 # seconds save_schedule = 30 # seconds -output_path = joinpath(@__DIR__, "output_periodic_tanh_interface_linear_background") +output_path = joinpath(@__DIR__, "periodic_tanh_interface_linear_background") simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, diff --git a/validation/single_interface_periodic_tanh_background.jl b/validation/single_interface_periodic_tanh_background.jl index 9ff5fa5..67bf684 100644 --- a/validation/single_interface_periodic_tanh_background.jl +++ b/validation/single_interface_periodic_tanh_background.jl @@ -23,7 +23,7 @@ sdns = StaircaseDNS(model_setup, interface_ics, nothing) Δt = 1e-1 stop_time = 4 * 60 * 60 # seconds save_schedule = 30 # seconds -output_path = joinpath(@__DIR__, "output_periodic_tanh_interface") +output_path = joinpath(@__DIR__, "periodic_tanh_interface") simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, output_path, max_Δt = 5) @@ -46,7 +46,7 @@ sdns = StaircaseDNS(model_setup, interface_ics, nothing) Δt = 1e-1 stop_time = 4 * 60 * 60 # seconds save_schedule = 30 # seconds -output_path = joinpath(@__DIR__, "output_periodic_tanh_interface_tanh_background") +output_path = joinpath(@__DIR__, "periodic_tanh_interface_tanh_background") simulation = SDNS_simulation_setup(sdns, stop_time, save_computed_output!, save_vertical_velocities!; Δt, save_schedule, output_path, max_Δt = 5) From e7558a8c62dad277c7e81cc2093769c86bd34197 Mon Sep 17 00:00:00 2001 From: jbisits Date: Sun, 24 Nov 2024 14:53:42 +1100 Subject: [PATCH 12/12] Fix validation setup --- src/interface_smoothing.jl | 6 +++--- validation/single_interface_periodic_linear_background.jl | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/src/interface_smoothing.jl b/src/interface_smoothing.jl index 44fb1d8..dfb451d 100644 --- a/src/interface_smoothing.jl +++ b/src/interface_smoothing.jl @@ -19,6 +19,6 @@ const Tanh = TanhInterfaceSmoothing{T} where {T} Base.summary(p::Type{<:TanhInterfaceSmoothing}) = "tanh smoothing" -struct NoSmoothing <: AbstractInterfaceSmoothing end -NoSmoothing() = nothing -Base.summary(ns::NoSmoothing) = "no interface smoothing" +struct NoInterfaceSmoothing <: AbstractInterfaceSmoothing end +const NoSmoothing = NoInterfaceSmoothing() = nothing +Base.summary(ns::NoInterfaceSmoothing) = "no interface smoothing" diff --git a/validation/single_interface_periodic_linear_background.jl b/validation/single_interface_periodic_linear_background.jl index 81f5ffd..898795c 100644 --- a/validation/single_interface_periodic_linear_background.jl +++ b/validation/single_interface_periodic_linear_background.jl @@ -13,7 +13,7 @@ salinity = [34.56, 34.70] temperature = [-1.5, 0.5] #### Background linear state with tanh smoothing -interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, +interface_ics = PeriodoicSingleInterfaceICs(eos, depth_of_interface, salinity, temperature, background_state = BackgroundLinear(), interface_smoothing = Tanh)