Castro
Castro copied to clipboard
AMReX-Astro
Remove grav_source_type and rot_source_type parameters
PR summary
The current default values (4 for both; the conservative energy formulation) are now the only option.
Closes #156
PR checklist
- [x] test suite needs to be run on this PR
- [x] this PR will change answers in the test suite
- [ ] all functions have docstrings as per the coding conventions
I would like to actually do this because it is my goal to implement fully conservative gravity and rotation options that are applied directly to the fluxes (for both momentum and energy). When I originally did the analysis of these source terms, I didn't see much benefit from doing it this way in terms of conservation properties, but I think it will be much better for code maintainability.
flame_wave uses grav_source_type = 2 currently. We need to see if the problem works well with 4. I vaguely remember there being issues with HSE, but that was a few years ago.
For constant gravity, the source term is very similar in both cases, except that the energy source term evaluates the time-centered (rho * v * g) at zone edges rather than zone centers.
The main case in which they would differ substantially is at a symmetry boundary, because the contribution on the boundary is zeroed out, as the fluxes are zeroed out.
Note that the MOL and SDC options still use an explicit source term, not the conservative flux-based approach. Also, it is not clear how to do this conservative approach to 4th order accuracy.
grav_source_type == 4 only matters for the corrector source, this shouldn't change anything for MOL/SDC.
agreed, but it should be possible to do this for MOL/SDC second-order, since we have the fluxes at the time we need the source.
so castro.grav_source_type=4 crashes flame_wave. Too many subcycles.
Is this conservative formulation axisymmetric-aware?
I suspect the problem here is the hydrostatic boundary. We need to figure out how to do an HSE boundary condition with conservative gravity.
Yes, the conservative formulation works for axisymmetric/spherical geometries.
Running flame_wave with the following causes the "too many subcycles" error.
Works fine with grav_source_type = 2.
inputs.noboost.500Hz.txt probin.noboost.500Hz.txt
I was able to reproduce this after 3234 steps with
jsrun -n 8 -a 1 -g 1 ./Castro2d.pgi.MPI.CUDA.ex inputs.noboost.500Hz castro.grav_source_type=4 amr.max_grid_size=1152
The code is crashing because of CFL violations on the fine grids. They are so large that the default number of retry subcycles cannot handle it, hence the abort.
Examining the plotfile at the end of the previous step shows that all of the velocities are reasonable. The large velocity is not present there on any level. Rather, the large velocity is being generated during a regrid. Unfortunately this is diffcult to avoid: even if the regrid step is extremum-preserving, the extrema preserved are the momenta and density, not velocity (which is not a state variable), and so it's still possible to generate a large velocity as a result of an unlucky combination of (rho u) and (rho). I verified that if I set amr.regrid_int=-1 and amr.regrid_on_restart=0 on the restart then the issue went away.
My recommendation is that you use amr.compute_new_dt_on_regrid=1 so that, if large velocities are generated by a regrid scheme, we catch this and drop the timestep to compensate for it. Actually, we do have issue #356 about making it on by default, so that would be another way to resolve the issue. I checked that this option resolves the crash for this case, and the timestep quickly returns to its original value thanks to change_max.
I don't know for sure why this happens only with grav.source_type = 4. But it's most likely just bad luck; this could happen regardless of other parameters, and it's the reason I have had compute_new_dt_on_regrid enabled in wdmerger for a while.
There are still issues with flame_wave. Using the above dt fix, the problem runs out longer but eventually hits a dt of 9e-15 (sometime after 27000 steps).
The removal of the hard_cfl_limit = 0 option in this run resolves that timestep crash (#723)
But the timestep is still about 4x smaller than the grav_source_type = 2 version.
tests: http://groot.astro.sunysb.edu/Castro/test-suite/gfortran/2022-05-27-002/index.html lots of gravity compilation failures
we still have the timestep issue with this change for flame_wave -- the timestep is about a factor of 2x smaller with castro.grav_source_type = 4 as compared with castro.grav_source_type = 2
Here's a comparison of the flame_wave with different castro.grav_source_type:
castro.grav_source_type = 2:
castro.grav_source_type = 4:
The biggest difference appears to be in the temperature structure at the top of the atmosphere. This is where the density is going toward 0, so it makes sense that the mass fluxes might be doing something odd here which changes the energy.
The castro.grav_source_type = 4 run has more retries and has them due to both "timestep validity check failed;" and "invalid density; proceeding to a retry.". The castro.grav_source_type = 2 run only has it due to timestep validity. And has fewer retries overall despite evolving for 4.8x longer (in simulation time)
Note this is with CFL = 0.8, so the timestep validity retries are not unexpected.