Overflow encountered in Multiply error for SWE IVP Sphere

[Rick Bonhof]

Hello all,

For my current MSc thesis I am working on a shallow water simulation of Jupiter's top atmospheric layer, using Dedalus3 and data from JWST for the zonal jet profile of Jupiter.

As a basis, I use the Spherical shallow water (sphere IVP) example script, with changes specified as per the attached script. However, outside of trying the test case of galewsky and one alternate test case with an opposite flow of galewsky's test case on the other hemisphere, all simulations experience "RuntimeWarning: overflow encountered in multiply np.multiply(arg0_exp_data, arg1_exp_data, out=out.data)", with ringing effects on the snapshot where the crash occurs like in the attached gif. This error tends to occur only 4 seconds into the simulation, with slight variation depending on restrictions of the data (for example, truncating the function or decreasing/increasing the resolution of the simulation).

I think the crash here occurs through Gibbs ringing based on similar-looking problems posed on this forum, but I have not been able to fix this yet. Currently I am working on testing the stable cases of opposing galewsky flows while changing the resolution and timesteps, to find a limit of the simulation where it still would remain stable to hopefully give insight into the cases where the Overflow error is encountered.

Are there any problems in my formulation or structure of my code? Any help would be greatly appreciated, as creating a stable version of the simulation has been a major roadblock up until now. If any further information about the simulation (e.g. results or plots of the solver) please let me know.

Thank you in advance.

Kind regards,

Rick Bonhof

[Geoffrey Vasil]

Hi Rick, 

Have you tried lowering your time-step size until (if) things run? I like that you’ve used common sense unit choices. But it’s hard to tell if the parameters are reasonable for 1sec time steps. You also have the comment "Diffusion matched at ell=32”. I assume you tuned this so diffusion is important for ell > 32. That seems reasonable if so. 

But given you have all the linear physics on the LHS, a lot comes down to the size of your initial condition. Without doing a lot of experimenting it’s hard to say what’s going on. Again I would try to lower the time step to try to make things work. The spatial resolutions seem pretty good. 

Cheers,

Geoff  

--
You received this message because you are subscribed to the Google Groups "Dedalus Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dedalus-user...@googlegroups.com.
To view this discussion visit https://groups.google.com/d/msgid/dedalus-users/943e2517-97e3-4d0d-b559-5d0a983292a5n%40googlegroups.com.
<jupiter_profile_truncated Upper Crash.gif><jupiter_fourier.py>

[Rick Bonhof]

Hi Geoff,

Thank you for your reply! It is good to hear that the unit choices are sensible, along with the choices for diffusion. For experimenting, I have tried different initial setups like increasing the resolution of the quadrature grid to double or reducing it to half, truncating the zonal jet profile from 70 degree positive and negative latitude to 40 degrees, with no success.

Currently I am able to create a far more stable simulation using multiple gaussian curves to approximate the zonal profile, but this is incredibly simplified (to 7 distinct flows instead of the complex zonal flows observed from Jupiter). As a result, the method with an analytical function constructed from the data with a DFT is more desirable as it would come closer to an actual shallow water simulation of the top atmospheric layer.

 

With a timestep of 1 millisecond, the simulation experiences the same Runtimewarning at approximately 20 milliseconds. Seeing as the end goal of the simulation is to run the jupiter zonal jet profile for hours, such that fluid elements can properly travel across the surface, going far smaller than the scale of milliseconds would create an extremely time-consuming simulation. I could attempt an even smaller size for timesteps, but would that be still within reason at this point? 

If you have any other ideas as to what could help avoid the issues I have encountered, they would be greatly appreciated.

Cheers,

Rick

[Rick Bonhof]

Hi everyone,

The gaussian case is still stable, but this still is a very rudimentary case of manually creating the profile instead of being able to adapt any data and return a profile (which would allow for more test cases rather than one manually fit for just jupiter).

Would there still be something that is incorrect as per the previous messages in this thread? The gaussian case right now seems to be incredibly stable, without perturbations it does not show any instabilities for timesteps from 10 seconds to 100 seconds, and from a resolution of 62x128 all the way to 1048x2048 (for 2048x4096 the HPC I use fails the job due to now being able to allocate enough memory to it). Current testing is done such that I can test extreme cases like near tidal locking and large scale height field perturbations, but smaller scale perturbations also do not appear to cause any instabilities.

Thank you in advance.

Cheers,

Rick