MAGMA2 SPH

[rennehan]

I've added the MAGMA2 SPH scheme with improved gradient estimation and 2nd order reconstruction into SWIFT. The implementation follows Rosswog 2020 extremely closely except for a few important changes.

• When the estimated G_ij vector is severely misaligned with the particle separation vector, the interaction falls back to a Gasoline2-like SPH to improve stability in highly anisotropic regions. That is set by the const_viscosity_cosine_limit compilation parameter in hydro_parameters.h.
• If the G_ij vector points in the opposite direction of the particle separation vector, it automatically falls-back to Gasoline2-like SPH.
• When the SPH falls-back to the Gasoline2-like SPH, the viscosity and conductivity are reduced by const_fallback_reduction_factor to prevent overly dissipative flows.
• There is a Balsara limiter for viscosity on top of the implementation of Rosswog 2020 or else the scheme will not work with cosmological simulations.
• There is a condition number limit on the matrix inversions from Rosswog 2020 set by const_condition_number_upper_limit. That number has been selected for stability in floating-point calculations and should not be changed from 60. If any of the matrices go above that condition number, the scheme falls-back to a Gasoline2-like SPH.
• The hydro calculations can be done in double or floating-point precision by compiling with hydro_props_use_double_precision or without it, respectively. It works fine in floating-point precision with the condition number above.
• The number of neighbors must be hard-coded into hydro_parameters.h for now because of the Van Leer slope limiting procedure in Rosswog 2020 and can be set with const_kernel_target_neighbours. I added in target_neighbours under the SPH section of the parameter file so that it may be set directly instead of resolution_eta.
• All of the extra debugging information (condition numbers, full matrices, etc.) can be output to disk by configuring SWIFT with --enable-magma2-debugging-checks.
• It is possible to use the improved gradient G_ij in an interaction but projected along the particle separation vector instead of the exact dot-products used in Rosswog 2020 by compiling with hydro_props_use_radial_artificial_terms in hydro_parameters.h. That should not be necessary, but is good for testing.
• All of the signal velocities are computed using the raw velocities and separation vectors, not any second order reconstructed quantities. It is very important that this remains the same, or the scheme is unstable.
• There is a time-step limiter that forces time steps to track at least 25% changes in specific internal energy due to conductivity only to ensure the differential equation is solved correctly.

The hydro tests work well and give good results in 3D and 2D. The results require at least 114 neighbours in 3D for the Wendland C2 kernel. I recommend values of const_viscosity_alpha=2 and const_viscosity_beta=2*const_viscosity_alpha for the viscosity. For conductivity, I recommend no higher than const_conductivity_alpha=0.075 or there will be severe spreading of shocks (see the SedovBlast_3D test).