diff --git a/Docs/source/FlowChart.rst b/Docs/source/FlowChart.rst index 419c58b877..e9c69ebf70 100644 --- a/Docs/source/FlowChart.rst +++ b/Docs/source/FlowChart.rst @@ -56,8 +56,7 @@ The time-integration method used is controlled by order integration implemented. At the moment, this does not support multilevel domains. Note: because of differences in the interfaces with the default Strang method, you must compile with ``USE_TRUE_SDC = TRUE`` for this - method to work (in particular, this defines ``EXTRA_THERMO`` which enables some - additional EOS derivatives). + method to work. * ``time_integration_method = 3``: this is the simplified SDC method described above that uses the CTU hydro advection and an ODE diff --git a/Exec/Make.Castro b/Exec/Make.Castro index 3b0b5cb2ec..4a645154ef 100644 --- a/Exec/Make.Castro +++ b/Exec/Make.Castro @@ -186,18 +186,6 @@ endif ifeq ($(USE_REACT), TRUE) Bdirs += Source/reactions DEFINES += -DREACTIONS - - ifeq ($(USE_TRUE_SDC), TRUE) - # we need the compositional derivatives for SDC - DEFINES += -DEXTRA_THERMO - endif - - ifeq ($(USE_SIMPLIFIED_SDC), TRUE) - # we need the compositional derivatives for SDC - DEFINES += -DEXTRA_THERMO - # we only implement this for C++ reactions - endif - endif ifeq ($(USE_REACT_SPARSE_JACOBIAN), TRUE) diff --git a/Exec/science/flame_wave/analysis/vis_3d/vol-xrb-abar.py b/Exec/science/flame_wave/analysis/vis_3d/vol-xrb-abar.py index 078cf554c3..c5faa43ac3 100644 --- a/Exec/science/flame_wave/analysis/vis_3d/vol-xrb-abar.py +++ b/Exec/science/flame_wave/analysis/vis_3d/vol-xrb-abar.py @@ -1,16 +1,19 @@ #!/usr/bin/env python -import matplotlib -matplotlib.use('agg') - import sys +import matplotlib +import numpy as np + import yt from yt.frontends.boxlib.api import CastroDataset -import numpy as np -#from yt.visualization.volume_rendering.render_source import VolumeSource -from yt.visualization.volume_rendering.api import create_volume_source, Scene from yt.units import cm +#from yt.visualization.volume_rendering.render_source import VolumeSource +from yt.visualization.volume_rendering.api import Scene, create_volume_source + +matplotlib.use('agg') + +resolution = (1920, 1080) # this is for the wdconvect problem @@ -56,7 +59,7 @@ def doit(plotfile): sc.get_source(0).transfer_function = tf cam = sc.add_camera(ds, lens_type="perspective") - cam.resolution = (1920, 1280) + cam.resolution = resolution # view 1 @@ -77,14 +80,24 @@ def doit(plotfile): cam.zoom(3.0) sc.camera = cam - sc.save_annotated("{}_abar_annotated_side.png".format(plotfile), - sigma_clip=3.0, + sc.save(f"{plotfile}_abar_noaxes_side.png", sigma_clip=3.0) + + sc.annotate_axes(alpha=0.005, thickness=6) + + sc.save(f"{plotfile}_abar_side.png", sigma_clip=3.0) + + sc.save_annotated(f"{plotfile}_abar_annotated_side.png", + sigma_clip=3.0, label_fmt="%.2f", text_annotate=[[(0.05, 0.05), f"t = {ds.current_time.d:7.5f} s", dict(horizontalalignment="left")]]) # view 2 + # remove the annotation source for now + print(list(sc.sources.keys())) + sc.sources.pop("source_01") + dx = ds.domain_right_edge[0] - ds.domain_left_edge[0] cam.position = [0.5*(ds.domain_left_edge[0] + ds.domain_right_edge[0]) + 0.0001 * dx, 0.5*(ds.domain_left_edge[1] + ds.domain_right_edge[1]), @@ -103,8 +116,13 @@ def doit(plotfile): cam.zoom(0.6) sc.camera = cam + sc.save(f"{plotfile}_abar_noaxes_top.png", sigma_clip=3.0) + + sc.annotate_axes(alpha=0.005, thickness=6) + + sc.save(f"{plotfile}_abar_top.png", sigma_clip=3.0) sc.save_annotated("{}_abar_annotated_top.png".format(plotfile), - sigma_clip=3.0, + sigma_clip=3.0, label_fmt="%.2f", text_annotate=[[(0.05, 0.05), f"t = {ds.current_time.d:7.5f} s", dict(horizontalalignment="left")]]) diff --git a/Exec/science/subchandra/GNUmakefile b/Exec/science/subchandra/GNUmakefile index ed57464fb5..bc66f50e7a 100644 --- a/Exec/science/subchandra/GNUmakefile +++ b/Exec/science/subchandra/GNUmakefile @@ -22,7 +22,7 @@ USE_SIMPLIFIED_SDC=TRUE EOS_DIR := helmholtz # This sets the network directory in $(MICROPHYSICS_HOME)/networks -NETWORK_DIR := subch_approx +NETWORK_DIR := subch_simple INTEGRATOR_DIR := VODE diff --git a/Exec/science/subchandra/inputs_2d.N14.coarse b/Exec/science/subchandra/inputs_2d.N14.coarse index 829526dac3..77781e0e3f 100644 --- a/Exec/science/subchandra/inputs_2d.N14.coarse +++ b/Exec/science/subchandra/inputs_2d.N14.coarse @@ -62,7 +62,7 @@ castro.sponge_timescale = 1.e-3 castro.cfl = 0.2 # cfl number for hyperbolic system castro.init_shrink = 0.05 # scale back initial timestep by this factor castro.change_max = 1.025 # factor by which dt is allowed to change each timestep -castro.sum_interval = 0 # timesteps between computing and printing volume averages +castro.sum_interval = 5 # timesteps between computing and printing volume averages castro.update_sources_after_reflux = 0 castro.time_integration_method = 3 @@ -147,12 +147,12 @@ integrator.rtol_spec = 1.e-5 integrator.atol_spec = 1.e-5 integrator.rtol_enuc = 1.e-5 integrator.atol_enuc = 1.e-5 -integrator.jacobian = 2 +integrator.jacobian = 1 integrator.X_reject_buffer = 4.0 # disable jacobian caching in VODE integrator.use_jacobian_caching = 0 -integrator.ode_max_steps = 500000 +integrator.ode_max_steps = 1000000 diff --git a/Source/driver/Castro.cpp b/Source/driver/Castro.cpp index 4401575242..807dada23c 100644 --- a/Source/driver/Castro.cpp +++ b/Source/driver/Castro.cpp @@ -39,7 +39,7 @@ #include #endif -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #include #endif @@ -1047,8 +1047,10 @@ Castro::initData () By_new.setVal(0.0); Bz_new.setVal(0.0); - for (MFIter mfi(S_new); mfi.isValid(); ++mfi) { - +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { auto geomdata = geom.data(); @@ -1093,9 +1095,12 @@ Castro::initData () #endif //MHD - for (MFIter mfi(S_new); mfi.isValid(); ++mfi) +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& box = mfi.validbox(); + const Box& box = mfi.tilebox(); auto s = S_new[mfi].array(); auto geomdata = geom.data(); @@ -1131,7 +1136,7 @@ Castro::initData () ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -1185,8 +1190,11 @@ Castro::initData () // Verify that the sum of (rho X)_i = rho at every cell - for (MFIter mfi(S_new); mfi.isValid(); ++mfi) { - const Box& bx = mfi.validbox(); +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { + const Box& bx = mfi.tilebox(); auto S_arr = S_new.array(mfi); @@ -1353,17 +1361,11 @@ Castro::initData () #ifdef RADIATION if (do_radiation) { - for (MFIter mfi(S_new); mfi.isValid(); ++mfi) { - int i = mfi.index(); - - if (radiation->verbose > 2) { - std::cout << "Calling RADINIT at level " << level << ", grid " - << i << std::endl; - } - - const Box& box = mfi.validbox(); - const int* lo = box.loVect(); - const int* hi = box.hiVect(); +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { + const Box& box = mfi.tilebox(); auto r = Rad_new[mfi].array(); auto geomdata = geom.data(); @@ -1388,9 +1390,7 @@ Castro::initData () // by a problem setup (defaults to an empty routine). problem_initialize_rad_data(i, j, k, r, xnu_pass, nugroup_pass, dnugroup_pass, geomdata); - }); - } } #endif // RADIATION @@ -2539,7 +2539,7 @@ Castro::advance_aux(Real time, Real dt) MultiFab& S_old = get_old_data(State_Type); MultiFab& S_new = get_new_data(State_Type); -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(S_old, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -2716,6 +2716,9 @@ Castro::reflux (int crse_level, int fine_level, bool in_post_timestep) // Now zero out any problematic flux corrections. +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif for (MFIter mfi(crse_state, TilingIfNotGPU()); mfi.isValid(); ++mfi) { const Box& bx = mfi.tilebox(); @@ -3078,7 +3081,7 @@ Castro::normalize_species (MultiFab& S_new, int ng) ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -3152,12 +3155,12 @@ Castro::enforce_consistent_e ( BL_PROFILE("Castro::enforce_consistent_e()"); -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(S, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& box = mfi.tilebox(); + const Box& box = mfi.tilebox(); auto S_arr = S.array(mfi); @@ -3216,15 +3219,13 @@ Castro::enforce_min_density (MultiFab& state_in, int ng) } -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(state_in, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& bx = mfi.growntilebox(ng); do_enforce_minimum_density(bx, state_in.array(mfi), verbose); - } if (print_update_diagnostics) @@ -3250,7 +3251,7 @@ Castro::check_for_negative_density () ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { @@ -3337,11 +3338,10 @@ Castro::enforce_speed_limit (MultiFab& state_in, int ng) return; } -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(state_in, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& bx = mfi.growntilebox(ng); auto u = state_in[mfi].array(); @@ -3467,13 +3467,13 @@ Castro::apply_problem_tags (TagBoxArray& tags, Real time) int lev = level; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif { - for (MFIter mfi(tags); mfi.isValid(); ++mfi) + for (MFIter mfi(tags, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& bx = mfi.validbox(); + const Box& bx = mfi.tilebox(); TagBox& tagfab = tags[mfi]; @@ -3532,10 +3532,10 @@ Castro::apply_tagging_restrictions(TagBoxArray& tags, [[maybe_unused]] Real time blocking_factor[dim] = parent->blockingFactor(lev)[dim]; } -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif - for (MFIter mfi(tags); mfi.isValid(); ++mfi) { + for (MFIter mfi(tags, TilingIfNotGPU()); mfi.isValid(); ++mfi) { const Box& bx = mfi.tilebox(); auto tag = tags[mfi].array(); @@ -3576,10 +3576,10 @@ Castro::apply_tagging_restrictions(TagBoxArray& tags, [[maybe_unused]] Real time auto geomdata = geom.data(); // Allow the user to limit tagging outside of some distance from the problem center. -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif - for (MFIter mfi(tags); mfi.isValid(); ++mfi) { + for (MFIter mfi(tags, TilingIfNotGPU()); mfi.isValid(); ++mfi) { const Box& bx = mfi.tilebox(); auto tag = tags[mfi].array(); @@ -3755,7 +3755,7 @@ Castro::reset_internal_energy( } // Ensure (rho e) isn't too small or negative -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(State, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -3792,12 +3792,12 @@ Castro::add_magnetic_e( MultiFab& Bx, MultiFab& State) { -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(State, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& box = mfi.tilebox(); + const Box& box = mfi.tilebox(); auto S_arr = State.array(mfi); auto Bx_arr = Bx.array(mfi); auto By_arr = By.array(mfi); @@ -3837,12 +3837,12 @@ Castro::check_div_B( MultiFab& Bx, using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(State, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& box = mfi.tilebox(); + const Box& box = mfi.tilebox(); auto Bx_arr = Bx.array(mfi); auto By_arr = By.array(mfi); auto Bz_arr = Bz.array(mfi); @@ -3973,12 +3973,11 @@ Castro::computeTemp( } #endif -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(State, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - int num_ghost = ng; #ifdef TRUE_SDC if (sdc_order == 4) { @@ -4038,7 +4037,6 @@ Castro::computeTemp( } }); } - } #ifdef TRUE_SDC @@ -4358,7 +4356,7 @@ Castro::getCPUTime() { int numCores = ParallelDescriptor::NProcs(); -#ifdef _OPENMP +#ifdef AMREX_USE_OMP numCores = numCores*omp_get_max_threads(); #endif diff --git a/Source/driver/Castro_advance.cpp b/Source/driver/Castro_advance.cpp index 210671e17c..00f078226a 100644 --- a/Source/driver/Castro_advance.cpp +++ b/Source/driver/Castro_advance.cpp @@ -432,9 +432,12 @@ Castro::initialize_advance(Real time, Real dt, int amr_iteration) amrex::Print() << "<<<<< drive initial convection reset >>>>" << std::endl; - for (MFIter mfi(S_old); mfi.isValid(); ++mfi) +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_old, TilingIfNotGPU()); mfi.isValid(); ++mfi) { - const Box& box = mfi.validbox(); + const Box& box = mfi.tilebox(); auto s = S_old[mfi].array(); auto geomdata = geom.data(); diff --git a/Source/driver/Castro_io.cpp b/Source/driver/Castro_io.cpp index f6cad5cb5c..ebbb328377 100644 --- a/Source/driver/Castro_io.cpp +++ b/Source/driver/Castro_io.cpp @@ -26,7 +26,7 @@ #include #endif -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #include #endif @@ -280,13 +280,14 @@ Castro::restart (Amr& papa, amrex::Abort(); } - for (MFIter mfi(S_new); mfi.isValid(); ++mfi) +#ifdef AMREX_USE_OMP +#pragma omp parallel +#endif + for (MFIter mfi(S_new, TilingIfNotGPU()); mfi.isValid(); ++mfi) { + const Box& bx = mfi.tilebox(); - const Box& bx = mfi.validbox(); - - if (! orig_domain.contains(bx)) { - + if (!orig_domain.contains(bx)) { auto s = S_new[mfi].array(); auto geomdata = geom.data(); @@ -297,7 +298,6 @@ Castro::restart (Amr& papa, // by a problem setup (defaults to an empty routine). problem_initialize_state_data(i, j, k, s, geomdata); }); - } } } @@ -554,7 +554,7 @@ Castro::writeJobInfo (const std::string& dir, const Real io_time) jobInfoFile << "inputs file: " << inputs_name << "\n\n"; jobInfoFile << "number of MPI processes: " << ParallelDescriptor::NProcs() << "\n"; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP jobInfoFile << "number of threads: " << omp_get_max_threads() << "\n"; #endif jobInfoFile << "\n"; diff --git a/Source/driver/sum_utils.cpp b/Source/driver/sum_utils.cpp index 1b8c3635f8..d84e03eb97 100644 --- a/Source/driver/sum_utils.cpp +++ b/Source/driver/sum_utils.cpp @@ -38,7 +38,7 @@ Castro::volWgtSum (const MultiFab& mf, int comp, bool local, bool finemask) ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(mf, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -101,7 +101,7 @@ Castro::locWgtSum (const MultiFab& mf, int comp, int idir, bool local) auto dx = geom.CellSizeArray(); auto problo = geom.ProbLoArray(); -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(mf, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -195,7 +195,7 @@ Castro::volProductSum (const MultiFab& mf1, ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(mf1, TilingIfNotGPU()); mfi.isValid(); ++mfi) @@ -250,7 +250,7 @@ Castro::locSquaredSum (const std::string& name, auto dx = geom.CellSizeArray(); auto problo = geom.ProbLoArray(); -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(*mf, TilingIfNotGPU()); mfi.isValid(); ++mfi) diff --git a/Source/driver/timestep.cpp b/Source/driver/timestep.cpp index 10f0655a1b..3f81ae798a 100644 --- a/Source/driver/timestep.cpp +++ b/Source/driver/timestep.cpp @@ -474,7 +474,7 @@ Castro::estdt_rad (int is_new) ReduceData reduce_data(reduce_op); using ReduceTuple = typename decltype(reduce_data)::Type; -#ifdef _OPENMP +#ifdef AMREX_USE_OMP #pragma omp parallel #endif for (MFIter mfi(stateMF, TilingIfNotGPU()); mfi.isValid(); ++mfi) diff --git a/Source/problems/hse_fill.cpp b/Source/problems/hse_fill.cpp index 9354b9ef0b..cbff091621 100644 --- a/Source/problems/hse_fill.cpp +++ b/Source/problems/hse_fill.cpp @@ -112,7 +112,7 @@ hse_fill(const Box& bx, Array4 const& adv, } } - bool converged_hse = false; + [[maybe_unused]] bool converged_hse = false; Real p_want; Real drho; @@ -304,7 +304,7 @@ hse_fill(const Box& bx, Array4 const& adv, } } - bool converged_hse = false; + [[maybe_unused]] bool converged_hse = false; Real p_want; Real drho; @@ -501,7 +501,7 @@ hse_fill(const Box& bx, Array4 const& adv, } } - bool converged_hse = false; + [[maybe_unused]] bool converged_hse = false; Real p_want; Real drho; @@ -694,7 +694,7 @@ hse_fill(const Box& bx, Array4 const& adv, } } - bool converged_hse = false; + [[maybe_unused]] bool converged_hse = false; Real p_want; Real drho; @@ -890,7 +890,7 @@ hse_fill(const Box& bx, Array4 const& adv, } } - bool converged_hse = false; + [[maybe_unused]] bool converged_hse = false; Real p_want; Real drho; diff --git a/Source/sdc/Castro_sdc_util.H b/Source/sdc/Castro_sdc_util.H index 7f84100c57..ef632056b5 100644 --- a/Source/sdc/Castro_sdc_util.H +++ b/Source/sdc/Castro_sdc_util.H @@ -108,7 +108,7 @@ f_sdc_jac(const Real dt_m, // returned by single_zone_react_source, since it is // more consistent T from e - eos_t eos_state; + eos_extra_t eos_state; eos_state.rho = U_full[URHO]; eos_state.T = T_old; // initial guess for (int n = 0; n < NumSpec; ++n) { diff --git a/Source/sdc/vode_rhs_true_sdc.H b/Source/sdc/vode_rhs_true_sdc.H index 5a1c18e86c..a01b11b9be 100644 --- a/Source/sdc/vode_rhs_true_sdc.H +++ b/Source/sdc/vode_rhs_true_sdc.H @@ -109,7 +109,7 @@ void jac (const Real time, burn_t& burn_state, DvodeT& vode_state, MatrixType& p // returned by single_zone_react_source, since it is // more consistent T from e - eos_t eos_state; + eos_extra_t eos_state; eos_state.rho = U_full[URHO]; eos_state.T = burn_state.T; // initial guess for (int n = 0; n < NumSpec; n++) { diff --git a/external/Microphysics b/external/Microphysics index bcf06921ae..37551260a0 160000 --- a/external/Microphysics +++ b/external/Microphysics @@ -1 +1 @@ -Subproject commit bcf06921ae6144e0f1d24ceb8181af43357e9f84 +Subproject commit 37551260a0c17a85dd6d08fc289a64c38dfa46a6 diff --git a/external/amrex b/external/amrex index 02ce32baa5..60b23d9d54 160000 --- a/external/amrex +++ b/external/amrex @@ -1 +1 @@ -Subproject commit 02ce32baa5f243bbf012ba2824cc4c44682c6748 +Subproject commit 60b23d9d54287d4b4616add69783ab19729664af