Creating a new doc example with drag on immersed boundaries

[tomchor]

Given that the doc examples are gonna run much faster due to https://github.com/CliMA/Oceananigans.jl/pull/3135 (I'm assuming it's gonna get merged soon) should we create an example with an immersed boundary with drag forces?

My main reasons for suggesting this are:

1. Probably there's a lot of user interest in this application, but no example showing this capability (https://github.com/CliMA/Oceananigans.jl/pull/3132 will add an ImmersedBoundaryGrid example, but with default boundary conditions);
2. Coming up with an intuitive way to impose flux boundary conditions on immersed boundaries is hard (see for example https://github.com/CliMA/Oceananigans.jl/pull/3142 and discussion in https://github.com/CliMA/Oceananigans.jl/pull/2437), so I think it'd be good to give users a guaranteed-to-work easy-to-find example to base themselves off of.
3. I also think having such examples probably makes counter-intuitive definitions and bugs more likely to be caught (and therefore fixed) early, since it makes the ImmersedBoundaryGrid API more visible to users, who will be more likely to use it for explorations. So my opinion is that including an example will increase our confidence in the IBG implementation (although, again, this is already being done with https://github.com/CliMA/Oceananigans.jl/pull/3132, albeit with default BCs).

If there's interest we can discuss what this example should look like.

Thoughts?

[navidcy]

Perhaps a side note: is there a benefit in adding drag to the internal tide example? I don't thing it would make any difference on the physics...

[glwagner]

I agree, it seems best simply to add drag to the internal tide example for the time being.

[glwagner]

The main rationale for not having an immersed boundary conditions example was because we regarded it as an experimental feature, and didn't want to trip new users up (which is the opposite of the motivation #3). We could change our philosophy though.

@tomchor you may not realize it, but you probably have the most experience with immersed boundary conditions. I think it's reasonable that nobody has added an example for these yet. I don't think there are any validation tests at all, or papers that use them, in constrast to much of the rest of the code. The fact that a major bug like the one found in #3142 still exists is a testament to the fact that these are unvalidated, and I think also validates our justification for not having an example yet.

[tomchor]

I guess imo it's a bit of a chicken and egg problem, since bugs are found when people start use the feature. And we don't wanna release the feature until bugs are found. So my thinking in suggesting an example is an attempt to break the cycle.

But it's a good point that the first users of a feature that's still unvalidated should probably more "experienced" ones.

Since now there's an example in the works with immersed boundaries, does that mean that that feature is not considered experimental anymore?

[glwagner]

I don't think there's a chicken and egg problem. Testing new features requires some initiative and also experience with / understanding of the relevant science applications. I think the "egg" is the validation test itself. A groomed, pedagogical docs example is the "chicken" that emerges from that egg. I also think that validation is not just about ensuring quantative accuracy, but about working out the user interface. We want the docs implementations to reflect relatively mature user interfaces, since the docs examples are relatively influential pieces of code.

I believe immersed boundaries are still experimental, not least in light of #3142. What are the quantitative validation tests that use immersed boundaries? The only one that I know if is the one I worked on, which shows that a stratified ocean remains at rest:

https://github.com/CliMA/Oceananigans.jl/blob/main/validation/immersed_boundaries/resting_stratified_bumpy_ocean.jl

We haven't documented those results well, though. The internal tide case could be adapted to validate immersed boundaries with the hydrostatic model, leveraging the vast literature on that subject. I don't think the example involves a comparison with theory or other published numerical results though. It wouldn't be that much work to turn that case into a validation test.

As for NonhydrostaticModel, by all accounts it seems that it should be experimental until we have a validated pressure solver. It could be interesting to see if we can validate the nonhydrostatic model for certain LES cases, somehow. Otherwise, we do know that pressure solver is in general incorrect (either impenetrability across immersed boundaries is not satisfied, or the velocity field is divergent close to the immersed boundary).

[tomchor]

I also think that validation is not just about ensuring quantative accuracy, but about working out the user interface. We want the docs implementations to reflect relatively mature user interfaces, since the docs examples are relatively influential pieces of code.

That's where the chicken and egg problem comes from IMO. It's harder to have a mature, friendly user interface without exposing it to users and course-correcting from input. But that's only my personal opinion and that's obviously up for debate.

I'll close this issue we apparently aren't ready for what I proposed.

[glwagner]

Do you want to close this now? We can have the internal tide example close it...

[glwagner]

It's harder to have a mature, friendly user interface without exposing it to users and course-correcting from input.

The process of setting up a validation test involves interacting with the interface as a user. So requiring a validation test does expose an interface to users.

[navidcy]

True.

@tomchor would you be keen to add drag in the internal tide example? If you just add the code I can deal w the rest.

[tomchor]

True.

@tomchor would you be keen to add drag in the internal tide example? If you just add the code I can deal w the rest.

I can, but I don't really see drag there adding value to the example. @navidcy, do you think it would be relevant to the physics? Or do we just want drag to be there as an example of using the interface?

Also, I'm not very familiar with these larger-scale, hydrostatic simulations. What's the physically-preferred way to add drag there? Does something like

\vec\tau = -C^D \vec u \sqrt{u^2+v^2}

(applied to bottom, west and east, with $C^D \approx$ 1e-3) look reasonable?

[navidcy]

Yeah, a quadratic drag line what you wrote looks good!

But what do you mean west/east? (The internal tide example is periodic in x, perhaps you thought it was bounded?)

[navidcy]

Indeed, I don't think drag is so relevant physically for this problem but it should demonstrate the interface.

[navidcy]

But it's a start? We can modify examples later on...

[simone-silvestri]

maybe apply it only to bottom, since for hydrostatic simulations, the spacing is strongly anisotropic (1:250 ratio in the internal tide example). The same drag law imposed on sides would result in a very large force.

[glwagner]

I don't think that's true --- side-wall drag is negligible for thin aspect ratios, because the horizontal spacing is huge.

I think we should illustrate imposing a drag law that would work for both thin and O(1) aspect ratio (drag on all boundaries, bottom and side).

[glwagner]

This might work (for example, for drag on the u component):

@inline u_drag(x, y, z, t, u, v, w, cd) = - cd * u * sqrt(u^2 + v^2 + w^2)
u_immersed_bc = FluxBoundaryCondition(u_drag, field_dependencies=(:u, :v, :w), parameters=cd)

[glwagner]

And right, for thin aspect ratios (H/L << 1), we typically neglect the dependence on w for two reasons: (1) w is small (eg w ~ u * H/L), and (2) w is non-zero only on side-walls, which make a negligible contribution to the overall flux, because side-wall areas are O(H/L) smaller than vertical-normal areas.

For this example, we can either use a hydrostatic-specific drag, or we can use a more general formulation and comment on how it could be simplified using an assumption of thin-aspect ratio (which is also the basis for the hydrostatic approximation).

In the context of Oceananigans, it turns out to be simpler to specify a 3D drag. So the "simplification" for thin-aspect ratio isn't as important as it is, for example, for analytical / theoretical calculations.