Generalize `TwoBandPhotosyntheticallyActiveRadiation` kernel to allow curvilinear grids #120

jagoosw · 2023-07-04T12:47:54Z

Not sure this is the most efficient way todo this, but since xnode is defiend differedntly for different grid types this is by far the easiest fix. The alternative would be to define xnode(i, grid, lx, ly, lz) for RectilinearGrids. Perhaps that solution would be better.

Closes #119

Not sure this is the most efficient way todo this, but since `xnode` is defiend differedntly for different grid types this is by far the easiest fix. The alternative would be to define `xnode(i, grid, lx, ly, lz)` for RectilinearGrids. Perhaps that solution would be better.

jagoosw · 2023-07-04T12:48:13Z

Not sure this is the most efficient way todo this, but since xnode is defiend differedntly for different grid
types this is by far the easiest fix. The alternative would be to define xnode(i, grid, lx, ly, lz) for RectilinearGrids.
Perhaps that solution would be better.

jagoosw · 2023-07-04T14:32:17Z

Probably should add some tests on Lat/lon grids

navidcy · 2023-07-04T21:40:57Z

What do you mean by

but since xnode is defined differently for different grid types this is by far the easiest fix

?

navidcy · 2023-07-04T21:45:10Z

xnode(i, j, k, grid, ℓx, ℓy, ℓz)

should work for all grids. If the above is not true then we should fix it in Oceananigans.

Note that xnode(s) will return the zonal coordinate (e.g. in meters), not the grid's native x-coordinate which can be degrees for some grids.

navidcy · 2023-07-04T21:47:48Z

src/Light/2band.jl

@@ -1,7 +1,7 @@
 @kernel function update_TwoBandPhotosyntheticallyActiveRadiation!(PAR, grid, P, surface_PAR, t, PAR_model) 
    i, j = @index(Global, NTuple)

-    x, y = xnode(i, grid, Center()), ynode(j, grid, Center())
+    x, y = @inbounds xnodes(grid, Center(), Center(), Center())[i], ynodes(grid, Center(), Center(), Center())[j]


not sure if k=1 is what we want here, but how about:

Suggested change

x, y = @inbounds xnodes(grid, Center(), Center(), Center())[i], ynodes(grid, Center(), Center(), Center())[j]

x = xnode(i, j, 1, grid, Center(), Center(), Center())

y = ynode(i, j, 1, grid, Center(), Center(), Center())

Note that for general curvilinear grids, the nodes depend on both i, j. So you probably want to use @navidcy's suggestion.

All current grids (and for the forseeable future) are independent of k so it's ok to take k=1. This makes a "thin spherical shell" approximation for the ocean.

jagoosw · 2023-07-05T20:14:49Z

Thank you both @navidcy @glwagner . I'm not sure exactly what I've got confused by, I think it was trying to do xnode(i, j, grid... rather than with i, j, k. I'll apply these suggestions and see if it works.

I was trying to run one of the ClimaOcean examples with a BGC model which initially wasn't working because the light just had xnode(i, grid..., and then I think I just fixed it wrong.

navidcy · 2023-07-05T20:23:09Z

Note that ClimaOcean might need some updating!

codecov · 2023-07-05T21:11:14Z

Codecov Report

Patch coverage: 100.00% and no project coverage change.

Comparison is base (e82d7f5) 61.46% compared to head (bfb86d9) 61.46%.

Additional details and impacted files

@@           Coverage Diff           @@
##             main     #120   +/-   ##
=======================================
  Coverage   61.46%   61.46%           
=======================================
  Files          26       26           
  Lines         968      968           
=======================================
  Hits          595      595           
  Misses        373      373

Impacted Files	Coverage Δ
src/Boundaries/Sediments/Sediments.jl	`61.53% <ø> (ø)`
src/Boundaries/Sediments/simple_multi_G.jl	`92.10% <ø> (ø)`
src/Light/Light.jl	`0.00% <ø> (ø)`
src/Light/2band.jl	`86.11% <100.00%> (ø)`
src/Models/Individuals/SLatissima.jl	`83.87% <100.00%> (ø)`
src/Particles/tracer_tendencies.jl	`100.00% <100.00%> (ø)`

☔ View full report in Codecov by Sentry.
📢 Do you have feedback about the report comment? Let us know in this issue.

jagoosw · 2023-07-05T21:11:31Z

It does work now

Note that ClimaOcean might need some updating!

I did have to remove some compat entries to install OceanBioME with it but it seems to be running fine with Oceananigans 0.83.0

src/Boundaries/Sediments/simple_multi_G.jl

navidcy · 2023-07-06T06:52:31Z

src/Models/Individuals/SLatissima.jl

shall we delete commented code? it's unclear why it's there..

#nodes, weights, normfactor = @inbounds get_nearest_nodes(x, y, z, grid, (Center(), Center(), Center())) #@unroll for (n, weight) in enumerate(weights)

Oh yeah, not sure why that's still there

jagoosw · 2023-07-06T10:17:30Z

I've just thought, it's probably more convenient for a user if the surface PAR is in the grid native coordinates, for example, if doing a global simulation data is going to be in lat/lon coordinates. Using Oceananigans.Grids.node would be the appropriate thing to use here then right?

navidcy · 2023-07-06T10:36:30Z

node(i, j, k, grid, ℓx, ℓy, ℓz)

will return a 1-, 2-, or 3- tuple with the native coordinates of the grid. (1-, 2-, or 3-tuple depending if any of the coordinates are Flat).

jagoosw · 2023-07-06T10:38:40Z

node(i, j, k, grid, ℓx, ℓy, ℓz)
will return a 1-, 2-, or 3- tuple with the native coordinates of the grid. (1-, 2-, or 3-tuple depending if any of the coordinates are Flat).

Oh yeah I mean to fix this. Although as I've done it now, if it returns anything other than a 3-tuple its going to fail right? Is there a way this can work for Flat dimensions?

navidcy · 2023-07-06T10:44:10Z

hm... I think what you did won't work... because you seem to assume the node will return a 3-tuple but it won't necessarily do that

navidcy · 2023-07-06T10:48:29Z

So I'm not sure what you are trying to do so I can't be very specific.
But how about this, is this helpful?

locs = node(i, j, k, grid, Center(), Center(), Center())

value = surface_PAR(loss...)

This way you don't unfold the tuple but just pass it to the surface_PAR as a tuple and use ... to unpack it then?

jagoosw · 2023-07-06T10:49:17Z

I just checked and it seems to work on flat dimensions always returning a 3-tuple, as long as the locations are all not nothing because that's what the different length results are dispatched on:

https://github.com/CliMA/Oceananigans.jl/blob/92791a962c9096746301cdb888a3050e32c4a58b/src/Grids/latitude_longitude_grid.jl#L624-L634

jagoosw · 2023-07-06T10:50:39Z

So I'm not sure what you are trying to do so I can't be very specific. But how about this, is this helpful?
locs = node(i, j, k, grid, Center(), Center(), Center())

value = surface_PAR(loss...)
This way you don't unfold the tuple but just pass it to the surface_PAR as a tuple and use ... to unpack it then?

I guess this would work although we currently have surface_PAR as a function of x and y only, but if this doesn't work we can change that to include z

navidcy · 2023-07-06T10:50:59Z

I just checked and it seems to work on flat dimensions always returning a 3-tuple, as long as the locations are all not nothing because that's what the different length results are dispatched on:

https://github.com/CliMA/Oceananigans.jl/blob/92791a962c9096746301cdb888a3050e32c4a58b/src/Grids/latitude_longitude_grid.jl#L624-L634

Hm... But this, e.g.,
https://github.com/CliMA/Oceananigans.jl/blob/92791a962c9096746301cdb888a3050e32c4a58b/src/Grids/latitude_longitude_grid.jl#L629
is a 2-tuple.

navidcy · 2023-07-06T10:53:32Z

So I'm not sure what you are trying to do so I can't be very specific. But how about this, is this helpful?
locs = node(i, j, k, grid, Center(), Center(), Center())

value = surface_PAR(loss...)
This way you don't unfold the tuple but just pass it to the surface_PAR as a tuple and use ... to unpack it then?
I guess this would work although we currently have surface_PAR as a function of x and y only, but if this doesn't work we can change that to include z

Sorry, ignore my suggestion then.

Perhaps I just don't know what you are trying to do and I was assuming things.

Let's restart: Yes, to get the coordinates of the grid at i, j, k then node(i,j, k, grid,...) will do the job. :)

jagoosw · 2023-07-06T10:56:30Z

I just checked and it seems to work on flat dimensions always returning a 3-tuple, as long as the locations are all not nothing because that's what the different length results are dispatched on:
https://github.com/CliMA/Oceananigans.jl/blob/92791a962c9096746301cdb888a3050e32c4a58b/src/Grids/latitude_longitude_grid.jl#L624-L634

Hm... But this, e.g., https://github.com/CliMA/Oceananigans.jl/blob/92791a962c9096746301cdb888a3050e32c4a58b/src/Grids/latitude_longitude_grid.jl#L629 is a 2-tuple.

Since we're always giving the locs as Center/Center/Center won't it always use the first method though? I guess it uses the others if you were to e.g. do locations(velocity field) which would return a nothing location on the flat dimension, and then use the function?

navidcy · 2023-07-06T10:58:08Z

Since we're always giving the locs as Center/Center/Center won't it always use the first method though? I guess it uses the others if you were to e.g. do locations(velocity field) which would return a nothing location on the flat dimension, and then use the function?

True! My bad...

navidcy · 2023-07-06T19:57:52Z

should we add a test?

jagoosw · 2023-07-13T11:43:35Z

I've changed the light test to test on different grids and models now, I don't think we need the other tests todo this since the others just test the properties of the bgc model components and their ability to work in different grids etc. should be being tested in the Oceananigans tests right?

I think maybe I should add a test that the particle tendencies and sediment models work across all model types. Although again their ability to function should be caught by the Oceananigans tests now. Perhaps I should test their combinations since that's an issue I found recently.

jagoosw · 2023-07-13T11:47:04Z

I forgot I hadn't pushed the sediment model fixes. I've just made this PR #121 which has them so I think we should merge these tests and fixes, then look at adding more sediment model tests in that PR.

navidcy · 2023-07-15T03:56:15Z

test/test_light.jl

+
+archs = (CPU(), )
+
+@testset "Light attenuaiton model" begin 


Suggested change

@testset "Light attenuaiton model" begin

@testset "Light attenuation model" begin

navidcy · 2023-07-15T03:57:32Z

test/test_light.jl

-    Pᵢ(x,y,z) = 2.5 + z
+    model = model_type(; grid, 
+                         biogeochemistry,
+                         tracers = unique((required_biogeochemical_tracers(biogeochemistry)..., :T, :S))) # because hydrostatic free surface will request T and S and some BGC models will too


Suggested change

tracers = unique((required_biogeochemical_tracers(biogeochemistry)..., :T, :S))) # because hydrostatic free surface will request T and S and some BGC models will too

tracers = unique((required_biogeochemical_tracers(biogeochemistry)..., :T, :S))) # because hydrostatic free surface model will request T and S and some BGC models will too

navidcy

lgtm

jagoosw added the bug Something isn't working label Jul 4, 2023

jagoosw requested a review from navidcy July 4, 2023 12:48

oops

448fde9

navidcy reviewed Jul 4, 2023

View reviewed changes

Applied review suggestions

757ae9f

navidcy changed the title ~~changed xnode use to xnodes~~ Generalize TwoBandPhotosyntheticallyActiveRadiation kernel to allow curvilinear grids Jul 5, 2023

navidcy added 3 commits July 6, 2023 16:46

code cleanup

dac9c2a

use user-facing znodes for default depth in SimpleMultiG

c3638bf

code cleanup

c530081

navidcy reviewed Jul 6, 2023

View reviewed changes

src/Boundaries/Sediments/simple_multi_G.jl Show resolved Hide resolved

navidcy reviewed Jul 6, 2023

View reviewed changes

jagoosw added 2 commits July 6, 2023 11:24

changed to native coodinates for surface PAR

9b4f2ce

tidied kelp model

a3b078b

corrected mistake

cf5c4c3

tests light over different grids and models

bfb86d9

jagoosw mentioned this pull request Jul 13, 2023

Sediment model fixes + simpler model #121

Merged

navidcy reviewed Jul 15, 2023

View reviewed changes

navidcy approved these changes Jul 15, 2023

View reviewed changes

jagoosw merged commit 1aa77d3 into main Jul 15, 2023
4 checks passed

jagoosw deleted the jsw/fix-light-for-lat-lon branch July 15, 2023 12:18

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Generalize `TwoBandPhotosyntheticallyActiveRadiation` kernel to allow curvilinear grids #120

Generalize `TwoBandPhotosyntheticallyActiveRadiation` kernel to allow curvilinear grids #120

jagoosw commented Jul 4, 2023

jagoosw commented Jul 4, 2023

jagoosw commented Jul 4, 2023

navidcy commented Jul 4, 2023

navidcy commented Jul 4, 2023

navidcy Jul 4, 2023

glwagner Jul 5, 2023

jagoosw commented Jul 5, 2023

navidcy commented Jul 5, 2023

codecov bot commented Jul 5, 2023 •

edited

Loading

jagoosw commented Jul 5, 2023

navidcy Jul 6, 2023

jagoosw Jul 6, 2023

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023 •

edited

Loading

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 13, 2023

jagoosw commented Jul 13, 2023

navidcy Jul 15, 2023

navidcy Jul 15, 2023

navidcy left a comment

	x, y = @inbounds xnodes(grid, Center(), Center(), Center())[i], ynodes(grid, Center(), Center(), Center())[j]
	x = xnode(i, j, 1, grid, Center(), Center(), Center())
	y = ynode(i, j, 1, grid, Center(), Center(), Center())

	@testset "Light attenuaiton model" begin
	@testset "Light attenuation model" begin

	tracers = unique((required_biogeochemical_tracers(biogeochemistry)..., :T, :S))) # because hydrostatic free surface will request T and S and some BGC models will too
	tracers = unique((required_biogeochemical_tracers(biogeochemistry)..., :T, :S))) # because hydrostatic free surface model will request T and S and some BGC models will too

Generalize TwoBandPhotosyntheticallyActiveRadiation kernel to allow curvilinear grids #120

Generalize TwoBandPhotosyntheticallyActiveRadiation kernel to allow curvilinear grids #120

Conversation

jagoosw commented Jul 4, 2023

jagoosw commented Jul 4, 2023

jagoosw commented Jul 4, 2023

navidcy commented Jul 4, 2023

navidcy commented Jul 4, 2023

navidcy Jul 4, 2023

Choose a reason for hiding this comment

glwagner Jul 5, 2023

Choose a reason for hiding this comment

jagoosw commented Jul 5, 2023

navidcy commented Jul 5, 2023

codecov bot commented Jul 5, 2023 • edited Loading

Codecov Report

jagoosw commented Jul 5, 2023

navidcy Jul 6, 2023

Choose a reason for hiding this comment

jagoosw Jul 6, 2023

Choose a reason for hiding this comment

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023 • edited Loading

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 6, 2023

navidcy commented Jul 6, 2023

navidcy commented Jul 6, 2023

jagoosw commented Jul 13, 2023

jagoosw commented Jul 13, 2023

navidcy Jul 15, 2023

Choose a reason for hiding this comment

navidcy Jul 15, 2023

Choose a reason for hiding this comment

navidcy left a comment

Choose a reason for hiding this comment

Generalize `TwoBandPhotosyntheticallyActiveRadiation` kernel to allow curvilinear grids #120

Generalize `TwoBandPhotosyntheticallyActiveRadiation` kernel to allow curvilinear grids #120

codecov bot commented Jul 5, 2023 •

edited

Loading

jagoosw commented Jul 6, 2023 •

edited

Loading