FourierFlows / GeophysicalFlows.jl

Geophysical fluid dynamics pseudospectral solvers with Julia and FourierFlows.jl.
https://fourierflows.github.io/GeophysicalFlowsDocumentation/stable/
MIT License
153 stars 31 forks source link

MultilayerQG example crashed at about 35000 steps #265

Closed szy21 closed 2 years ago

szy21 commented 2 years ago

I just run examples/multilayerqg_2layer.jl on a fresh clone on Julia 1.6.0, but with nsteps = 40000. The simulation crashes at about 35000 steps. Here is the plot of upper and lower q in the last output and there seem to be some numerical instability.

I think adding some hyperviscosity might fix it, but just want to know if this is as expected? image

szy21 commented 2 years ago

I tried using the default de-aliasing and it didn't crash, but cfl is still larger than 1 towards the end of the simulation (and I guess the simulation would crash if I run it for a longer time period). I also tried adding a hyperviscosity, and a large enough hyperviscosity can stabilize it, but I don't know what an appropriate hyperviscosity should be. I wonder if there is a recommended setup?

glwagner commented 2 years ago

Ah the never-ending quest for "recommended hyperviscosity"!

There's no such thing (@navidcy calls this an "open research question")...

But, if you can estimate a time-scale (eg an eddy turnover time or forcing time scale, using the input parameters of the problem) then you can write the hyperviscosity as

C = 1e-2
nu4 = C * dx^4 / time_scale

where dx is the grid spacing and C is a free parameter. In my experience we typically choose a pretty small number for C. This fixes the rate of decay of motions at the grid scale as some fraction of the energy injection / nonlinear time-scale.

Another scaling argument might use dt as the time-scale. The scaling ends up similar because dt is also set by the eddy turnover time-scale (ideally).

If you don't care about energy budgets (or analytical PDEs) then you can use a filtered time-stepper like FilteredRK4, which is a kind of LES scheme that mops up energy at the smallest scales, thereby stabilizing your simulation if the energy cascade isn't too dramatic.

szy21 commented 2 years ago

Thanks @glwagner! In the example, dealiasing is turned off by setting aliased_fraction=0. I wonder if this is what I should do? Or is it better to use the default dealiasing?

navidcy commented 2 years ago

Is it a Filtered timestep? If so aliasing won't make almost any difference.

Try reducing the time-step? Perhaps the example was not run long enough to equilibrate and thus flow speed keeps rising?

szy21 commented 2 years ago

Yes, I was using FilteredRK4. I'll try to reduce the time step, thanks!

szy21 commented 2 years ago

I reduced the time step by half and the simulation is stable. I will close this issue then, thanks @glwagner and @navidcy!

navidcy commented 2 years ago

I'll leave it open just as a reminder to reduce the tilmestep in the example. :)