How to implement multiple sponge layers to the same variable

[Dahaipang6561]

Hi, as the title suggests, I am struggling to implement multiple sponge layers to the same variables u and v to make them damped near the edges of x and y. I have got a reference for implementation from a previous discussion here, and the code is running. The problem now is that it is taking forever to run. I am wondering if I have done anything wrong in implementing those sponge layers and if there are better things I can do for the same purpose.

damping_rate = 1/100
xstart_mask = GaussianMask{:x}(center=x[1], width=grid.Lx/32)
xend_mask = GaussianMask{:x}(center=x[end], width=grid.Lx/32)

ystart_mask = GaussianMask{:y}(center=y[1], width=grid.Ly/32)
yend_mask = GaussianMask{:y}(center=y[end], width=grid.Ly/32)

xstart_sponge = Relaxation(rate=damping_rate, mask=xstart_mask, target=0)
xend_sponge = Relaxation(rate=damping_rate, mask=xend_mask, target=0)
ystart_sponge = Relaxation(rate=damping_rate, mask=ystart_mask, target=0)
yend_sponge = Relaxation(rate=damping_rate, mask=yend_mask, target=0)

uv_forcings =  (xstart_sponge, xend_sponge, ystart_sponge, yend_sponge)

forcing = (; u = uv_forcings, v = uv_forcings)



#set up the governing equations
model = NonhydrostaticModel(; grid,
                            buoyancy=BuoyancyTracer(),
                            coriolis=FPlane(f=1),
                            advection=UpwindBiased(order=5),
                            tracers=(:P, :b),
                            forcing=forcing)

[glwagner]

Does this setup compile on GPU? Wondering if there's an inference issue.

PS for code blocks use triple backticks not just one (```)

[glwagner]

I'd like to change this into an issue, since I think you have found a real limitation that motivates changing the source code. Let me know if you're alright with that. And thanks for raising it!

[Dahaipang6561]

Of course I am alright with it!

[Dahaipang6561]

For 1,2,3 forcings on CPU the code runs and there are results, but it takes roughly 10 times longer compared to without the sponge layers, so on CPU it could be just it is very computationally heavy that I am not getting anything even after an hour. I have no luck with GPU runnings so far, it doesn't give me anything.

[glwagner]

Ok, good to know. Let us know if writing ContinuousForcing manually works.

[Dahaipang6561]

Hi I was just trying to get it run on both gpu and cpu to see what happens. In the end I used the following codes(and a similar version for cpu)

x0 = CUDA.@allowscalar x[1]
xN = CUDA.@allowscalar x[end]
y0 = CUDA.@allowscalar y[1]
yN = CUDA.@allowscalar y[end]

const parameters = (dx = grid.Lx / 32, x0 = x0, xN = xN, dy = grid.Ly / 32, y0 = y0, yN = yN, r = 1 / 100)

@inline function mask(x::Float64, y::Float64, p)
    return exp(-((x - p.x0)^2) / (2 * p.dx^2)) +
           exp(-((x - p.xN)^2) / (2 * p.dx^2)) +
           exp(-((y - p.y0)^2) / (2 * p.dy^2)) +
           exp(-((y - p.yN)^2) / (2 * p.dy^2))
end

@inline function u_sponge(x::Float64, y::Float64, z, t, u, p)
    return -p.r * mask(x, y, p) * u
end

@inline function v_sponge(x::Float64, y::Float64, z, t, v, p)
    return -p.r * mask(x, y, p) * v
end

u_forcing = Forcing(u_sponge, field_dependencies = :u, parameters = parameters)
v_forcing = Forcing(v_sponge, field_dependencies = :v, parameters = parameters)

forcing = (; u = u_forcing, v = v_forcing)

Both of them actually run and shows outputs for simulation steps, but it creates out of memory issue on HPC, and roughly speaking, the time required is tripled/quadrupled. For the reference, I am using a 25625632 grid with 0.01 timestep and the total simulation time is 5 minutes.

I am wondering if that is an issue with my code? Thank you for your help and patience:)