Merge branch 'development' into add_cubic_interp

.github/workflows/dependencies_hip.sh

@@ -1,17 +1,26 @@
 #!/usr/bin/env bash
 #
-# Copyright 2020-2022 The AMReX Community
+# Copyright 2020 The AMReX Community
 #
 # License: BSD-3-Clause-LBNL
 # Authors: Axel Huebl
 
+# search recursive inside a folder if a file contains tabs
+#
+# @result 0 if no files are found, else 1
+#
+
 set -eu -o pipefail
 
+# `man apt.conf`:
+#   Number of retries to perform. If this is non-zero APT will retry
+#   failed files the given number of times.
+echo 'Acquire::Retries "3";' | sudo tee /etc/apt/apt.conf.d/80-retries
 
 # Ref.: https://rocmdocs.amd.com/en/latest/Installation_Guide/Installation-Guide.html#ubuntu
 curl -O https://repo.radeon.com/rocm/rocm.gpg.key
 sudo apt-key add rocm.gpg.key
-echo 'deb [arch=amd64] https://repo.radeon.com/rocm/apt/debian/ ubuntu main' \
+echo "deb [arch=amd64] https://repo.radeon.com/rocm/apt/${1-debian}/ ubuntu main" \
   | sudo tee /etc/apt/sources.list.d/rocm.list
 echo 'export PATH=/opt/rocm/llvm/bin:/opt/rocm/bin:/opt/rocm/profiler/bin:/opt/rocm/opencl/bin:$PATH' \
   | sudo tee -a /etc/profile.d/rocm.sh
@@ -34,12 +43,14 @@ sudo apt-get install -y --no-install-recommends \
     roctracer-dev   \
     rocprofiler-dev \
     rocrand-dev     \
-    rocprim-dev
+    rocprim-dev     \
+    hiprand-dev
 
 # activate
 #
 source /etc/profile.d/rocm.sh
 hipcc --version
+hipconfig --full
 which clang
 which clang++
 which flang

CHANGES.md

@@ -1,3 +1,7 @@
+# 23.12
+
+  * The radiation solver port to C++ has been completed (#2638, #2648)
+
 # 23.11
 
   * Problem GNUmakefiles have been standardized and now allow for the

Docs/source/io.rst

@@ -562,7 +562,15 @@ By default, 4 output files are created:
 
        The species masses are given in units of solar masses.
 
-Some problems have custom versions of the diagnostics with additional information.
+``Castro/Util/scripts/diag_parser.py`` contains Python code for parsing
+these output files into Numpy arrays.  Usage instructions are included
+in the file, along with an example script at
+``Castro/Util/scripts/plot_species.py``.  This reads a
+``species_diag.out`` file provided on the command line and makes a plot
+of the total mass fractions over time.
+
+Some problems have custom versions of the diagnostics with additional
+information.  These are not currently supported by the Python parser.
 
 
 .. _sec:parallel_io:

Exec/Make.Castro

@@ -270,6 +270,10 @@ ifeq ($(USE_ROTATION), TRUE)
   DEFINES += -DROTATION
 endif
 
+ifeq ($(USE_SPECIES_SOURCES), TRUE)
+  DEFINES += -DCONS_SPECIES_HAVE_SOURCES
+endif
+
 ifeq ($(USE_PARTICLES), TRUE)
   Bdirs += Source/particles
 endif

Exec/reacting_tests/nse_test/problem_initialize_state_data.H

@@ -5,6 +5,7 @@
 #include <eos.H>
 #include <network.H>
 #if defined(NSE_TABLE)
+#include <nse_table_type.H>
 #include <nse_table.H>
 #elif defined(NSE_NET)
 #include <nse_solver.H>
@@ -64,14 +65,15 @@ void problem_initialize_state_data (int i, int j, int k,
     burn_state.y_e = ye;
 
 #if defined(NSE_TABLE)
-    Real abar;
-    Real dq;
-    Real dyedt;
-    Real dabardt;
-    Real dbeadt;
-    Real e_nu;
-
-    nse_interp(T, problem::rho0, ye, abar, dq, dyedt, dabardt, dbeadt, e_nu, burn_state.xn);
+    nse_table_t nse_state;
+    nse_state.T = T;
+    nse_state.rho = problem::rho0;
+    nse_state.Ye = ye;
+    nse_interp(nse_state);
+
+    for (int n = 0; n < NumSpec; ++n) {
+        burn_state.xn[n] = nse_state.X[n];
+    }
 #elif defined(NSE_NET)
     Real eps = 1.e-10_rt;
     bool input_ye_is_valid = true;
@@ -105,16 +107,16 @@ void problem_initialize_state_data (int i, int j, int k,
     for (int n = 0; n < NumSpec; n++) {
         burn_state.xn[n] /= sumX;
     }
-    
+
 #ifdef NSE_NET
     state(i,j,k,UMUP) = burn_state.mu_p;
     state(i,j,k,UMUN) = burn_state.mu_n;
 #endif
-    
+
 #ifdef AUX_THERMO
     burn_state.aux[AuxZero::iye] = ye;
-    burn_state.aux[AuxZero::iabar] = abar;
-    burn_state.aux[AuxZero::ibea] = dq;
+    burn_state.aux[AuxZero::iabar] = nse_state.abar;
+    burn_state.aux[AuxZero::ibea] = nse_state.bea;
 #endif
 
     eos(eos_input_rt, burn_state);

Exec/science/massive_star/analysis/massive_star_multi.py

@@ -1,22 +1,24 @@
 #!/usr/bin/env python3
 
+import argparse
+import os
+
 import matplotlib
-matplotlib.use('agg')
+import matplotlib.pyplot as plt
 
-import os
-import sys
 import yt
-import matplotlib.pyplot as plt
-import numpy as np
-from functools import reduce
+from yt.frontends.boxlib.api import CastroDataset
+# assume that our data is in CGS
+from yt.units import cm
+
+matplotlib.use('agg')
 
-from mpl_toolkits.axes_grid1 import ImageGrid
 
-# assume that our data is in CGS
-from yt.units import cm, amu
-from yt.frontends.boxlib.api import CastroDataset
 
-def make_plot(plotfile, prefix="plot", size=(19.2, 10.8), cbar_location="right"):
+
+def make_plot(plotfile, prefix="plot", width_frac=0.1,
+              layout=(1, 4),
+              size=(19.2, 10.8), cbar_location="right"):
 
     ds = CastroDataset(plotfile)
 
@@ -36,8 +38,6 @@ def make_plot(plotfile, prefix="plot", size=(19.2, 10.8), cbar_location="right")
     fig = plt.figure()
     fig.set_size_inches(size)
 
-    width_frac = 0.1
-
     fields = ["MachNumber", "magvort", "abar", "enuc"]
 
     sp = yt.SlicePlot(ds, "theta", fields,
@@ -75,17 +75,39 @@ def make_plot(plotfile, prefix="plot", size=(19.2, 10.8), cbar_location="right")
 
     sp.set_axes_unit("cm")
 
-    fig = sp.export_to_mpl_figure((1, len(fields)), cbar_location=cbar_location, cbar_pad="5%")
+    fig = sp.export_to_mpl_figure((layout[0], layout[1]), axes_pad=(1.0, 0.4),
+                                  cbar_location=cbar_location, cbar_pad="2%")
 
     fig.subplots_adjust(left=0.05, right=0.95, bottom=0.025, top=0.975)
-    fig.text(0.02, 0.02, "time = {:8.5f} s".format(float(ds.current_time)), transform=fig.transFigure)
+    fig.text(0.02, 0.02, f"time = {float(ds.current_time):8.2f} s",
+             transform=fig.transFigure)
     fig.set_size_inches(size)
+    fig.tight_layout()
+    extra = ""
+    if layout[0] >= layout[1]:
+        extra = "_vertical"
 
-    fig.savefig(f"{prefix}_{os.path.basename(plotfile)}_slice.png", pad_inches=0)
+    fig.savefig(f"{prefix}_{os.path.basename(plotfile)}_w{width_frac:04.2f}{extra}.pdf", pad_inches=0)
 
 
 if __name__ == "__main__":
 
-    plotfile = sys.argv[1]
-
-    make_plot(plotfile, "all")
+    p = argparse.ArgumentParser()
+    p.add_argument("--vertical", action="store_true",
+                   help="plot 2x2 or 1x4")
+    p.add_argument("--width_fraction", type=float, default=0.1,
+                   help="fraction of domain to show")
+    p.add_argument("plotfile", type=str, nargs=1,
+                   help="plotfile to plot")
+
+    args = p.parse_args()
+    plotfile = args.plotfile[0]
+
+    if args.vertical:
+        size = (7.5, 11.0)
+        layout = (2, 2)
+    else:
+        size = (19.2, 8.5)
+        layout = (1, 4)
+
+    make_plot(plotfile, "all", layout=layout, width_frac=args.width_fraction, size=size)

Exec/science/massive_star/problem_initialize_state_data.H

@@ -5,6 +5,7 @@
 #include <eos.H>
 #include <model_parser.H>
 #ifdef NSE_TABLE
+#include <nse_table_type.H>
 #include <nse_table_check.H>
 #endif
 
@@ -80,21 +81,16 @@ void problem_initialize_state_data (int i, int j, int k,
 
     if (nse_check) {
 
-        Real abar;
-        Real dq;
-        Real dyedt;
-        Real dabardt;
-        Real dbeadt;
-        Real e_nu;
-        Real xn[NumSpec];
+        nse_table_t nse_state;
+        nse_state.T = state(i,j,k,UTEMP);
+        nse_state.rho = state(i,j,k,URHO);
+        nse_state.Ye = state(i,j,k,UFX+AuxZero::iye);
+        nse_interp(nse_state);
 
-        nse_interp(state(i,j,k,UTEMP), state(i,j,k,URHO),
-                   state(i,j,k,UFX+AuxZero::iye), abar, dq, dyedt, dabardt, dbeadt, e_nu, xn);
-
-        state(i,j,k,UFX+AuxZero::iabar) = abar;
-        state(i,j,k,UFX+AuxZero::ibea) = dq;
+        state(i,j,k,UFX+AuxZero::iabar) = nse_state.abar;
+        state(i,j,k,UFX+AuxZero::ibea) = nse_state.bea;
         for (int n = 0; n < NumSpec; n++) {
-            state(i,j,k,UFS+n) = amrex::max(xn[n], small_x);
+            state(i,j,k,UFS+n) = amrex::max(nse_state.X[n], small_x);
         }
 
         // we just got new X's, so we need to re-renormalize them

Source/driver/Castro.cpp

@@ -1260,8 +1260,7 @@ Castro::initData ()
              {
                const Box& box = mfi.validbox();
 
-               tmp.resize(box, 1);
-               Elixir elix_tmp = tmp.elixir();
+               tmp.resize(box, 1, The_Async_Arena());
                auto tmp_arr = tmp.array();
 
                make_fourth_in_place(box, Sborder.array(mfi), tmp_arr, domain_lo, domain_hi);
@@ -1287,8 +1286,7 @@ Castro::initData ()
            {
              const Box& box = mfi.growntilebox(2);
 
-             tmp.resize(box, 1);
-             Elixir elix_tmp = tmp.elixir();
+             tmp.resize(box, 1, The_Async_Arena());
              auto tmp_arr = tmp.array();
 
              make_cell_center_in_place(box, Sborder.array(mfi), tmp_arr, domain_lo, domain_hi);
@@ -1336,8 +1334,7 @@ Castro::initData ()
            {
              const Box& box = mfi.validbox();
 
-             tmp.resize(box, 1);
-             Elixir elix_tmp = tmp.elixir();
+             tmp.resize(box, 1, The_Async_Arena());
              auto tmp_arr = tmp.array();
 
              make_fourth_in_place(box, Sborder.array(mfi), tmp_arr, domain_lo, domain_hi);
@@ -3950,8 +3947,7 @@ Castro::computeTemp(
       compute_lap_term(bx0, Stemp.array(mfi), Eint_lap.array(mfi), UEINT,
                        domain_lo, domain_hi);
 
-      tmp.resize(bx, 1);
-      Elixir elix_tmp = tmp.elixir();
+      tmp.resize(bx, 1, The_Async_Arena());
       auto tmp_arr = tmp.array();
 
       make_cell_center_in_place(bx, Stemp.array(mfi), tmp_arr, domain_lo, domain_hi);
@@ -4067,8 +4063,7 @@ Castro::computeTemp(
 
       const Box& bx = mfi.tilebox();
 
-      tmp.resize(bx, 1);
-      Elixir elix_tmp = tmp.elixir();
+      tmp.resize(bx, 1, The_Async_Arena());
       auto tmp_arr = tmp.array();
 
       // only temperature

Source/driver/Castro_advance_sdc.cpp

@@ -122,8 +122,7 @@ Castro::do_advance_sdc (Real time,
           for (MFIter mfi(S_new); mfi.isValid(); ++mfi) {
             const Box& bx = mfi.tilebox();
 
-            tmp.resize(bx, 1);
-            Elixir elix_tmp = tmp.elixir();
+            tmp.resize(bx, 1, The_Async_Arena());
             auto tmp_arr = tmp.array();
 
             make_fourth_in_place(bx, old_source.array(mfi), tmp_arr, domain_lo, domain_hi);
@@ -305,8 +304,7 @@ Castro::do_advance_sdc (Real time,
       // pass in the reaction source at centers (Sburn_arr), including
       // one ghost cell and derive everything that is needed including
       // 1 ghost cell
-      R_center.resize(obx, R_new.nComp());
-      Elixir elix_r_center = R_center.elixir();
+      R_center.resize(obx, R_new.nComp(), The_Async_Arena());
       auto const R_center_arr = R_center.array();
 
       Array4<const Real> const Sburn_arr = Sburn.array(mfi);
@@ -315,8 +313,7 @@ Castro::do_advance_sdc (Real time,
       ca_store_reaction_state(obx, Sburn_arr, R_center_arr);
 
       // convert R_new from centers to averages in place
-      tmp.resize(bx, 1);
-      Elixir elix_tmp = tmp.elixir();
+      tmp.resize(bx, 1, The_Async_Arena());
       auto const tmp_arr = tmp.array();
 
       make_fourth_in_place(bx, R_center_arr, tmp_arr, domain_lo, domain_hi);

Source/driver/Derive.cpp

@@ -352,8 +352,7 @@ extern "C"
       const Box& obx = amrex::grow(bx, 1);
 
       FArrayBox coeff_cc;
-      coeff_cc.resize(obx, 1);
-      Elixir elix_coeff_cc = coeff_cc.elixir();
+      coeff_cc.resize(obx, 1, The_Async_Arena());
       Array4<Real> const coeff_arr = coeff_cc.array();
 
       auto const dat = datfab.array();

Source/driver/_cpp_parameters

@@ -92,9 +92,11 @@ hybrid_hydro                 int           0
 # reconstruction type:
 # 0: piecewise linear;
 # 1: classic Colella \& Woodward ppm;
-# 2: extrema-preserving ppm
 ppm_type                     int           1
 
+# do we limit the ppm parabola?
+ppm_do_limiting              int           1
+
 # For MHD + PLM, do we limit on characteristic or primitive variables
 mhd_limit_characteristic     int           1
 

Source/driver/_variables

@@ -16,9 +16,9 @@
    energy-density         UEDEN     NSRC                                  1                 None
    internal-energy        UEINT     NSRC                                  1                 None
    temperature            UTEMP     NSRC                                  1                 None
-   advected               UFA       None                                  NumAdv            None
-   species                UFS       None                                  NumSpec           None
-   auxiliary              UFX       None                                  NumAux            None
+   advected               UFA       [(NSRC, CONS_SPECIES_HAVE_SOURCES)]   NumAdv            None
+   species                UFS       [(NSRC, CONS_SPECIES_HAVE_SOURCES)]   NumSpec           None
+   auxiliary              UFX       [(NSRC, CONS_SPECIES_HAVE_SOURCES)]   NumAux            None
    shock                  USHK      None                                  1                 SHOCK_VAR
    mu_p                   UMUP      None                                  1                 NSE_NET
    mu_n                   UMUN      None                                  1                 NSE_NET