Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Question: What is maximum_sinking_velocity(bgc) expected behaviour? #229

Open
radka-j opened this issue Nov 8, 2024 · 4 comments
Open
Labels
bug Something isn't working

Comments

@radka-j
Copy link

radka-j commented Nov 8, 2024

Hi!

I'm working on implementing a biogeochemistry model using OceanBioME and have found it incredibly useful, so thank you!

I have a question about themaximum_sinking_velocity(bgc) method. I have noticed that NPZD and LOBSTER both define it but not PISCES. it is also implemented in the main Biogeochemistry class where it is set to return 0 and this is the method that seems to be called in the code.

Additionally, if I use the example in the README and try to call this method, it works if called on the biogeochemistry class (returning a 0) but gives an error if called on biogeochemistry.underlying_biogeochemistry (which is NPZD in this case):

julia> OceanBioME.maximum_sinking_velocity(biogeochemistry)
0.0
julia> OceanBioME.maximum_sinking_velocity(biogeochemistry.underlying_biogeochemistry)
ERROR: type Field has no field w

Is the expected behaviour that the returned value is always 0? What is this used for?

@jagoosw
Copy link
Collaborator

jagoosw commented Nov 8, 2024

Hi @radka-j,

The intention of this function was always to return the maximum sinking velocity that a bgc model imposed, but it seems like it was never tested and then forgotten. They could then be used by the code you linked in a timestep wizard-like:

wizard = TimeStepWizard(cfl = 0.15, diffusive_cfl = 0.15, max_change = 2.0, min_change = 0.5, cell_diffusion_timescale = column_diffusion_timescale, cell_advection_timescale = column_advection_timescale)

maximum_sinking_velocity(bgc) is hitting the fallback:

@inline maximum_sinking_velocity(bgc) = 0.0

because I never implemented:

maximum_sinking_velocity(bgc::CompleteBiogeochemistry) = maximum_sinking_velocity(bgc.underlying_biogeochemistry)

and at some point I forgot to update the NPZD and LOBSTER methods which should be like:

 @inline maximum_sinking_velocity(bgc::LOBSTER) = maximum(abs, bgc.sinking_velocities.bPOM) 

not

@inline maximum_sinking_velocity(bgc::LOBSTER) = maximum(abs, bgc.sinking_velocities.bPOM.w)

I can fix this at some point, or if you'd like to use it feel free to PR the changes if I don't get round to it soon enough!

@jagoosw jagoosw added the bug Something isn't working label Nov 8, 2024
@glwagner
Copy link
Collaborator

glwagner commented Nov 8, 2024

The intention of this function was always to return the maximum sinking velocity that a bgc model imposed, but it seems like it was never tested and then forgotten. They could then be used by the code you linked in a timestep wizard-like:

That'd be appropriate to implement in Oceananigans if it is extending TimeStepWizard.

@jagoosw
Copy link
Collaborator

jagoosw commented Nov 8, 2024

The intention of this function was always to return the maximum sinking velocity that a bgc model imposed, but it seems like it was never tested and then forgotten. They could then be used by the code you linked in a timestep wizard-like:

That'd be appropriate to implement in Oceananigans if it is extending TimeStepWizard.

That's probably a much better approach since it should just be able to return something like maximum(map(n->maximum(drift_velocities[n].w), 1:length(drift_velocities)), I guess this never actually needed to be defined for all models but we can just add methods to make it more efficient if we want to

@glwagner
Copy link
Collaborator

glwagner commented Nov 8, 2024

I think the right solution is to use total_velocities to compute the CFL here:

https://github.com/CliMA/Oceananigans.jl/blob/5dd01ddfe73dd03515b52e6499b27c1b61c71364/src/Models/NonhydrostaticModels/nonhydrostatic_tendency_kernel_functions.jl#L259-L261

This is guaranteed to work universally I think and wouldn't require maintaining any code in OceanBioME specifically?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug Something isn't working
Projects
None yet
Development

No branches or pull requests

3 participants