Add internal tide example in Docs #3132

navidcy · 2023-06-03T23:54:50Z

This is a simple example with immersed boundary.

barotropic_tide.mp4

navidcy · 2023-06-03T23:56:45Z

This PR requires #3094. At the moment it includes the changes there but ideally we want to conclude the discussion in #3094 and merge that one before this.

Or perhaps use the split-explicit free surface?

navidcy · 2023-06-03T23:57:56Z

cc @djlikesdjs

tomchor · 2023-06-04T00:00:29Z

@navidcy can you please include the link to the docs preview page so that we see how it looks before approving?

navidcy · 2023-06-04T00:25:29Z

Could you build the docs locally?

Comment out all examples here

Oceananigans.jl/docs/make.jl

Lines 30 to 43 in 72a2399

    
           examples = [ 
        
               "One-dimensional diffusion"        => "one_dimensional_diffusion", 
        
               "Two-dimensional turbulence"       => "two_dimensional_turbulence", 
        
               "Internal wave"                    => "internal_wave", 
        
               "Barotropic tide"                  => "barotropic_tide", 
        
               "Convecting plankton"              => "convecting_plankton", 
        
               "Ocean wind mixing and convection" => "ocean_wind_mixing_and_convection", 
        
               "Langmuir turbulence"              => "langmuir_turbulence", 
        
               "Baroclinic adjustment"            => "baroclinic_adjustment", 
        
               "Kelvin-Helmholtz instability"     => "kelvin_helmholtz_instability", 
        
               "Shallow water Bickley jet"        => "shallow_water_Bickley_jet", 
        
               "Horizontal convection"            => "horizontal_convection", 
        
               "Tilted bottom boundary layer"     => "tilted_bottom_boundary_layer" 
        
           ]

except the tide one and run:

$ rm -rf docs/build; rm -rf docs/src/generated; julia --project -e 'using Pkg; Pkg.instantiate()'; julia --project=docs/ -e 'using Pkg; Pkg.precompile()'; JULIA_DEBUG=Documenter julia --project=docs/ docs/make.jl

and then open, e.g., docs/build/generated/barotropic_tide.html.

If you even change

Oceananigans.jl/docs/make.jl

Lines 141 to 151 in 72a2399

    
           makedocs(bib, 
        
             sitename = "Oceananigans.jl", 
        
              authors = "Climate Modeling Alliance and contributors", 
        
               format = format, 
        
                pages = pages, 
        
              modules = [Oceananigans], 
        
              doctest = true, 
        
               strict = true, 
        
                clean = true, 
        
            checkdocs = :exports 
        
           )

to doctest = false, strict = false and checkdocs = :none, it'll speed things even more.

(It took ~8 mins on my laptop.)

tomchor · 2023-06-04T18:30:08Z

Could you build the docs locally?

I know how to do it, but I didn't try. I thought the best practice was to set push_preview=true temporarily for docs changes like this one as per your first comment in #2870.

But fair enough. I'll build and view them locally for review.

simone-silvestri · 2023-06-04T19:14:50Z

examples/barotropic_tide.jl

+hill(x) = h₀ * exp(-x^2 / 2width^2)
+bottom(x, y) = - H + hill(x)
+
+grid = ImmersedBoundaryGrid(underlying_grid, GridFittedBottom(bottom))


This example is probably the perfect test for PartialCellBottom

Yes, or even better for shaved cells!
I'd be very interested to see if these streaks of w-velocities in the initial condition diminish with shaved cells.

cc @siddharthabishnu

Why don't we use PartialCellBottom here?

When I tried to change

grid = ImmersedBoundaryGrid(underlying_grid, GridFittedBottom(bottom))

to

using Oceananigans.ImmersedBoundaries: PartialCellBottom grid = ImmersedBoundaryGrid(underlying_grid, PartialCellBottom(bottom))

the simulation NaNed immediately. Even when I reduced the max_Δt from 6 minutes down to 2 minutes..

navidcy · 2024-03-07T15:38:25Z

Can't we just plot a fixed-color shape on top on the mount on all panels? Makie.poly seems to be able to do that.

We can.

We can also mask the output after we load it, e.g, via

using something like

using Oceananigans.ImmersedBoundaries: mask_immersed_field!

function mask_and_get_interior(φ_t, n; value=NaN)
    mask_immersed_field!(φ_t[n], value)
    return interior(φ_t[n], :, 1, :)
end

u′ₙ = @lift mask_and_get_interior(u′_t, $n)

which gives

internal_tide.mp4

But either of these solutions complicate the example a bit. Ideally, mask_immersed_field! should not be user-facing. And plotting a mountain on top of the mountain is a bit of a hack and would require a bit of explaining and justifying why we do that.

navidcy · 2024-03-07T15:45:45Z

But I wouldn't mind adding the post-processing masking step with the hope we can eliminate it when masking is enabled for OutputWriters; see #3092

simone-silvestri · 2024-03-07T15:49:23Z

I am in favor of not doing 979c2a7. It hardcodes a change in the source code that is only required for aesthetic reasons. I would wait for #3092 and in the meantime mask the output in the example

navidcy · 2024-03-07T16:06:24Z

I am in favor of not doing 979c2a7. It hardcodes a change in the source code that is only required for aesthetic reasons. I would wait for #3092 and in the meantime mask the output in the example

Deal, done

francispoulin · 2024-03-08T12:52:04Z

The movie looks great and I'm happy to see this example coming to life!

One question I have is in the plot for stratification, at the very top. It looks like there is a black line and then there is a thin layer right below that. Anyone know what that is?

glwagner · 2024-03-08T20:28:58Z

I am in favor of not doing 979c2a7. It hardcodes a change in the source code that is only required for aesthetic reasons. I would wait for #3092 and in the meantime mask the output in the example

I agree with @simone-silvestri that we want this feature, but that it should be done properly

glwagner · 2024-03-08T20:35:46Z

Can't we just plot a fixed-color shape on top on the mount on all panels? Makie.poly seems to be able to do that.

We can.

We can also mask the output after we load it, e.g, via

using something like
using Oceananigans.ImmersedBoundaries: mask_immersed_field!

function mask_and_get_interior(φ_t, n; value=NaN)
    mask_immersed_field!(φ_t[n], value)
    return interior(φ_t[n], :, 1, :)
end

u′ₙ = @lift mask_and_get_interior(u′_t, $n)
which gives

internal_tide.mp4
But either of these solutions complicate the example a bit. Ideally, mask_immersed_field! should not be user-facing. And plotting a mountain on top of the mountain is a bit of a hack and would require a bit of explaining and justifying why we do that.

I think we should mask the solution during visualization in the script. This is actually what users will need to do currently to make decent visualizations, so it is good to illustrate how to do it --- even if, hopefully, we will eventually have a better solution.

As for plotting a shape on top, I'm inclined to encourage visualization that directly represents the domain / data for the purposes of an Oceananigans example. I think its ok if people want to incorporate visualizations like that in their own work but we may not want to promote it as "the best" way to visualize complex domains in this example.

navidcy · 2024-03-09T08:04:16Z

I also agree. That’s how the example is now. If you approve then please approve?

navidcy · 2024-03-09T12:24:36Z

Masking is done:

Oceananigans.jl/examples/internal_tide.jl

Lines 213 to 217 in ca43714

    
           using Oceananigans.ImmersedBoundaries: mask_immersed_field! 
        
           for φ_t in (u′_t, w_t, N²_t), n in 1:length(times) 
        
               mask_immersed_field!(φ_t[n], NaN) 
        
           end

navidcy · 2024-03-09T12:45:52Z

One question I have is in the plot for stratification, at the very top. It looks like there is a black line and then there is a thin layer right below that. Anyone know what that is?

@francispoulin nice catch!
You are correct... it's because b has no flux boundary conditions and therefore the first halo point outside the interior (top and bottom) have values exactly equal to the values just inside (so that the derivative is zero). And when we compute N² = ∂z(b) the top and bottom values are zero.

using Oceananigans

Nx, Nz = 25, 12

underlying_grid = RectilinearGrid(size = (Nx, Nz), x = (-1000e3, 1000e3), z = (-2000, 0),
                                  halo = (4, 4), topology = (Periodic, Flat, Bounded))

h₀, width = 250, 20e3
bottom(x) = - H + h₀ * exp(-x^2 / 2width^2)

grid = ImmersedBoundaryGrid(underlying_grid, GridFittedBottom(bottom))

model = HydrostaticFreeSurfaceModel(; grid, buoyancy = BuoyancyTracer(), tracers = :b)

Nᵢ² = 1e-4
bᵢ(x, z) = Nᵢ² * z
set!(model, b=bᵢ)

b = model.tracers.b
N² = Field(∂z(b))
compute!(N²)

using GLMakie

xb, _, zb = nodes(b)
fig = Figure()
ax = Axis(fig[1, 1])
hm = heatmap!(ax, xb, zb, interior(b, :, 1, :))
Colorbar(fig[1, 2], hm)
save("b.png", current_figure())

xN², _, zN² = nodes(N²)

fig = Figure()
ax = Axis(fig[1, 1])
hm = heatmap!(ax, xN², zN², interior(N², :, 1, :), colorrange = (0, 2Nᵢ²))
Colorbar(fig[1, 2], hm)
save("N2.png", current_figure())

I don't know what would be the solution to this. I think it's physical in our domain to have no flux condition of b and therefore ∂z(b) = 0 at top and bottom, right?

francispoulin · 2024-03-09T12:53:53Z

Makes sense.

In this case would it be better to specify a non zero flux on the buoyancy, the ambient stratification? It won't change the nice results you are getting over the topography but makes sense to me. No need to make this change.

navidcy · 2024-03-09T13:01:14Z

I think the proper way to do this is to use BackgroundField for the ambient stratification. I’ll give it a go

navidcy · 2024-03-09T13:10:57Z

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

glwagner · 2024-03-09T17:07:57Z

Masking is done:

Oceananigans.jl/examples/internal_tide.jl

Lines 213 to 217 in ca43714

using Oceananigans.ImmersedBoundaries: mask_immersed_field!

for φ_t in (u′_t, w_t, N²_t), n in 1:length(times)

mask_immersed_field!(φ_t[n], NaN)

end

That's great.

glwagner · 2024-03-09T17:10:50Z

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

Right, we do not support that right now. I think there are some niche uses but, for the most part, when we go hydrostatic we are often doing large scale and a bit more complicated modeling with parameterizations, etc, where getting background fields right would be more of a pain than with NonhydrostaticModel (which still does not work properly with nonlinear subgrid-scale turbulence closures).

I think its better to go more mainstream for this example. I think the use cases that this example serves as a starting point for are probably not going to use background fields, typically.

navidcy · 2024-03-09T17:14:25Z

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

Right, we do not support that right now. I think there are some niche uses but, for the most part, when we go hydrostatic we are often doing large scale and a bit more complicated modeling with parameterizations, etc, where getting background fields right would be more of a pain than with NonhydrostaticModel (which still does not work properly with nonlinear subgrid-scale turbulence closures).

I think its better to go more mainstream for this example. I think the use cases that this example serves as a starting point for are probably not going to use background fields, typically.

True true.

Any idea how to deal with the dz(b) that vanishes at the top? Or should we live with it?

glwagner · 2024-03-09T17:56:24Z

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

Hm.... HydrostaticFreeSurfaceModel does not support BackgroundFields....

Right, we do not support that right now. I think there are some niche uses but, for the most part, when we go hydrostatic we are often doing large scale and a bit more complicated modeling with parameterizations, etc, where getting background fields right would be more of a pain than with NonhydrostaticModel (which still does not work properly with nonlinear subgrid-scale turbulence closures).
I think its better to go more mainstream for this example. I think the use cases that this example serves as a starting point for are probably not going to use background fields, typically.

True true.

Any idea how to deal with the dz(b) that vanishes at the top? Or should we live with it?

Because of the no-flux boundary condition? I think that's correct, ok to have. You could have a mixed layer at the top too!

navidcy · 2024-03-10T09:21:52Z

hm... I don't know about a mixed layer...

let's just ignore the dz(b)=0 at the top&bottom or we can add a remark in the example about it?

glwagner · 2024-03-10T17:33:42Z

hm... I don't know about a mixed layer...

let's just ignore the dz(b)=0 at the top&bottom or we can add a remark in the example about it?

Yes definitely ignore.

What we want to focus on is providing an example that is a useful starting point for other, new simulations. This is the key priority.

This is why it makes sense to prioritize relatively simple examples (and also simple visualization) that can be quickly digested and are amenable to generalization.

We want to deprioritize bespoke or niche cases that are a dead-end. For example, if we use GradientBoundaryCondition at the top in this case, then most users will have to delete / change that line in order to generalize this example into the case they are interested in using.

glwagner

This is a great start and I think in the future we can think about adding bottom drag using ImmersedBoundaryCondition

navidcy added 4 commits June 4, 2023 09:43

don't hardcode 3 halo points for pcg

a32f127

add barotropic tide example

e6a12ba

add barotropic tide example

d122aaf

run for 4 days

292b4e8

navidcy added the documentation 📜 The sacred scrolls label Jun 3, 2023

navidcy requested review from siddharthabishnu, tomchor, glwagner and simone-silvestri June 3, 2023 23:54

tidal forcing 0 at t=0

72a2399

navidcy added 14 commits June 4, 2023 10:29

fix latex rendering

29a44fb

adds remarks about saved output

d2258d8

add more explanation regarding the tidal forcing

bbd8ee4

add remarks about IBG

900db08

more remarks about IBG

e067cea

(i) -> (ii)

8cbe588

first load, then plot

ba801a2

missing _2

b3f39a4

some more remarks

803d21b

fix bottom redefinition + beautifications

dc0ea63

simpler panel titles

fbde03b

simpler plotting

373e5f5

z [m] -> km

9f6590d

simplify

b8c0efb

simone-silvestri reviewed Jun 4, 2023

View reviewed changes

navidcy added 4 commits March 7, 2024 18:00

revert change

06e20a5

add post-process masking step

da4502a

revert changes

84bebe8

find max abs values before masking

cf18546

navidcy added 2 commits March 8, 2024 10:37

Merge branch 'main' into ncc/tide

f8c0dea

plot domain with band

cf1ad7b

simpler masking

ca43714

fix missing .jl

6777dd7

glwagner approved these changes Mar 11, 2024

View reviewed changes

navidcy merged commit f1aacf3 into main Mar 11, 2024
48 checks passed

navidcy deleted the ncc/tide branch March 11, 2024 15:49

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Add internal tide example in Docs #3132

Add internal tide example in Docs #3132

navidcy commented Jun 3, 2023 •

edited

Loading

navidcy commented Jun 3, 2023

navidcy commented Jun 3, 2023

tomchor commented Jun 4, 2023

navidcy commented Jun 4, 2023

tomchor commented Jun 4, 2023 •

edited

Loading

simone-silvestri Jun 4, 2023

navidcy Jun 4, 2023

navidcy Jun 4, 2023

glwagner Jun 4, 2023

navidcy Jun 5, 2023

navidcy commented Mar 7, 2024 •

edited

Loading

navidcy commented Mar 7, 2024

simone-silvestri commented Mar 7, 2024

navidcy commented Mar 7, 2024

francispoulin commented Mar 8, 2024

glwagner commented Mar 8, 2024

glwagner commented Mar 8, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024 •

edited

Loading

francispoulin commented Mar 9, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024

glwagner commented Mar 9, 2024

glwagner commented Mar 9, 2024 •

edited by navidcy

Loading

navidcy commented Mar 9, 2024

glwagner commented Mar 9, 2024

navidcy commented Mar 10, 2024

glwagner commented Mar 10, 2024

glwagner left a comment

Add internal tide example in Docs #3132

Add internal tide example in Docs #3132

Conversation

navidcy commented Jun 3, 2023 • edited Loading

navidcy commented Jun 3, 2023

navidcy commented Jun 3, 2023

tomchor commented Jun 4, 2023

navidcy commented Jun 4, 2023

tomchor commented Jun 4, 2023 • edited Loading

simone-silvestri Jun 4, 2023

Choose a reason for hiding this comment

navidcy Jun 4, 2023

Choose a reason for hiding this comment

navidcy Jun 4, 2023

Choose a reason for hiding this comment

glwagner Jun 4, 2023

Choose a reason for hiding this comment

navidcy Jun 5, 2023

Choose a reason for hiding this comment

navidcy commented Mar 7, 2024 • edited Loading

navidcy commented Mar 7, 2024

simone-silvestri commented Mar 7, 2024

navidcy commented Mar 7, 2024

francispoulin commented Mar 8, 2024

glwagner commented Mar 8, 2024

glwagner commented Mar 8, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024 • edited Loading

francispoulin commented Mar 9, 2024

navidcy commented Mar 9, 2024

navidcy commented Mar 9, 2024

glwagner commented Mar 9, 2024

glwagner commented Mar 9, 2024 • edited by navidcy Loading

navidcy commented Mar 9, 2024

glwagner commented Mar 9, 2024

navidcy commented Mar 10, 2024

glwagner commented Mar 10, 2024

glwagner left a comment

Choose a reason for hiding this comment

navidcy commented Jun 3, 2023 •

edited

Loading

tomchor commented Jun 4, 2023 •

edited

Loading

navidcy commented Mar 7, 2024 •

edited

Loading

navidcy commented Mar 9, 2024 •

edited

Loading

glwagner commented Mar 9, 2024 •

edited by navidcy

Loading