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

New and more correct turbulence closures #383

Merged
merged 31 commits into from
Aug 29, 2019
Merged

Conversation

glwagner
Copy link
Member

@glwagner glwagner commented Aug 29, 2019

This PR makes a few major contributions:

  • Fixes a bug in Anisotropic Minimum Dissipation (AMD) type closures that led to (very) wrong subgrid-scale diffusivities.
  • Adds a new flavor of AMD proposed by Verstappen (2018) and described in detail by Vreugdenhil and Taylor (2018). This new flavor of AMD, dubbed the "Verstappen" AMD, is now the default name bound to AnisotropicMinimumDissipation
  • The closure formerly known as AnisotropicMinimumDissipation is now called RozemaAnisotropicMinimumDissipation (it was first described by Rozema et al. (2015)). It is not exported, but can be accessed via Oceananigans.RozemaAnisotropicMinimumDissipation.
  • There is new Smagorinsky closure variant called BlasiusSmagorinsky which is untested and unexported. It can be accessed via Oceananigans.BlasiusSmagorinsky.

…roduces more accurate approximation for near-boundary diffusivity in amd
Merge branch 'master' into new-closures
@glwagner glwagner added science 🌊 Sometimes addictive, sometimes allergenic feature 🌟 Something new and shiny bug 🐞 Even a perfect program still has bugs labels Aug 29, 2019
@codecov
Copy link

codecov bot commented Aug 29, 2019

Codecov Report

Merging #383 into master will decrease coverage by 8.21%.
The diff coverage is 52.87%.

Impacted file tree graph

@@            Coverage Diff             @@
##           master     #383      +/-   ##
==========================================
- Coverage   78.37%   70.15%   -8.22%     
==========================================
  Files          22       23       +1     
  Lines        1244     1528     +284     
==========================================
+ Hits          975     1072      +97     
- Misses        269      456     +187
Impacted Files Coverage Δ
src/boundary_conditions.jl 64.28% <ø> (-5.72%) ⬇️
src/models.jl 92.68% <ø> (ø) ⬆️
src/turbulence_closures/TurbulenceClosures.jl 56.25% <0%> (-27.09%) ⬇️
src/turbulence_closures/closure_operators.jl 64.44% <0%> (-0.73%) ⬇️
...rbulence_closures/constant_diffusivity_closures.jl 100% <100%> (+40.54%) ⬆️
src/time_steppers.jl 75.97% <100%> (-0.61%) ⬇️
...closures/rozema_anisotropic_minimum_dissipation.jl 23.95% <23.95%> (ø)
...ures/verstappen_anisotropic_minimum_dissipation.jl 53.33% <53.33%> (ø)
src/turbulence_closures/smagorinsky.jl 67.79% <67.79%> (ø)
...c/turbulence_closures/velocity_tracer_gradients.jl 69.72% <69.72%> (ø)
... and 10 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 5a12ab6...7e7449e. Read the comment docs.

@codecov
Copy link

codecov bot commented Aug 29, 2019

Codecov Report

Merging #383 into master will decrease coverage by 8.17%.
The diff coverage is 52.87%.

Impacted file tree graph

@@            Coverage Diff            @@
##           master    #383      +/-   ##
=========================================
- Coverage   78.37%   70.2%   -8.18%     
=========================================
  Files          22      23       +1     
  Lines        1244    1527     +283     
=========================================
+ Hits          975    1072      +97     
- Misses        269     455     +186
Impacted Files Coverage Δ
src/boundary_conditions.jl 64.28% <ø> (-5.72%) ⬇️
src/models.jl 92.68% <ø> (ø) ⬆️
src/turbulence_closures/TurbulenceClosures.jl 56.25% <0%> (-27.09%) ⬇️
src/turbulence_closures/closure_operators.jl 64.44% <0%> (-0.73%) ⬇️
...rbulence_closures/constant_diffusivity_closures.jl 100% <100%> (+40.54%) ⬆️
src/time_steppers.jl 75.97% <100%> (-0.61%) ⬇️
...closures/rozema_anisotropic_minimum_dissipation.jl 23.95% <23.95%> (ø)
...ures/verstappen_anisotropic_minimum_dissipation.jl 53.33% <53.33%> (ø)
src/turbulence_closures/smagorinsky.jl 67.79% <67.79%> (ø)
...c/turbulence_closures/velocity_tracer_gradients.jl 69.72% <69.72%> (ø)
... and 10 more

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 5a12ab6...be136c3. Read the comment docs.

Copy link
Member

@ali-ramadhan ali-ramadhan left a comment

Choose a reason for hiding this comment

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

Nice! Just left some minor comments but otherwise looks good to merge. Smag test isn't skipped anymore so we should be back to 0 broken tests now.

We should release v0.10.0 once this is merged.

Just some comments:

Do we have an idea of what the Poincaré constant C should be in AMD? Or maybe it's something that needs to be investigated depending on the numerical methods used? Does the advection scheme and temporal discretization matter or is it mainly the numerical method (e.g. finite volume vs. spectral)?

In particular, I notice that Rozema AMD uses C=0.33 and Versteppen AMD uses C=1/12 which are pretty different.

Presumably the same question applies to Smag, although @jm-c would know what constant to use for constant Smag combined with our numerical methods as it has been used in the MITgcm.

General concerns:

  1. We've already discussed whether we want to maintain multiple versions of AMD And Smag but maybe what to keep will become clear with time. Maybe they have different use applications so have multiple versions could be a nice feature? Presumably a lot of code is shared between the two AMD closures which, if merged, would make maintaining both easy.
  2. Still not sure why TurbulenceClosures uses it's own differentiation and interpolation operators when we have Oceananigans.Operators. They'll have to be rewritten in finite volume form when we move to a stretched grid (PRs Finite volume operators through multiple dispatch #283 ) so maybe that'll be an opportunity to unify both sets of operators.

src/time_steppers.jl Outdated Show resolved Hide resolved
src/turbulence_closures/TurbulenceClosures.jl Show resolved Hide resolved
# at 20ᵒC and 35 psu, according to Sharqawy et al., "Thermophysical
# properties of seawater: A review of existing correlations and data" (2010).
const ν₀ = 1.05e-6
const κ₀ = 1.46e-7
Copy link
Member

Choose a reason for hiding this comment

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

Would a subscript m be better to denote that these are molecular values?

const νₘ = 1.05e-6
const κₘ = 1.46e-7

Copy link
Member Author

Choose a reason for hiding this comment

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

We can do that --- I would keep the 0, however, to emphasize that these values are chosen for a particular 'state', e.g. a temperature of 20C and a salinity of 35 psu. If you agree with this logic I will make these changes.

Copy link
Member

Choose a reason for hiding this comment

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

Ah ok if 0 refers to the reference state then that makes sense. I might include a comment explaining this but not really necessary.

Copy link
Member Author

Choose a reason for hiding this comment

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

@ali-ramadhan I have included the comment

# Approximate viscosities and thermal diffusivities for seawater
# at 20ᵒC and 35 psu, according to Sharqawy et al., "Thermophysical 
# properties of seawater: A review of existing correlations and data" (2010).

What do you think?

Copy link
Member

Choose a reason for hiding this comment

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

Yeah I saw that. I guess the 0 subscript was bothering me because it could imply many things. I might be more explicit about 0 implying reference state, but this is a minor comment. How many people are going to read through that file?

#####

#####
##### the heinous
Copy link
Member

Choose a reason for hiding this comment

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

Yup

@glwagner
Copy link
Member Author

glwagner commented Aug 29, 2019

@ali-ramadhan the constants have slightly different meaning, which creates confusion --- let me explain:

  • In the Rozema version of AMD, the 'filter width' δ is defined as the grid spacing. Thus on an isotropic grid, the eddy viscosity becomes

ν = -C_roz * Δ^2 * r / q,

where r and q are defined in equations (9) and (11) in Rozema et al. (2015)

  • In the Verstappen version of AMD (and as described in Vreugdenhil and Taylor 2018), the 'filter width' δ is defined as twice the grid spacing for finite differences (and 3/2 for spectral discretization). Thus on an isotropic grid and with finite differences the Verstappen version of AMD reduces to

ν = -C_ver * δ^2 * r / q = - 4 * C_ver * Δ^2 * r / q

Thus in this case C_ver = 1/12 is equivalent to C_roz = 4/12 = 0.33.

On anisotropic grids, things can be slightly different.

I'm open to changing the implementation, but I'm not sure what the best route is. Possibly we will phase out the Rozema version... ?

@glwagner
Copy link
Member Author

@jm-c would know what constant to use for constant Smag combined with our numerical methods as it has been used in the MITgcm.

We were told that Smagorinsky has not been tested in MITgcm so I'm not sure that's true.

Still not sure why TurbulenceClosures uses it's own differentiation and interpolation operators when we have Oceananigans.Operators.

Migrating to the other operators would have taken some time and effort and possibly involved bugs so I chose to keep the ones I was familiar with. I find the new operators a little more difficult to read than the ones we currently use, and I didn't want to have to worry about that when implementing complex new functions involved with the VerstappenAnisotropicMinimumDissipation.

@ali-ramadhan
Copy link
Member

Migrating to the other operators would have taken some time and effort and possibly involved bugs so I chose to keep the ones I was familiar with. I find the new operators a little more difficult to read than the ones we currently use, and I didn't want to have to worry about that when implementing complex new functions involved with the VerstappenAnisotropicMinimumDissipation.

Fair enough, I realize this was a decision made back when AMD was first implemented. Hopefully we'll find a way to merge the two in the future.

@ali-ramadhan
Copy link
Member

Thanks for explaining why the C values are different, didn't realize that.

Not sure where it comes from but shouldn't it be C_roz = 4/24 = 0.167 or C_roz = 8/24 = 0.33?

I'm open to changing the implementation, but I'm not sure what the best route is. Possibly we will phase out the Rozema version... ?

Not sure how much code overlaps between the two, but if both are deemed important enough to maintain, then the shared functionality can go into amd_core.jl or something. Although probably a bad idea to do this refactor before each AMD closure has it's own regression test.

@glwagner
Copy link
Member Author

Not sure where it comes from but shouldn't it be C_roz = 4/24 = 0.167 or C_roz = 8/24 = 0.33?

Typo. 4 * 1/12 = 4/12 = 0.33.

@ali-ramadhan
Copy link
Member

Looks good to merge!

@glwagner
Copy link
Member Author

Not sure how much code overlaps between the two, but if both are deemed important enough to maintain, then the shared functionality can go into amd_core.jl or something.

I considered this when I started developing the Verstappen closure. Almost none of the code is shared, which is why there are two files, unlike the Smagorinsky closures which have a lot of shared code. In principle it might be possible to refactor both implementations to share (a bit) more code but that would be a lot of work, and would also probably make both implementations more complex...

@ali-ramadhan
Copy link
Member

I considered this when I started developing the Verstappen closure. Almost none of the code is shared, which is why there are two files, unlike the Smagorinsky closures which have a lot of shared code. In principle it might be possible to refactor both implementations to share (a bit) more code but that would be a lot of work, and would also probably make both implementations more complex...

Didn't realize they were that different. In that case, maybe makes more sense to keep two separate files. More complex implementation of AMD sounds like a bad idea.

@glwagner
Copy link
Member Author

Didn't realize they were that different.

The difference between them is that both derivatives and velocity fields are normalized by Verstappen, while only derivatives are normalized by Rozema. To implement the Verstappen functions, I redefined functions that compute velocity gradients and tracer gradients by prefixing them with norm_. In the Rozema version, the normalization of derivatives is easily factored out by dividing the largest expression (for 27 terms) into 3 parts, so I don't need special operators and just use ordinary velocity gradient operators there.

@glwagner
Copy link
Member Author

@glwagner
Copy link
Member Author

@ali-ramadhan feel free to hit the red button when you are ready.

@ali-ramadhan ali-ramadhan merged commit 3f5739e into master Aug 29, 2019
@ali-ramadhan ali-ramadhan deleted the new-closures-no-wall branch August 29, 2019 15:12
arcavaliere pushed a commit to arcavaliere/Oceananigans.jl that referenced this pull request Nov 6, 2019
New and more correct turbulence closures

Former-commit-id: 3f5739e
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug 🐞 Even a perfect program still has bugs feature 🌟 Something new and shiny science 🌊 Sometimes addictive, sometimes allergenic
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants