Sponge layer not behaving as expected #1488
-
I set-up a sponge layer for my simulation following the examples here, but I'm getting results that are different from what I expected. Here's how I'm setting my grid, model and sponge:
From the docs, I expected this to
which leads me to believe that there's some a sharp cut-off in the function mask so that velocities above the range However, when I plot the vertical velocities (and other variables with zero mean) I get the results below, which indicate that the sponge layer is active in the approximate range Am I missing something here? Thanks! |
Beta Was this translation helpful? Give feedback.
Replies: 25 comments 1 reply
-
Can you plot bottom_mask to make sure that looks like what you expect it to be? |
Beta Was this translation helpful? Give feedback.
-
I honestly don't know how to do that. It's type is |
Beta Was this translation helpful? Give feedback.
-
Ah not sure we can plot
might not be the best description of I guess if you want to reduce the spatial effect of the sponge you could reduce the |
Beta Was this translation helpful? Give feedback.
-
Sorry @tomchor for my suggestion. It might be nice to be able to plot the mask so the user knows what they are actually doing. Do you people agree? If yes I can create an issue asking whether we want to do this or not. |
Beta Was this translation helpful? Give feedback.
-
From relaxation.jl, I see that GaussianMask in the z direction is defined as follows:
I see that your center is Maybe you can define this function with replacing z with the znodes and then plot it to make sure everything looks right? |
Beta Was this translation helpful? Give feedback.
-
That could be a pretty neat feature (and helpful for debugging). Might be easiest to do by adding a plotting recipe (#1115)? |
Beta Was this translation helpful? Give feedback.
-
This behavior might be expected. Within 100 m of the bottom boundary, the Gaussian has a value julia> exp(-0.5 * 1^2)
0.6065306597126334 so I think we expect the solution to be strongly damped. within 200 m, it falls to 14% julia> exp(-0.5 * 2^2)
0.1353352832366127 Which is still significant. At 300 m its about 1%, so it probably has a negligible impact above -200 m. |
Beta Was this translation helpful? Give feedback.
-
I couldn't agree more. It's amazing how much you can see when you look. |
Beta Was this translation helpful? Give feedback.
-
Thanks for the help guys. I also agree with @francispoulin that plotting the mask would make it very helpful! Here are some thought that I had that might be worth discussing: 1 - Sponge layers (afaik) should only alter the NS equations in a finite part of the domain, right? Otherwise, even far from the relaxed edge, you're not really solving the NS equations --- it's being modified by an extra unphysical forcing. So one question is: is this really the behavior you want your sponge layer to have? (And here I'm assuming you're using relaxation as a synonym for sponge, because that's what the docs imply.)
Hopefully that makes sense. |
Beta Was this translation helpful? Give feedback.
-
@tomchor, indeed you are correct on 1. above, but let me remark that sponge layers could be a bit tricky. I explain: What you want sponge layers to do is to allow flow propagation into the sponge region but where drag would dissipate any fluid motion there. The "tricky" part is to make the sponge layer transition smooth enough so that the fluid does not experience it as a "wall". If the transition region is very short then the fluid "sees" the sponge as a wall and waves are reflected back into the fluid. But how "short" is short-enough depends on the scales of the fluid motion... In my experience always there is some fiddling to be done for every problem or even for same problem at different parameter regimes. |
Beta Was this translation helpful? Give feedback.
-
Oceananigans.jl/src/Forcings/relaxation.jl Lines 81 to 82 in 35f749e The snippet shows that So, sponge layers are simply one application of @tomchor, perhaps what you're saying is that a Gaussian is not appropriate for a sponge layer |
Beta Was this translation helpful? Give feedback.
-
I admit I haven't read @tomchor script above in detail before making my comment... From my experience, sponge layers are usually applied as linear drag with a variable linear coefficient. E.g., ∂u/∂t = ... + μ(z) u where μ(z) is a function that is 0 everywhere and becomes 1 in the sponge with some transitional region of certain width. As, @glwagner points out, |
Beta Was this translation helpful? Give feedback.
-
@glwagner yes, that was more along the lines of what I was saying. From the docs it seemed to me using a Gaussian mask in the relaxation should be equivalent to a sponge layer, which is why I expected the gaussian function to be set-up with a sharp cut-off. But I think we're on the same page about what a sponge layer should be. So I think it may be a matter of maybe making the docs clearer? Maybe also, like you said, constructing an example that uses sponge layers (I think there's an issue about there somewhere, right?). The Langmuir example might a good one. In the original paper by Jim they use a radiation boundary condition at the bottom which doesn't have an analog in the Langmuir example at the moment. |
Beta Was this translation helpful? Give feedback.
-
Based on the answers here I tried to set my own mask function. It appears to run fine on the CPU but when I try to run it on the GPU I get this error when running the simulation (with a few subsequent errors after these lines, but this is the first one):
Am I missing something here? Like I said, it runs fine on the CPU. Thanks! |
Beta Was this translation helpful? Give feedback.
-
When running GPU kernels, KernelAbstractions.jl uses Casette.jl (which defines
I assume there was a a function it failed to convert to be GPU-friendly. But the only function I see in there is your @tomchor Can you post or share your full script? Might help with debugging. |
Beta Was this translation helpful? Give feedback.
-
@ali-ramadhan Thanks for the info. Here's the my code: https://pastebin.com/ACKJhEj7 Thanks! |
Beta Was this translation helpful? Give feedback.
-
Thanks @tomchor! I think the issue was that the GPU was not able to access the In general, for functions that are to be embedded within GPU kernel (like your So with this change I was able to run your script on the GPU: # Set-up sponge layer
heaviside(x) = ifelse(x < 0, zero(x), one(x))
const H = grid.Lz
sponge_one = minimum(oc.Grids.znodes(Face, grid))
sponge_zero = sponge_one + grid.Lz/10
function bottom_mask_func(x, y, z)
sponge_one = -H/2
sponge_zero = sponge_one + H/10
return heaviside(-(z-sponge_zero)) * (z - sponge_zero)^2 / (sponge_one - sponge_zero)^2
end
#bottom_mask = GaussianMask{:z}(center=-grid.Lz/2, width=grid.Lz/40)
mom_sponge = Relaxation(rate=1/20, mask=bottom_mask_func, target=0) You could also declare |
Beta Was this translation helpful? Give feedback.
-
Thanks! Indeed that seems to work 👍 Apparently running on GPUs has a few caveats even after you guys did most of the hard work in Oceananigans |
Beta Was this translation helpful? Give feedback.
-
Awesome! Yeah there are some GPU idiosyncrasies (biggest one being that GPU errors can be very obscure) with user-defined forcing functions and boundary conditions. Should this issue be closed or do you want to modify the title to reflect that it's about improving the language in the documentation? |
Beta Was this translation helpful? Give feedback.
-
Might be good to clarify some of this discussion: a "sponge layer" refers to any region in a simulation domain that dissipates energy. A wide variety of forcing functions are "valid" sponge layers, including the example in the documentation. The Gaussian mask appears to be behaving correctly, and as one expects given the documentation, source code, and known/measurable/plottable properties of Gaussians. Using a Gaussian to define the masked region is valid strategy used successfully in research. Of course, this does not mean that other masking functions may be more appropriate, or that a "smooth step" function might be more useful to users, since the width of the sponge, and width of the sponge layer transition region may be modulated independently. This is not the case for the Gaussian, which has only one width parameter. A sharp heaviside is not a good choice for a sponge layer when the sponge is intended to absorb radiating internal waves. The reason is that waves can reflect off a sponge layer if the transition between undamped and damping regions is too abrupt (as pointed out by @navidcy). was under the impression this might even be an issue for problems without waves... ? The requirement that sponge layers involve smooth transitions does mean that sponge layers can take up a substantial portion of the computational domain. This is often cited as a downside of sponge layers (eg Klemp and Durran 1993). |
Beta Was this translation helpful? Give feedback.
-
@tomchor can you clarify what you mean when you say you expected the "gaussian function to be set-up with a sharp cut-off"? Do you mean you expected the width of the Gaussian to be narrower? We can change the example in the docs to have a narrower width; perhaps
To clarify: the implementation in the documentation is a valid implementation of a sponge layer, correct? |
Beta Was this translation helpful? Give feedback.
-
@tomchor the reason your code did not compile is because this function function bottom_mask_func(x, y, z)
sponge_one = -grid.Lz/2; sponge_zero = sponge_one + grid.Lz/10
return heaviside(-(z-sponge_zero)) * (z - sponge_zero)^2 / (sponge_one - sponge_zero)^2
end references One strategy to avoid having to use It might be nice to add a |
Beta Was this translation helpful? Give feedback.
-
@glwagner Thanks for the clarifying the global-var-limitation for GPUs! This mat be something to add to the docs at some point (I didn't see it there!) (I'll answer the sponge layer questions on a separate answer in a minute!) |
Beta Was this translation helpful? Give feedback.
-
No, I mean I expected the sponge layer to be confined to a finite part of the domain. That is, I expected the sponge layer to be identically zero for most of the domain. Us not "being on the same page" here (for the lack of a better description! haha) may be due to different philosophies of what a sponge layer should be. So, due to my background, a sponge layer should be a function that only acts in a finite part of the domain. So, for example, in my code I set-up the sponge layer mask as
with H(z) being the Heaviside function. Comparing this to the Gaussian mask, this has a disadvantage of not being as smooth (although it should be smooth enough because of the (z - z₀)² / (z₁ - z₀)² term, and we can always increase the order of the exponents). However it has a very important advantage: I know exactly when it stops influencing my solution. So I know that in the range I was always taught (and I still agree) that setting up a sponge layer mask as something that only reaches zero at infinity (like a Gaussian) is not recommended since you don't know when do start believing your solution. In practice if you're 4σs away from the center of the Gaussian you're pretty sure that the sponge layer influences are minimal (but they're not exactly zero!). But a harder question is: how far away is far enough? A note here should be that if my mask was simply So, again, I think it boils down to your background and modeling philosophy! But it would be nice to have options for finite-range masks, such as the one I imposed, as well. (Although there are many other functions that are used and mine is probably the easiest one you can come up with...) I hope I made the issue clearer! And sorry about any confusion. |
Beta Was this translation helpful? Give feedback.
-
@tomchor Not sure if this "issue" was resolved but seemed more suitable as a discussion so I converted it to one. |
Beta Was this translation helpful? Give feedback.
@tomchor the reason your code did not compile is because this function
references
grid.Lz
as a global variable. Global variables can only be used inside functions on the GPU if they are declaredconst
as @ali-ramadhan has done.One strategy to avoid having to use
const
is to build callable objects likeGaussianMask
that store their parameters. This is why a parameterizedGaussianMask
cal…