Skip to content

Commit

Permalink
Merge pull request #9 from jagoosw/jsw/change-to-discrete-bcs
Browse files Browse the repository at this point in the history
(0.5.3) Changes everything to be discrete boundary conditions to make compatible with flattend dimensions
  • Loading branch information
jagoosw authored Aug 9, 2024
2 parents 183dba8 + e89f271 commit 75848da
Show file tree
Hide file tree
Showing 11 changed files with 192 additions and 235 deletions.
233 changes: 50 additions & 183 deletions Manifest.toml

Large diffs are not rendered by default.

6 changes: 3 additions & 3 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "Walrus"
uuid = "bec5164f-1c7a-4c32-8768-be9354254d6a"
authors = ["Jago Stong-Wright <[email protected]>"]
version = "0.5.2"
version = "0.5.3"

[deps]
Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e"
Expand All @@ -14,8 +14,8 @@ Roots = "f2b01f46-fcfa-551c-844a-d8ac1e96c665"
Adapt = "3.7"
CUDA = "5"
KernelAbstractions = "0.9.22"
OceanBioME = "0.10.3"
Oceananigans = "0.91.3"
OceanBioME = "0.10.5"
Oceananigans = "0.91.7"
Roots = "2"
Test = "1.10"
julia = "1.10"
Expand Down
13 changes: 1 addition & 12 deletions src/Walrus.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,18 +18,7 @@ using Adapt: adapt

import Adapt: adapt_structure

# This is how we will standarise function vs values, I will have to think of a way to deal with discrete vs continuous at some point
struct ReturnValue{FT}
value :: FT
end

adapt_structure(to, rv::ReturnValue) = ReturnValue(adapt(to, rv.value))

@inline (return_value::ReturnValue)(args...) = return_value.value

display_input(return_value::ReturnValue) = return_value.value
display_input(::Function) = "Function"

include("get_value.jl")
include("interpolations.jl")
include("wall_model.jl")
include("radiative_transfer/radiative_transfer.jl")
Expand Down
35 changes: 35 additions & 0 deletions src/get_value.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,35 @@
using Oceananigans.Fields: AbstractField, Center
using Oceananigans.Grids: xnode, ynode

# fallback
@inline get_value(f, args...) = f

struct ContinuousSurfaceFunction{F}
func :: F
end

@inline function get_value(f::ContinuousSurfaceFunction, i, j, grid, clock, args...)
t = clock.time

x = xnode(i, j, grid.Nz, grid, Center(), Center(), Center())
y = ynode(i, j, grid.Nz, grid, Center(), Center(), Center())

return f.func(x, y, t, args...)
end

@inline get_value(f::ContinuousSurfaceFunction, ::Nothing, ::Nothing, t, args...) = f(0, 0, t)

struct DiscreteSurfaceFuncton{F}
func :: F
end

@inline get_value(f::DiscreteSurfaceFuncton, i, j, grid, clock, args...) = f.func(i, j, grid, clock, args...)

# GPU compatible, mainly for fields
@inline get_value(f::AbstractArray{<:Any, 2}, i, j, grid, clock, args...) = @inbounds f[i, j]
@inline get_value(f::AbstractArray{<:Any, 3}, i, j, grid, clock, args...) = @inbounds f[i, j, grid.Nz]

# fallback
normalise_surface_function(f; kwargs...) = f

normalise_surface_function(f::Function; discrete_form = false) = discrete_form ? DiscreteSurfaceFuncton(f) : ContinuousSurfaceFunction(f)
3 changes: 2 additions & 1 deletion src/radiative_transfer/radiative_transfer.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,6 @@ module RadiativeTransfer
export HomogeneousBodyHeating, PARModelHeating

using KernelAbstractions
using Walrus: ReturnValue

using Adapt: adapt
using Oceananigans.Architectures: architecture
Expand All @@ -20,6 +19,8 @@ using Oceananigans.Utils: launch!

using KernelAbstractions.Extras: @unroll

using Walrus: normalise_surface_function

import Oceananigans.Biogeochemistry: update_biogeochemical_state!, update_tendencies!, AbstractBiogeochemistry

import Adapt: adapt_structure
Expand Down
32 changes: 19 additions & 13 deletions src/surface_heating.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@ using Adapt: adapt

using Oceananigans.BoundaryConditions: FluxBoundaryCondition

using Walrus: ReturnValue, display_input
using Walrus: get_value, normalise_surface_function
using Walrus.WindStressModel: WindStress,
LogarithmicNeutralWind,
find_velocity_roughness_length
Expand Down Expand Up @@ -119,10 +119,10 @@ julia> using Walrus
julia> using Oceananigans
julia> wind_stress = WindStress(; reference_wind_speed = 0., reference_wind_direction = 90.)
(::WindStress{Walrus.ReturnValue{Float64}, Walrus.ReturnValue{Float64}, LogarithmicNeutralWind{Float64, Nothing}, Float64}) (generic function with 2 methods)
(::WindStress{Float64, Float64, LogarithmicNeutralWind{Float64, Nothing}, Float64}) (generic function with 2 methods)
julia> surface_heat_exchange = SurfaceHeatExchange(; wind_stress)
(::SurfaceHeatExchange{WindStress{Walrus.ReturnValue{Float64}, Walrus.ReturnValue{Float64}, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Walrus.ReturnValue{Int64}, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64}) (generic function with 1 method)
(::SurfaceHeatExchange{WindStress{Float64, Float64, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Int64, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64}) (generic function with 1 method)
julia> boundary_conditions = (; T = FieldBoundaryConditions(top = FluxBoundaryCondition(surface_heat_exchange, field_dependencies = (:T, :u, :v))))
(T = Oceananigans.FieldBoundaryConditions, with boundary conditions
Expand All @@ -131,7 +131,7 @@ julia> boundary_conditions = (; T = FieldBoundaryConditions(top = FluxBoundaryCo
├── south: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── north: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── bottom: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing)
├── top: FluxBoundaryCondition: ContinuousBoundaryFunction (::SurfaceHeatExchange{WindStress{Walrus.ReturnValue{Float64}, Walrus.ReturnValue{Float64}, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Walrus.ReturnValue{Int64}, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64}) at (Nothing, Nothing, Nothing)
├── top: FluxBoundaryCondition: ContinuousBoundaryFunction (::SurfaceHeatExchange{WindStress{Float64, Float64, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Int64, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64}) at (Nothing, Nothing, Nothing)
└── immersed: DefaultBoundaryCondition (FluxBoundaryCondition: Nothing),)
```
Expand All @@ -147,7 +147,7 @@ function SurfaceHeatExchange(; wind_stress,
air_water_mixing_ratio = 0.001, # kg / kgs
stephan_boltzman_constant = 5.670374419e-8) # W / K⁴

isa(air_temperature, Function) || (air_temperature = ReturnValue(air_temperature))
air_temperature = normalise_surface_function(air_temperature)

return SurfaceHeatExchange(wind_stress, air_temperature,
latent_heat_vaporisation, vapour_pressure,
Expand Down Expand Up @@ -190,10 +190,10 @@ Example
julia> using Walrus
julia> wind_stress = WindStress(; reference_wind_speed = 0., reference_wind_direction = 90.)
(::WindStress{Walrus.ReturnValue{Float64}, Walrus.ReturnValue{Float64}, LogarithmicNeutralWind{Float64, Nothing}, Float64}) (generic function with 2 methods)
(::WindStress{Float64, Float64, LogarithmicNeutralWind{Float64, Nothing}, Float64}) (generic function with 2 methods)
julia> surface_heat_exchange = SurfaceHeatExchangeBoundaryCondition(; wind_stress)
FluxBoundaryCondition: ContinuousBoundaryFunction (::SurfaceHeatExchange{WindStress{Walrus.ReturnValue{Float64}, Walrus.ReturnValue{Float64}, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Walrus.ReturnValue{Int64}, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64}) at (Nothing, Nothing, Nothing)
FluxBoundaryCondition: DiscreteBoundaryFunction with (::SurfaceHeatExchange{WindStress{Float64, Float64, LogarithmicNeutralWind{Float64, Nothing}, Float64}, Int64, Walrus.SurfaceHeatingModel.EmpiricalLatentHeatVaporisation{Float64}, Walrus.SurfaceHeatingModel.AugustRocheMagnusVapourPressure{Float64}, Float64})
```
"""
Expand All @@ -208,15 +208,15 @@ function SurfaceHeatExchangeBoundaryCondition(; wind_stress,
air_water_mixing_ratio = 0.001, # kg / kg
stephan_boltzman_constant = 5.670374419e-8) # W / K⁴

isa(air_temperature, Function) || (air_temperature = ReturnValue(air_temperature))
air_temperature = normalise_surface_function(air_temperature)

surface_heat_exchange = SurfaceHeatExchange(wind_stress, air_temperature,
latent_heat_vaporisation, vapour_pressure,
water_specific_heat_capacity, water_density,
air_specific_heat_capacity, air_density, air_water_mixing_ratio,
stephan_boltzman_constant)

return FluxBoundaryCondition(surface_heat_exchange, field_dependencies = (:T, :u, :v))
return FluxBoundaryCondition(surface_heat_exchange, discrete_form=true)
end

@inline function (drag_coefficient::LogarithmicNeutralWind, wind_speed)
Expand Down Expand Up @@ -268,7 +268,7 @@ end
##### Cʰ parameterisations
#####

@inline function (interface::SurfaceHeatExchange)(x, y, t, T, u, v)
@inline function (interface::SurfaceHeatExchange)(i, j, grid, clock, model_fields)
σ = interface.stephan_boltzman_constant
ρᵃ = interface.air_density
cₚᵃ = interface.air_specific_heat_capacity
Expand All @@ -277,10 +277,16 @@ end
ρʷ = interface.water_density
cₚʷ = interface.water_specific_heat_capacity

air_temperature = interface.air_temperature(x, y, t)
t = clock.time

wind_speed = interface.wind_stress.reference_wind_speed(x, y, t)
wind_direction = interface.wind_stress.reference_wind_direction(x, y, t)
u = @inbounds model_fields.u[i, j, grid.Nz]
v = @inbounds model_fields.v[i, j, grid.Nz]
T = @inbounds model_fields.T[i, j, grid.Nz]

air_temperature = get_value(interface.air_temperature, i, j, grid, clock)

wind_speed = get_value(interface.wind_stress.reference_wind_speed, i, j, grid, clock)
wind_direction = get_value(interface.wind_stress.reference_wind_direction, i, j, grid, clock)

= - wind_speed * sind(wind_direction)
= - wind_speed * cosd(wind_direction)
Expand Down
1 change: 1 addition & 0 deletions src/wall_model.jl
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@ using Oceananigans.Grids: znode

using Oceananigans.Architectures: on_architecture, CPU, architecture

using Walrus: get_value, normalise_surface_function
using Walrus.Interpolations: SimpleInterpolation

import Adapt: adapt_structure
Expand Down
3 changes: 2 additions & 1 deletion src/wind_driven_stokes.jl
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ using Oceananigans.Architectures: on_architecture, CPU, architecture
using Oceananigans.BuoyancyModels: g_Earth
using Oceananigans.StokesDrifts: UniformStokesDrift

using Walrus: get_value
using Walrus.Interpolations: SimpleInterpolation

import Adapt: adapt_structure
Expand Down Expand Up @@ -107,7 +108,7 @@ end
@inline function surface_drift_velocity(z, t, params)
wind = params.wind

= abs(wind.reference_wind_speed(0, 0, t))
= abs(get_value(wind.reference_wind_speed, nothing, nothing, nothing, t))

dc = wind.drag_coefficient

Expand Down
Loading

2 comments on commit 75848da

@jagoosw
Copy link
Owner Author

@jagoosw jagoosw commented on 75848da Aug 9, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@JuliaRegistrator
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Registration pull request created: JuliaRegistries/General/112776

Tip: Release Notes

Did you know you can add release notes too? Just add markdown formatted text underneath the comment after the text
"Release notes:" and it will be added to the registry PR, and if TagBot is installed it will also be added to the
release that TagBot creates. i.e.

@JuliaRegistrator register

Release notes:

## Breaking changes

- blah

To add them here just re-invoke and the PR will be updated.

Tagging

After the above pull request is merged, it is recommended that a tag is created on this repository for the registered package version.

This will be done automatically if the Julia TagBot GitHub Action is installed, or can be done manually through the github interface, or via:

git tag -a v0.5.3 -m "<description of version>" 75848daf5e22ab8b161281aed209f8ccb1251d9f
git push origin v0.5.3

Also, note the warning: Version 0.5.3 skips over 0.5.1
This can be safely ignored. However, if you want to fix this you can do so. Call register() again after making the fix. This will update the Pull request.

Please sign in to comment.