Merge branch 'development' into new_shock_flag_call

.github/workflows/docs-test.yml

@@ -9,6 +9,10 @@ on:
     branches:
       - development
 
+env:
+  # enable color output from Sphinx
+  FORCE_COLOR: "1"
+
 jobs:
   docs:
     runs-on: ubuntu-latest
@@ -33,7 +37,10 @@ jobs:
       - name: Build docs
         run: |
           cd Docs/
-          make SPHINXOPTS=-v NO_DOXYGEN=TRUE html
+          # remove missing TOC entries that would otherwise be generated by
+          # doxygen to avoid warnings about missing references
+          sed -i -e 's/^   filelist$//; s/^   classlist$//' source/index.rst
+          make SPHINXOPTS='-v -W --keep-going' NO_DOXYGEN=TRUE html
 
       - name: Check links
         run: |

.gitignore

@@ -72,6 +72,8 @@ Docs/source/namespacelist.rst
 Docs/source/runtime_parameters.rst
 Docs/source/*_files.rst
 Docs/source/preprocessed_files
+Docs/source/yt_example.rst
+Docs/source/yt_example_files/
 
 
 amr_diag.out

Docs/rp.py

@@ -7,17 +7,17 @@
 from more_itertools import unique_everseen
 
 MAIN_HEADER = """
-+--------------------------------------------+-------------------------------------------------------------+---------------+
-| parameter                                  | description                                                 | default value |
-+============================================+=============================================================+===============+
++--------------------------------------------+-------------------------------------------------------------+-----------------------------+
+| parameter                                  | description                                                 | default value               |
++============================================+=============================================================+=============================+
 """
 
 SEPARATOR = """
-+--------------------------------------------+-------------------------------------------------------------+---------------+
++--------------------------------------------+-------------------------------------------------------------+-----------------------------+
 """
 
 ENTRY = """
-| {:42} | {:59} | {:13} |
+| {:42} | {:59} | {:27} |
 """
 
 WRAP_LEN = 59

Docs/source/FlowChart.rst

@@ -309,7 +309,7 @@ In the code, the objective is to evolve the state from the old time,
    A. Create ``Sborder``, initialized from ``S_old``
 
    B. Call ``clean_state()`` to make sure the thermodynamics are in
-     sync, in particular, compute the temperature.
+      sync, in particular, compute the temperature.
 
    C. [``SHOCK_VAR``] zero out the shock flag.
 

Docs/source/Hydrodynamics.rst

@@ -342,6 +342,14 @@ accounted for in **Steps 1** and **6**. The source terms are:
    S_{{\rm ext},\rho Y_k}
    \end{array}\right)^n.
 
+.. index:: USE_SPECIES_SOURCES
+
+.. note:: To reduce memory usage, we do not include source terms for the
+   advected quantities, species, and auxiliary variables in the conserved
+   state vector by default. If your application needs external source terms for
+   these variables, set ``USE_SPECIES_SOURCES=TRUE`` when compiling so that space
+   will be allocated for them.
+
 Primitive Forms
 ===============
 
@@ -585,9 +593,6 @@ runtime parameters for hydrodynamics:
 
    See :ref:`sponge_section` for more details on the sponge.
 
--  ``castro.normalize_species``: enforce that :math:`\sum_i X_i = 1`
-   (0 or 1; default: 0)
-
 .. index:: castro.small_dens, castro.small_temp, castro.small_pres
 
 Several floors are imposed on the thermodynamic quantities to prevet unphysical

Docs/source/docutils.conf

@@ -0,0 +1,2 @@
+[parsers]
+line_length_limit = 1000000

Docs/source/index.rst

@@ -60,7 +60,7 @@ https://github.com/amrex-astro/Castro
 
    filelist
    classlist
-   .. namespacelist
+.. namespacelist
 
 .. toctree::
    :caption: References

Docs/source/io.rst

@@ -271,8 +271,6 @@ radiation quantities).
 | (where X is any of the species    | :math:`\omegadot_k = DX_k/Dt`                     |                                      |
 | defined in the network)           |                                                   |                                      |
 +-----------------------------------+---------------------------------------------------+--------------------------------------+
-| ``enuc``                          | Nuclear energy generation rate / gram             | :math:`{\rm erg~g^{-1}~s^{-1}}`      |
-+-----------------------------------+---------------------------------------------------+--------------------------------------+
 | ``rho_enuc``                      | Nuclear energy generation rate density            | :math:`{\rm erg~cm^{-3}~s^{-1}}`     |
 +-----------------------------------+---------------------------------------------------+--------------------------------------+
 | ``phiGrav``                       | Gravitational potential                           | :math:`{\rm erg~g^{-1}}`             |
@@ -393,7 +391,8 @@ Derived variables
 | ``y_velocity``,                   | :math:`\ub = (\rho \ub)/\rho`                     |                             |                                         |
 | ``z_velocity``                    |                                                   |                             |                                         |
 +-----------------------------------+---------------------------------------------------+-----------------------------+-----------------------------------------+
-
+| ``enuc``                          | Nuclear energy generation rate / gram             | ``derenuc``                 | :math:`{\rm erg~g^{-1}~s^{-1}}`         |
++-----------------------------------+---------------------------------------------------+-----------------------------+-----------------------------------------+
 
 problem-specific plotfile variables
 ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Exec/hydro_tests/oddeven/inputs.2d

@@ -22,6 +22,7 @@ castro.do_hydro = 1
 castro.do_react = 0
 castro.ppm_type = 0
 castro.hybrid_riemann = 1
+castro.shock_detection_threshold = 0.5
 
 # TIME STEP CONTROL
 

Exec/reacting_tests/nse_test/GNUmakefile

@@ -11,6 +11,7 @@ USE_MPI    = TRUE
 
 USE_REACT  = TRUE
 
+USE_SIMPLIFIED_SDC := TRUE
 
 CASTRO_HOME ?= ../../..
 

Exec/reacting_tests/nse_test/README.md

@@ -1,4 +1,4 @@
-# nse_test
+# `nse_test`
 
 This is a simple test problem designed to explore how well hydro and
 reactions are coupled when a system enters NSE.
@@ -7,7 +7,10 @@ This version is based on ``reacting_convergence`` (which is in turn
 based on ``acoustic_pulse_general``), but using the ``aprox19``
 network with the NSE table enabled.
 
-You can run the Strang convergence test with the script convergence_strang.sh
+You can run the simplified-SDC convergence test with the script
+`convergence_simplified_sdc_w_vel.sh`
+
+The script `create_pretty_tables.py` will take the 2 (or 3) output
+files and make a single LaTeX-formatted table of the results.  Use
+the `--simple` argument to format in plaintext.
 
-The script create_pretty_tables.py will take the 2 output files and
-make a single LaTeX-formatted table of the results.

Exec/reacting_tests/nse_test/convergence_simplified_sdc_w_vel.sh

@@ -10,17 +10,17 @@ problem.u0=1.e8
 problem.v0=1.e8
 "
 
-mpiexec -n 8 ${EXEC} inputs.32 ${RUNPARAMS} >& /dev/null
-mpiexec -n 8 ${EXEC} inputs.64 ${RUNPARAMS} >& /dev/null
-mpiexec -n 8 ${EXEC} inputs.128 ${RUNPARAMS} >& /dev/null
+mpiexec -n 8 ${EXEC} inputs.32 ${RUNPARAMS} >& 32.out
+mpiexec -n 8 ${EXEC} inputs.64 ${RUNPARAMS} >& 64.out
+mpiexec -n 8 ${EXEC} inputs.128 ${RUNPARAMS} >& 128.out
 
 RichardsonConvergenceTest2d.gnu.ex coarFile=nse_test_32_plt00080 mediFile=nse_test_64_plt00160 fineFile=nse_test_128_plt00320  >& nse_convergence_simple_sdc_vlo.out
 
-mpiexec -n 8 ${EXEC} inputs.256 ${RUNPARAMS} >& /dev/null
+mpiexec -n 8 ${EXEC} inputs.256 ${RUNPARAMS} >& 256.out
 
 RichardsonConvergenceTest2d.gnu.ex coarFile=nse_test_64_plt00160 mediFile=nse_test_128_plt00320 fineFile=nse_test_256_plt00640  >& nse_convergence_simple_sdc_lo.out
 
-mpiexec -n 8 ${EXEC} inputs.512 ${RUNPARAMS} >& /dev/null
+mpiexec -n 8 ${EXEC} inputs.512 ${RUNPARAMS} >& 512.out
 
 RichardsonConvergenceTest2d.gnu.ex coarFile=nse_test_128_plt00320 mediFile=nse_test_256_plt00640 fineFile=nse_test_512_plt01280 >& nse_convergence_simple_sdc_hi.out
 

Exec/reacting_tests/nse_test/problem_initialize_state_data.H

@@ -16,7 +16,6 @@ void problem_initialize_state_data (int i, int j, int k,
                                     Array4<Real> const& state,
                                     const GeometryData& geomdata)
 {
-    int coord_type = geomdata.Coord();
 
     const Real* dx = geomdata.CellSize();
     const Real* problo = geomdata.ProbLo();
@@ -101,11 +100,11 @@ void problem_initialize_state_data (int i, int j, int k,
     // since the species are interpolated, normalize them
 
     Real sumX = 0.0_rt;
-    for (int n = 0; n < NumSpec; n++) {
-        sumX += burn_state.xn[n];
+    for (auto X : burn_state.xn) {
+        sumX += X;
     }
-    for (int n = 0; n < NumSpec; n++) {
-        burn_state.xn[n] /= sumX;
+    for (auto & X : burn_state.xn) {
+        X /= sumX;
     }
 
 #ifdef NSE_NET

Exec/science/massive_star/analysis/initial_model.py

@@ -10,9 +10,8 @@ def find_r_for_rho(r, rho, rho_want):
     return r[idx]
 
 
-#file = "../15m_500_sec.aprox19.hse.6400"
-file = "../15m_500_sec.aprox19.hse.20.0km"
-Lx = 1.6384e10
+file = "../15m_500_sec.aprox19.hse.5.00km"
+Lx = 8.192e9
 
 data = np.loadtxt(file)
 
@@ -65,6 +64,10 @@ def find_r_for_rho(r, rho, rho_want):
 
 ax.set_xlabel("r [cm]")
 ax.set_ylabel(r"$\rho~[\rm{g/cm^3}]$, $T~[K]$")
+
+ax.set_xlim(None, 5.e10)
+ax.grid(color="b", alpha=0.5, ls=":")
+
 ax2 = ax.twinx()
 
 ax2.set_ylabel(r"$Y_e$")
@@ -74,8 +77,9 @@ def find_r_for_rho(r, rho, rho_want):
 lns = l1 + l2 + l3
 labs = [l.get_label() for l in lns]
 
-ax.legend(lns, labs, frameon=False)
+ax.legend(lns, labs, frameon=False, loc=6)
 
+# species plot
 
 ax = fig.add_subplot(212)
 
@@ -92,27 +96,34 @@ def find_r_for_rho(r, rho, rho_want):
 
     if Xmax > threshold:
         if Xmax > 0.5:
-            lw = 2
+            lw = 3
+            ls = "-"
         else:
             lw = 1
-        ax.plot(data[:,0], data[:,n], label=var, lw=lw)
+            ls = "--"
+        ax.plot(data[:,0], data[:,n], label=var, lw=lw, ls=ls)
 
 
 ax.axvline(r_refine, color="0.25", ls=":")
 
-ax.axvline(Lx, color="0.25", ls="-")
+ax.axvline(Lx, color="0.5", ls="-")
 
 ax.set_xscale("log")
 ax.set_yscale("log")
 
 ax.set_xlabel("r [cm]")
 ax.set_ylabel("mass fraction")
 
-ax.legend(frameon=False, fontsize="small", ncol=2)
+ax.set_xlim(None, 5.e10)
+ax.set_ylim(1.e-6, None)
+
+ax.grid(color="b", alpha=0.5, ls=":")
+
+ax.legend(frameon=True, edgecolor="w", ncol=1, framealpha=0.5)
 
 fig.set_size_inches((8, 12))
 
 fig.tight_layout()
 
-fig.savefig("initial_model.png")
+fig.savefig("initial_model.pdf")
 

Exec/science/massive_star/analysis/massive_star_multi.py

@@ -22,6 +22,13 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1,
 
     ds = CastroDataset(plotfile)
 
+    t_drive = 0.0
+    if "[*] castro.drive_initial_convection_tmax" in ds.parameters:
+        t_drive = ds.parameters["[*] castro.drive_initial_convection_tmax"]
+    elif "castro.drive_initial_convection_tmax" in ds.parameters:
+        t_drive = ds.parameters["castro.drive_initial_convection_tmax"]
+    print(t_drive)
+
     xmin = ds.domain_left_edge[0]
     xmax = ds.domain_right_edge[0]
 
@@ -70,32 +77,33 @@ def make_plot(plotfile, prefix="plot", width_frac=0.1,
 
         if f == "enuc":
             # now do a contour of NSE
-            sp.annotate_contour("in_nse", ncont=1, clim=(0.5, 0.5), take_log=False,
-                                plot_args={"colors": "k", "linewidths": 2})
+            if ("boxlib", "in_nse") in ds.derived_field_list:
+                sp.annotate_contour("in_nse", ncont=1, clim=(0.5, 0.5), take_log=False,
+                                    plot_args={"colors": "k", "linewidths": 2})
 
     sp.set_axes_unit("cm")
 
     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, f"time = {float(ds.current_time):8.2f} s",
+    fig.text(0.02, 0.02, f"$t - \\tau_\\mathrm{{drive}}$ = {float(ds.current_time) - t_drive:6.1f} 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)}_w{width_frac:04.2f}{extra}.pdf", pad_inches=0)
+    fig.savefig(f"{prefix}_{os.path.basename(plotfile)}_w{width_frac:04.2f}{extra}.pdf", pad_inches=0.1, bbox_inches="tight")
 
 
 if __name__ == "__main__":
 
     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,
+    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")