FINALLY reverse the k index #462

ali-ramadhan · 2019-10-11T20:16:55Z

This PR finally reverses the k index although it's still a work in progress as a few tests fail. In particular, it seems that incompressibility is not satisfied even though I managed to get the "recomputing w from continuity" test to pass.

Need to think about this some more.

Resolves #90
Resolves #468
Resolves #480

glwagner · 2019-10-11T20:26:38Z

Are we integrating bottom-up or top-down now?

Ultimately it might be nice to eliminate that step and explicitly time-step the vertical velocity, since thats valid for general boundary conditions...

src/time_steppers.jl

ali-ramadhan · 2019-10-15T13:46:03Z

For w is

k = Nz now the surface implying that k=0 is the bottom and a halo cell?
or is k = 1 the bottom implying that k = Nz+1 is the surface and a halo cell?

Seems that we're actually implementing the second.

Might be nice to have direct access to the surface via k = Nz but this would break the pattern that i = 1 for u and j = 1 for v coincide with the boundary whereas i = Nx and j = Ny do not.

glwagner · 2019-10-15T13:52:44Z

The convention we use for the other variables stored on faces is that our "interior" data spans from i=1 to i=N. Thus we call the i=N+1 a "halo point" for fields on faces, even though it really is not (its on the boundary, see #455).

In addition, we use a convention that face indices are "left" of cell indices. Thus iFace=1 corresponds to the "left" boundary of the domain. iCell=N is the rightmost cell in the domain; and iFace=N+1 corresponds to the rightmost boundary.

Why do you want "direct" access to the surface at k=Nz? Can you explain why you want this? This is confusing to me and I don't see why we would want to have an odd and surprising convention for the vertical coordinate.

The fact is simply that the face at i=N is a very concrete and real place. It lies in the interior of the domain, to the left of the cell point i=N. This is a fact of our grid and must be internalized by anyone who wants to work with the discrete data and index fields directly by raw index.

If/when we use named axis conventions for fields (see #457), we would be able to access the surface via u[Z(At(0))] (or some syntax like that. We can also create an alias Surface() = Z(At(0)).

ali-ramadhan · 2019-10-15T13:59:16Z

Why do you want "direct" access to the surface at k=Nz? Can you explain why you want this?

We've mostly been concerned with surface ocean physics so naively it would seem nice to have the surface not be a halo cell.

Might have been an issue in the past when we didn't have vertical halos because then you don't even have a cell for the surface.

Out of curiousity then, would an implementation of surface wave effects (#443) have to touch the halos?

glwagner · 2019-10-15T14:03:43Z

We've mostly been concerned with surface ocean physics so naively it would seem nice to have the surface not be a halo cell.

Your suggestion seems more to change indexing behavior rather than convention for halos --- but maybe thats not true given that we use the halos more as ghost points now, rather than overlap regions for MPI stuff, etc.

Presumably, we could simply define the i=N+1 point for fields on faces as "within" the domain, and define the points outside of i=1:N+1 as halo points. This is a convention of our halo regions?

codecov · 2019-10-16T00:58:56Z

Codecov Report

Merging #462 into master will not change coverage.
The diff coverage is 86.45%.

@@           Coverage Diff           @@
##           master     #462   +/-   ##
=======================================
  Coverage   73.34%   73.34%           
=======================================
  Files          27       27           
  Lines        1508     1508           
=======================================
  Hits         1106     1106           
  Misses        402      402

Impacted Files	Coverage Δ
src/Oceananigans.jl	`62.5% <ø> (ø)`	⬆️
src/output_writers.jl	`80.33% <ø> (ø)`	⬆️
src/utils.jl	`81.55% <0%> (ø)`	⬆️
src/grids.jl	`93.54% <100%> (ø)`	⬆️
src/TurbulenceClosures/closure_operators.jl	`64.44% <100%> (ø)`	⬆️
src/diagnostics.jl	`79.31% <50%> (ø)`	⬆️
src/boundary_conditions.jl	`80.59% <72.72%> (ø)`	⬆️
src/halo_regions.jl	`90.27% <87.23%> (+0.13%)`	⬆️
src/Operators/ops_regular_cartesian_grid.jl	`84.21% <92.3%> (ø)`	⬆️
src/time_steppers.jl	`76.08% <92.85%> (-0.13%)`	⬇️

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 5114a22...419f176. Read the comment docs.

ali-ramadhan · 2019-10-16T01:00:09Z

Fixed a few more things and I've updated recompute_w_from_continuity! so that the integral is computed starting from the bottom upwards.

Regression tests fail so I'll see what's wrong there tomorrow.

We should probably go in this order:

Merge LES regression tests (PR Regression tests for models with LES closures #479)
Merge arbitrary tracers (PR Arbitrary tracers #452)
Make sure this PR passes the LES regression tests and works with arbitrary tracers.

ali-ramadhan · 2019-10-18T01:21:47Z

Regression tests all pass on my laptop!

Next I'll merge in master to make sure reverse-k-index works with arbitrary tracers and can pass the LES regression tests.

Then we can regenerate the regression data with a reversed k index.

Differences are all very small so I'm confident we're producing the same numbers as before (up to machine precision and truncation error):

Running regression tests...
  Testing thermal bubble regression [CPU]
[ Info: Δu: min=-3.01224e-15, max=3.01289e-15, mean=-3.94216e-22, absmean=5.13638e-16, std=7.36011e-16
[ Info: Δv: min=-3.01208e-15, max=3.00899e-15, mean=-1.19521e-21, absmean=5.13326e-16, std=7.35825e-16
[ Info: Δw: min=-7.64756e-15, max=5.88619e-15, mean=1.35968e-21, absmean=1.03765e-15, std=1.48909e-15
[ Info: ΔT: min=-1.77636e-15, max=1.77636e-15, mean=7.80626e-18, absmean=3.55618e-17, std=2.51247e-16
[ Info: ΔS: min=0, max=0, mean=0, absmean=0, std=0
  Testing Rayleigh–Bénard tracer regression [CPU]
[ Info: Δu: min=-3.65818e-14, max=3.21965e-14, mean=-1.92406e-18, absmean=3.3529e-15, std=4.83978e-15
[ Info: Δv: min=-2.85535e-14, max=3.72202e-14, mean=1.62974e-19, absmean=3.29709e-15, std=4.79574e-15
[ Info: Δw: min=-4.56857e-14, max=4.32987e-14, mean=2.54957e-18, absmean=2.97761e-15, std=4.68161e-15
[ Info: ΔT: min=-5.19029e-14, max=3.71925e-14, mean=8.223e-18, absmean=2.55428e-15, std=4.14735e-15
[ Info: ΔS: min=-2.37588e-14, max=1.60982e-14, mean=1.15112e-18, absmean=1.13258e-15, std=1.79417e-15

But the tracer values at the top and bottom boundaries are wrong after one time step, so must be a boundary issue...

ali-ramadhan · 2019-10-21T11:29:15Z

After reversing the k index u, v, w, and S are all correct if you take one time step for the LES regression tests (non-LES regression tests all pass). But T values at the top and bottom grid cells are wrong.

Still thinking of what could be effecting T but not S (which also runs with AMD and flux boundary conditions)...

T is the only field with a gradient boundary condition, but enforcing gradient boundary conditions shouldn't have changed and looks fine to me...

I think this is the last bug, if fixed this should hopefully become ready to merge .

src/halo_regions.jl

ali-ramadhan · 2019-10-21T22:10:03Z

Thanks @glwagner! All four regression tests pass now!

I'll regenerate the regression data in a subsequent PR which also removes all the reverse statements.

ali-ramadhan · 2019-10-21T22:20:57Z

Differences are close to machine epsilon so pretty sure this branch produces the same numbers as before, up to truncation error (and we do integrate from the bottom up to enforce continuity now).

Running regression tests...
  Testing thermal bubble regression [CPU]
[ Info: Δu: min=-3.01224e-15, max=3.01289e-15, mean=-3.94216e-22, absmean=5.13638e-16, std=7.36011e-16
[ Info: Δv: min=-3.01208e-15, max=3.00899e-15, mean=-1.19521e-21, absmean=5.13326e-16, std=7.35825e-16
[ Info: Δw: min=-7.64756e-15, max=5.88619e-15, mean=1.35968e-21, absmean=1.03765e-15, std=1.48909e-15
[ Info: ΔT: min=-1.77636e-15, max=1.77636e-15, mean=7.80626e-18, absmean=3.55618e-17, std=2.51247e-16
[ Info: ΔS: min=0, max=0, mean=0, absmean=0, std=0
  Testing Rayleigh–Bénard tracer regression [CPU]
[ Info: Δu: min=-3.65818e-14, max=3.21965e-14, mean=-1.92406e-18, absmean=3.3529e-15, std=4.83978e-15
[ Info: Δv: min=-2.85535e-14, max=3.72202e-14, mean=1.62974e-19, absmean=3.29709e-15, std=4.79574e-15
[ Info: Δw: min=-4.56857e-14, max=4.32987e-14, mean=2.54957e-18, absmean=2.97761e-15, std=4.68161e-15
[ Info: ΔT: min=-5.19029e-14, max=3.71925e-14, mean=8.223e-18, absmean=2.55428e-15, std=4.14735e-15
[ Info: ΔS: min=-2.37588e-14, max=1.60982e-14, mean=1.15112e-18, absmean=1.13258e-15, std=1.79417e-15
  Testing oceanic large eddy simulation regression [VerstappenAnisotropicMinimumDissipation, CPU]
[ Info: Δu: min=-8.1532e-16, max=7.47666e-16, mean=3.10996e-19, absmean=1.34093e-16, std=1.74514e-16
[ Info: Δv: min=-7.94503e-16, max=7.32053e-16, mean=1.62668e-19, absmean=1.30967e-16, std=1.70133e-16
[ Info: Δw: min=-2.21871e-15, max=2.4928e-15, mean=-4.54454e-18, absmean=2.67538e-16, std=4.09625e-16
[ Info: ΔT: min=-1.06581e-14, max=1.06581e-14, mean=-1.76661e-17, absmean=1.34445e-15, std=2.4561e-15
[ Info: ΔS: min=-2.13163e-14, max=2.13163e-14, mean=-9.38131e-17, absmean=2.61336e-15, std=4.82786e-15
  Testing oceanic large eddy simulation regression [SmagorinskyLilly, CPU]
[ Info: Δu: min=-7.65013e-16, max=8.1532e-16, mean=-4.16382e-20, absmean=1.23114e-16, std=1.62955e-16
[ Info: Δv: min=-7.56339e-16, max=8.43076e-16, mean=-1.08495e-19, absmean=1.21821e-16, std=1.61572e-16
[ Info: Δw: min=-2.19616e-15, max=1.9585e-15, mean=6.95974e-18, absmean=2.40883e-16, std=3.76408e-16
[ Info: ΔT: min=-1.06581e-14, max=1.06581e-14, mean=-9.13764e-18, absmean=1.30546e-15, std=2.43045e-15
[ Info: ΔS: min=-2.13163e-14, max=2.13163e-14, mean=-1.34019e-17, absmean=2.589e-15, std=4.83055e-15

glwagner

Looks great. Also happy about the house cleaning. Looks like some NetCDF stuff snuck in here too; not sure how that happened, but it seems minor.

glwagner · 2019-10-22T14:07:27Z

src/Oceananigans.jl

@@ -49,7 +49,7 @@ export

    # Model output writers
    Checkpointer, restore_from_checkpoint, read_output,
-    JLD2OutputWriter, FieldOutput, FieldOutputs,
+    JLD2OutputWriter, NetCDFOutputWriter, FieldOutput, FieldOutputs,


Yeah this was an issue after I merged master in, I guess we forgot to export it in PR #433.

Bad practice on my part to include changes not relevant to reversing the k index.

glwagner · 2019-10-22T14:10:30Z

src/boundary_conditions.jl

@@ -214,7 +211,7 @@ function ChannelBCs(;  north = BoundaryCondition(Flux, nothing),

    x = PeriodicBCs()
    y = CoordinateBoundaryConditions(south, north)
-    z = CoordinateBoundaryConditions(top, bottom)
+    z = CoordinateBoundaryConditions(bottom, top)


Very satisfying.

glwagner · 2019-10-22T14:15:43Z

test/test_turbulence_closures.jl

@@ -136,9 +135,8 @@ function test_function_differentiation(T=Float64)
    ∂y_ϕ_f = ϕ²[2, 2, 2] - ϕ²[2, 1, 2]
    ∂y_ϕ_c = ϕ²[2, 3, 2] - ϕ²[2, 2, 2]

-    # Note reverse indexing here!


ali-ramadhan · 2019-10-22T14:31:04Z

I was thinking of releasing v0.13 once this is merged (need to sort out v0.12 tag first).

Would be good to tag a version and be able to say "pre-v0.13" is bad k-indexing and "v0.13+" is good k-indexing as this is a pretty big change and model output is now different (some analysis and plotting scripts will have to be modified).

FINALLY reverse the k index Former-commit-id: 41a2b55

ali-ramadhan added 5 commits September 20, 2019 07:34

Reverse index for operators

ed46bb2

Reverse index for turbulence operators

0a929bd

Reverse index for halos

826d302

Reverse index for compute_w_from_continuity

9ec27d3

Reverse zC and zF

534c4c1

ali-ramadhan added cleanup 🧹 Paying off technical debt numerics 🧮 So things don't blow up and boil the lobsters alive labels Oct 11, 2019

glwagner reviewed Oct 11, 2019

View reviewed changes

src/time_steppers.jl Outdated Show resolved Hide resolved

glwagner mentioned this pull request Oct 12, 2019

AbstractOperations for diagnostics (and beyond...) #463

Merged

4 tasks

ali-ramadhan mentioned this pull request Oct 12, 2019

Strategy for implementing and testing models with vertically stretched grids #471

Closed

ali-ramadhan added 2 commits October 15, 2019 08:28

Fix grid tests

d362be6

Model is incompressible again

6f47019

ali-ramadhan added 3 commits October 15, 2019 17:29

Update k indexing on turbulence closure tests

8e26cf8

Add test for grid range lengths

2343bdf

Update variable names in ancient grid tests

79d055d

ali-ramadhan mentioned this pull request Oct 15, 2019

Round-off error in RegularCartesianGrid produces ranges of wrong length #480

Closed

ali-ramadhan added 7 commits October 15, 2019 18:09

Add test for roundoff error in ranges that fails

43bd7e9

Fix grid ranges.

b08f6fa

Rename Pearson vortex to Taylor-Green vortex

52d85ea

Update halo filling for no-penetration boundary conditions.

05fb3ac

Integrate from bottom upwards to recompute w from continuity.

1377c34

Fix indexing on turbulence closure tests.

d1c906f

Perform buoyancy integral from the surface downwards.

7c1fa1c

ali-ramadhan added 6 commits October 17, 2019 19:12

Correct buoyancy integral

9a8aa2f

Reverse output for rising thermal bubble regression test.

bc049e1

Remapping :bottom -> :left and :top -> :right.

20be6f4

Reverse k index and remap top/bottom for filling VBC and GBC halos

0445cd4

Reverse output for Rayleigh-Benard convection regression test.

5ea35a3

Only compare k = 2:Nz as w[:, :, 1] = 0.

0de2765

ali-ramadhan added 7 commits October 17, 2019 22:22

Fix bottom/top alias when using CoordinateBoundaryConditions.

fb8d906

Fix flux boundary conditions after reversing k index.

8f09f58

Merge branch 'master' into reverse-k-index

f0a91d6

Update old regression tests and fix post-merge typo.

a8b5a9a

Reverse output for LES regression tests.

5dff125

But the tracer values at the top and bottom boundaries are wrong after one time step, so must be a boundary issue...

Fix bottom/top mixup in applying flux boundary conditions.

6ee988b

return -> return nothing

ed452ff

glwagner reviewed Oct 21, 2019

View reviewed changes

src/halo_regions.jl Outdated Show resolved Hide resolved

Switch _fill_bottom_halo! and _fill_top_halo! where appropriate.

861a0ea

Merge branch 'master' into reverse-k-index

4170bbd

ali-ramadhan requested a review from glwagner October 21, 2019 22:16

ali-ramadhan added 2 commits October 21, 2019 22:38

Need to export NetCDFOutputWriter.

41575ad

Fix thermal bubble regression test by directly using NCDatasets.jl

419f176

glwagner approved these changes Oct 22, 2019

View reviewed changes

ali-ramadhan changed the title ~~[WIP] FINALLY reverse the k index~~ FINALLY reverse the k index Oct 22, 2019

ali-ramadhan merged commit 41a2b55 into master Oct 22, 2019

ali-ramadhan deleted the reverse-k-index branch October 22, 2019 21:11

ali-ramadhan mentioned this pull request Nov 4, 2019

Regenerate regression test data following reversal of k index #526

Merged

arcavaliere pushed a commit to arcavaliere/Oceananigans.jl that referenced this pull request Nov 6, 2019

FINALLY reverse the k index (CliMA#462)

4c468ce

FINALLY reverse the k index Former-commit-id: 41a2b55

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

FINALLY reverse the k index #462

FINALLY reverse the k index #462

ali-ramadhan commented Oct 11, 2019 •

edited

Loading

glwagner commented Oct 11, 2019

ali-ramadhan commented Oct 15, 2019 •

edited

Loading

glwagner commented Oct 15, 2019 •

edited

Loading

ali-ramadhan commented Oct 15, 2019

glwagner commented Oct 15, 2019 •

edited

Loading

codecov bot commented Oct 16, 2019 •

edited

Loading

ali-ramadhan commented Oct 16, 2019 •

edited

Loading

ali-ramadhan commented Oct 18, 2019 •

edited

Loading

ali-ramadhan commented Oct 21, 2019 •

edited

Loading

ali-ramadhan commented Oct 21, 2019 •

edited

Loading

ali-ramadhan commented Oct 21, 2019

glwagner left a comment

glwagner Oct 22, 2019

ali-ramadhan Oct 22, 2019

glwagner Oct 22, 2019

glwagner Oct 22, 2019

ali-ramadhan commented Oct 22, 2019 •

edited

Loading

FINALLY reverse the k index #462

FINALLY reverse the k index #462

Conversation

ali-ramadhan commented Oct 11, 2019 • edited Loading

glwagner commented Oct 11, 2019

ali-ramadhan commented Oct 15, 2019 • edited Loading

glwagner commented Oct 15, 2019 • edited Loading

ali-ramadhan commented Oct 15, 2019

glwagner commented Oct 15, 2019 • edited Loading

codecov bot commented Oct 16, 2019 • edited Loading

Codecov Report

ali-ramadhan commented Oct 16, 2019 • edited Loading

ali-ramadhan commented Oct 18, 2019 • edited Loading

ali-ramadhan commented Oct 21, 2019 • edited Loading

ali-ramadhan commented Oct 21, 2019 • edited Loading

ali-ramadhan commented Oct 21, 2019

glwagner left a comment

Choose a reason for hiding this comment

glwagner Oct 22, 2019

Choose a reason for hiding this comment

ali-ramadhan Oct 22, 2019

Choose a reason for hiding this comment

glwagner Oct 22, 2019

Choose a reason for hiding this comment

glwagner Oct 22, 2019

Choose a reason for hiding this comment

ali-ramadhan commented Oct 22, 2019 • edited Loading

ali-ramadhan commented Oct 11, 2019 •

edited

Loading

ali-ramadhan commented Oct 15, 2019 •

edited

Loading

glwagner commented Oct 15, 2019 •

edited

Loading

glwagner commented Oct 15, 2019 •

edited

Loading

codecov bot commented Oct 16, 2019 •

edited

Loading

ali-ramadhan commented Oct 16, 2019 •

edited

Loading

ali-ramadhan commented Oct 18, 2019 •

edited

Loading

ali-ramadhan commented Oct 21, 2019 •

edited

Loading

ali-ramadhan commented Oct 21, 2019 •

edited

Loading

ali-ramadhan commented Oct 22, 2019 •

edited

Loading