CliMA / Oceananigans.jl

🌊 Julia software for fast, friendly, flexible, ocean-flavored fluid dynamics on CPUs and GPUs
https://clima.github.io/OceananigansDocumentation/stable
MIT License
965 stars 190 forks source link

Weighted WENO stencils for stretched grids #1704

Closed glwagner closed 2 years ago

glwagner commented 3 years ago

In A simple algorithm to improve the performance of the WENO scheme on non-uniform grids a simple method for generalizing WENO stencils for stretched grids is described. This scheme is easy to implement and limits to what we already have implemented in the case of a regular grid.

glwagner commented 3 years ago

To implement this scheme we first need to reformulate all the stencils in https://github.com/CliMA/Oceananigans.jl/blob/master/src/Advection/weno_fifth_order.jl in terms of interpolation operators.

The next step is to use weighted interpolation operators on stretched grids according to this definition:

image

where h is the grid spacing.

"That's it."

tomchor commented 3 years ago

Any clue as to what amount of error we're introducing by using WENO5+StretchedGrids? Maybe enough to put an approx upper bound to the level of stretching we should use? I'm assuming this only becomes important when the grid stretching is too aggressive but who knows...

Also I wonder how this plays out with https://github.com/CliMA/Oceananigans.jl/issues/1705

glwagner commented 3 years ago

No clue, but the reference linked in the original post tests some cases and includes a discussion.

francispoulin commented 3 years ago

If we wanted to quantify the difference, we can simply find an exact solution and simulate it and compute the error. There are lots of examples to choose from.

glwagner commented 3 years ago

Here's a couple results picked from the paper:

Stationary vortex on a grid with "random" spacing

image

Stationary vortex on a curvilinear distorted grid

image

Contour plots of a perfectly circular vortex are good at highlighting the differences. Their "production" case of flow around an airfoil is more subtle:

image

glwagner commented 3 years ago

In the above, WENO-JS is the scheme we have now, and WENO-NM is the "weighted" scheme they propose.

The weighted scheme is extremely simple fairly simple. It's just a little algebra to get it to work for us. On reading the paper again, I realize that the smoothness indicators are a bit annoying. But doable.

francispoulin commented 3 years ago

Thanks @glwagner , Certainly Figure 5 look sdifferent and that is good inspiration to consider how to impliment WENO-NM.

tomchor commented 2 years ago

@glwagner I've been thinking about the WENO5 algorithm based on our discussion in https://github.com/CliMA/Oceananigans.jl/discussions/2054, and I just wanted to ask a quick question to make sure I understand the problem in this issue.

If I understand correctly, the WENO5 scheme calculates the advection with different advection schemes (all of which are correct regardless of the grid) and then does a weighted sum of these solutions according so some criterion. So the issue with WENO5 and StretchedGrids is not that the advection calculated is "wrong", it's just that the weighting of the solutions is done in a non-optimal manner, and therefore introduces errors that could be avoided, right?

Basically I'm asking the question: is it the case the shouldn't use WENO5 with stretched grids because the answer is wrong, period? Or we can use it, it's just not gonna be as accurate as it potentially could be?

glwagner commented 2 years ago

The results in A simple algorithm to improve the performance of the WENO scheme on non-uniform grids show that WENO reverts to 2nd-order accuracy on stretched grids. So its not "incorrect", just less accurate.

However, for any given problem, the convergence rate is only one aspect of the accuracy of a solution. It's possible that WENO5 is still more accurate than any other numerical scheme, even if the solution only converges to the true solution at a 2nd order rate.

It also seems likely that this is true for any reconstruction stencil. I'm not sure how much attention is paid to this issue in other codes.