Sponge layer implementation

[suhanidalal]

Hi!

I'm trying to implement a sponge layer to dampen tracers approaching the north/south latitude boundary of my near-global simulation. Referencing this discussion, I tried creating the mask like this:

const sponge_size = 2.0    
const slope = 1.0         
const ymin = minimum(y)
const ymax = maximum(y)
x, y, z = nodes(grid, (Center(), Center(), Center()))


@inline mask_func(x,y,z) = ((
     tanh((y-(ymax-sponge_size))/slope)
    *tanh((y-(ymin+sponge_size))/slope)
)+1)/2

rate = 1 / 60

sponge_SC = Relaxation(rate=rate, mask=mask_func, target=0)

I also tried to implement the mask with LinearlyTaperedPolarMask, but that didn't work either:

mask = LinearlyTaperedPolarMask(southern=(-70, -65), northern=(65, 70));

The error I'm getting seems to be related to the GPU compiler and part of it states:

InvalidIRError: compiling MethodInstance for Oceananigans.Models.HydrostaticFreeSurfaceModels.gpu_compute_hydrostatic_free_surface_Gc!(::KernelAbstractions.
...
Reason: unsupported call to an unknown function (call to ijl_get_nth_field_checked)

Do you know what I'm doing wrong? Is there a different way to create the mask? Thank you!

[glwagner]

Does it run on CPU? That is the first thing to check.

It may be helpful to see your whole script as well as the entire stack trace.

[suhanidalal]

That was definitely my problem – it was an OffsetVector. I just changed it and it's working now. Thank you! Also, the simulation still crashes at DatasetRestoring on GPU, but it works fine on CPU. Any thoughts on how to resolve that? Thanks again.

[glwagner]

For the future I suggest using Field rather than raw arrays, because then the data will be automatically allocated correctly on the architecture you're using.

As for DatasetRestoring --- I don't see any problem at first glance. Is the GPU error a failure to compile compute_Gc? Does it produce the expected result on CPU? (Eg not merely running, but running correctly.) You can try launching julia with julia --check-bounds=yes to see if there is some issue there. Did you try changing the inpainting?

[suhanidalal]

Got it!

The GPU error is:

Scalar indexing is disallowed.
Invocation of getindex resulted in scalar indexing of a GPU array.
This is typically caused by calling an iterating implementation of a method.
Such implementations *do not* execute on the GPU, but very slowly on the CPU,
and therefore should be avoided.

It does produce the expected result on CPU and I did try changing the inpainting.

[glwagner]

I'm at a loss then. There is a time-stepping test here:

ClimaOcean.jl/test/runtests_setup.jl

Line 169 in 4216d36

T_forcing = DatasetRestoring(T_meta, arch; inpainting, rate=1/1000)

which I think runs on GPU. Not sure why it doesn't work for you.

@navidcy

[simone-silvestri]

I think we are hitting now the same error in the tests https://buildkite.com/clima/climaocean-ci/builds/4379#01980802-5ac2-4dae-9e32-a3f2274dc5d1