Address review comments

docs/cluster.md

@@ -184,15 +184,15 @@ Pertubations are controlled by the following parameters:
 # for the velocity field
 sigma_v = 0.0      # volume weighted RMS |v| in code velocity; default: 0.0 (disabled)
 l_peak_v = ???     # lengthscale (in code length units) where the velocity spectrum peaks. No default value.
-k_peak_v = ???     # (alternative to l_peak_v): wavenumber in normalized units where the velocity spectrum peaks. No default value.
+k_peak_v = ???     # (exclusive alternative to l_peak_v): wavenumber in normalized units where the velocity spectrum peaks. No default value.
 num_modes_v = 40   # (optional) number of wavemodes in spectral space; default: 40
 sol_weight_v = 1.0 # (optional) power in solenoidal (rotational) modes of the perturbation. Range between 0 (fully compressive) and 1.0 (default, fully solenoidal).
 rseed_v = 1        # (optional) integer seed for RNG for wavenumbers and amplitudes
 
 # for the magnetic field
 sigma_b = 0.0      # volume weighted RMS |B| in code magnetic; default: 0.0 (disabled)
 l_peak_b = ???     # lengthscale (in code length units) where the magnetic field spectrum peaks. No default value.
-k_peak_b = ???     # (alternative to l_peak_b): wavenumber in normalized units where the magnetic field spectrum peaks. No default value.
+k_peak_b = ???     # (exclusive alternative to l_peak_b): wavenumber in normalized units where the magnetic field spectrum peaks. No default value.
 num_modes_b = 40   # (optional) number of wavemodes in spectral space; default: 40
 rseed_b = 2        # (optional) integer seed for RNG for wavenumbers and amplitudes
 ```

src/hydro/prolongation/custom_ops.hpp

@@ -41,6 +41,10 @@ using parthenon::ParArray6D;
 using parthenon::TE;
 using parthenon::TopologicalElement;
 
+// Multi-dimensional, limited prolongation:
+// Multi-dim stencil corresponds to eq (5) in Stone et al. (2020).
+// Limiting based on implementation in AMReX (see https://github.com/AMReX-Codes/amrex/blob/735c3513153f1d06f783e64f455816be85fb3602/Src/AmrCore/AMReX_MFInterp_3D_C.H#L89)
+// to preserve extrema.
 struct ProlongateCellMinModMultiD {
   static constexpr bool OperationRequired(TopologicalElement fel,
                                           TopologicalElement cel) {
@@ -77,7 +81,7 @@ struct ProlongateCellMinModMultiD {
     Real dx1fm = 0;
     Real dx1fp = 0;
     Real gx1c = 0;
-    if constexpr (INCLUDE_X2) {
+    if constexpr (INCLUDE_X1) {
       Real dx1m, dx1p;
       GetGridSpacings<1, el>(coords, coarse_coords, cib, ib, i, fi, &dx1m, &dx1p, &dx1fm,
                              &dx1fp);

src/pgen/cluster/magnetic_tower.hpp

@@ -201,9 +201,6 @@ class MagneticTower {
       PARTHENON_REQUIRE_THROWS(offset_ <= 0.0 && thickness_ <= 0.0,
                                "Please disable (set to zero) tower offset and thickness "
                                "for the Li tower model");
-    } else {
-      PARTHENON_FAIL(
-          "Unknown potential for magnetic tower. Current options are: donut, li")
     }
 
     // Vector potential is only locally used, so no need to

src/pgen/cluster/stellar_feedback.cpp

@@ -140,9 +140,7 @@ void StellarFeedback::FeedbackSrcTerm(parthenon::MeshData<parthenon::Real> *md,
         AddDensityToConsAtFixedVelTemp(cell_delta_rho, cons, prim, eos.GetGamma(), k, j,
                                        i);
         //  Then add thermal energy
-        const auto cell_delta_mass = -cell_delta_rho * coords.CellVolume(k, j, i);
-        const auto cell_delta_energy_density =
-            mass_to_energy * cell_delta_mass / coords.CellVolume(k, j, i);
+        const auto cell_delta_energy_density = -mass_to_energy * cell_delta_rho;
         PARTHENON_REQUIRE(cell_delta_energy_density > 0.0,
                           "Sanity check failed. Added thermal energy should be positive.")
         cons(IEN, k, j, i) += cell_delta_energy_density;