From fe81e5d363dea3c5eaaf05abe28f9b243b1c25c2 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Wed, 26 Jun 2024 07:52:24 -0400 Subject: [PATCH 01/10] fix some namespace issues + add missing headers to runtime_parameters.H this will make it easier to use this for non-castro things (like exact_riemann) --- Source/driver/runtime_parameters.H | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/Source/driver/runtime_parameters.H b/Source/driver/runtime_parameters.H index 0619070811..4881668f2b 100644 --- a/Source/driver/runtime_parameters.H +++ b/Source/driver/runtime_parameters.H @@ -1,6 +1,7 @@ #ifndef RUNTIME_PARAMETERS_H #define RUNTIME_PARAMETERS_H +#include #include #include @@ -25,34 +26,34 @@ initialize_cpp_runparams() { { - ParmParse pp("castro"); + amrex::ParmParse pp("castro"); #include } #ifdef AMREX_PARTICLES { - ParmParse pp("particles"); + amrex::ParmParse pp("particles"); #include } #endif #ifdef DIFFUSION { - ParmParse pp("diffusion"); + amrex::ParmParse pp("diffusion"); #include } #endif #ifdef GRAVITY { - ParmParse pp("gravity"); + amrex::ParmParse pp("gravity"); #include } #endif #ifdef RADIATION { - ParmParse pp("radiation"); + amrex::ParmParse pp("radiation"); #include } #endif @@ -102,14 +103,14 @@ validate_runparams() for (const auto& nm: check_namespaces) { // "castro" - if (ParmParse::hasUnusedInputs(nm)) { + if (amrex::ParmParse::hasUnusedInputs(nm)) { amrex::Print() << "Warning: the following " + nm + ".* parameters are ignored\n"; - auto unused = ParmParse::getUnusedInputs(nm); + auto unused = amrex::ParmParse::getUnusedInputs(nm); for (const auto& p: unused) { amrex::Print() << p << "\n"; } amrex::Print() << std::endl; - if (abort_on_invalid_params) { + if (castro::abort_on_invalid_params) { amrex::Error("Error: invalid parameters"); } } From 000768fca4a2332245d320f32b90a57d4d882954 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 10:45:59 -0400 Subject: [PATCH 02/10] more work on getting castro params into exact Riemann --- Exec/Make.auto_source | 8 ++++- Source/driver/parse_castro_params.py | 13 +++++--- Util/exact_riemann/GNUmakefile | 7 ++++ Util/exact_riemann/main.cpp | 49 +++++++++++++++++++++------- 4 files changed, 60 insertions(+), 17 deletions(-) diff --git a/Exec/Make.auto_source b/Exec/Make.auto_source index e9e93fc284..dc3dfe194b 100644 --- a/Exec/Make.auto_source +++ b/Exec/Make.auto_source @@ -72,9 +72,15 @@ CPP_PARAMETERS := $(TOP)/Source/driver/_cpp_parameters $(CASTRO_AUTO_SOURCE_DIR)/runtime_params.cpp: $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H +STRUCT_USE_CASTRO_CLASS ?= TRUE +PARSE_ARGS := +ifneq ($(STRUCT_USE_CASTRO_CLASS), TRUE) + PARSE_ARGS += --without-castro-class +endif + $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H: $(CPP_PARAMETERS) @if [ ! -d $(CASTRO_AUTO_SOURCE_DIR) ]; then mkdir -p $(CASTRO_AUTO_SOURCE_DIR); fi - PYTHONPATH=$(MICROPHYSICS_HOME)/util/build_scripts $(TOP)/Source/driver/parse_castro_params.py -o $(CASTRO_AUTO_SOURCE_DIR) $(CPP_PARAMETERS) + PYTHONPATH=$(MICROPHYSICS_HOME)/util/build_scripts $(TOP)/Source/driver/parse_castro_params.py $(PARSE_ARGS) -o $(CASTRO_AUTO_SOURCE_DIR) $(CPP_PARAMETERS) # for debugging test_cxx_params: $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H diff --git a/Source/driver/parse_castro_params.py b/Source/driver/parse_castro_params.py index 89f91a8ea7..4d869d09dd 100755 --- a/Source/driver/parse_castro_params.py +++ b/Source/driver/parse_castro_params.py @@ -129,7 +129,7 @@ def read_param_file(infile): return params -def write_headers_and_source(params, out_directory, struct_name): +def write_headers_and_source(params, out_directory, struct_name, without_castro_class): # output @@ -176,19 +176,23 @@ def write_headers_and_source(params, out_directory, struct_name): cq.write(CWARNING) + class_name = "Castro" + if without_castro_class: + class_name = None + for ifdef in ifdefs: if ifdef is None: for p in [q for q in params_nm if q.ifdef is None]: cq.write(p.get_default_string()) cq.write(p.get_query_string()) - cq.write(p.get_query_struct_string(struct_name=struct_name, class_name="Castro")) + cq.write(p.get_query_struct_string(struct_name=struct_name, class_name=class_name)) cq.write("\n") else: cq.write(f"#ifdef {ifdef}\n") for p in [q for q in params_nm if q.ifdef == ifdef]: cq.write(p.get_default_string()) cq.write(p.get_query_string()) - cq.write(p.get_query_struct_string(struct_name=struct_name, class_name="Castro")) + cq.write(p.get_query_struct_string(struct_name=struct_name, class_name=class_name)) cq.write("\n") cq.write("#endif\n") cq.write("\n") @@ -305,13 +309,14 @@ def main(): help="output directory for the generated files") parser.add_argument("-s", type=str, default="params", help="name for the name struct that will hold the parameters") + parser.add_argument("--without-castro-class", action="store_true", help="don't include Castro:: in the struct namespace") parser.add_argument("input_file", type=str, nargs=1, help="input file containing the list of parameters we will define") args = parser.parse_args() p = read_param_file(args.input_file[0]) - write_headers_and_source(p, args.o, args.s) + write_headers_and_source(p, args.o, args.s, args.without_castro_class) if __name__ == "__main__": main() diff --git a/Util/exact_riemann/GNUmakefile b/Util/exact_riemann/GNUmakefile index 76bc8e8a9f..8be2709e1a 100644 --- a/Util/exact_riemann/GNUmakefile +++ b/Util/exact_riemann/GNUmakefile @@ -23,9 +23,16 @@ EOS_DIR := helmholtz NETWORK_DIR := general_null NETWORK_INPUTS = ignition.net +# for the Castro runtime parameters, we don't want to use Castro:: +STRUCT_USE_CASTRO_CLASS := FALSE + EXTERN_SEARCH += . Bpack := ./Make.package Blocs := . +# we explicitly want runtime_parameters.H so we have access to +# Castro's runtime parameter system +Blocs += $(CASTRO_HOME)/Source/driver + include $(CASTRO_HOME)/Exec/Make.Castro diff --git a/Util/exact_riemann/main.cpp b/Util/exact_riemann/main.cpp index c2b3ed296d..6e5a5f7150 100644 --- a/Util/exact_riemann/main.cpp +++ b/Util/exact_riemann/main.cpp @@ -13,27 +13,52 @@ #include -int main(int argc, char *argv[]) { +#include +#include +#include - amrex::Initialize(argc, argv); +// create a dummy Castro class just to hold the runtime parameters struct - std::cout << "starting the exact Riemann solver..." << std::endl; +class Castro { + + public: + + static params_t params; - // initialize the external runtime parameters in C++ + void doit() { - init_prob_parameters(); + initialize_cpp_runparams(); - init_extern_parameters(); + // initialize the external runtime parameters in C++ - // now initialize the C++ Microphysics + init_prob_parameters(); + + init_extern_parameters(); + + // now initialize the C++ Microphysics #ifdef REACTIONS - network_init(); + network_init(); #endif - amrex::Real small_temp = 1.e-200; - amrex::Real small_dens = 1.e-200; - eos_init(small_temp, small_dens); + amrex::Real small_temp = 1.e-200; + amrex::Real small_dens = 1.e-200; + eos_init(small_temp, small_dens); + + exact_riemann(); + + } + +}; + + +int main(int argc, char *argv[]) { + + amrex::Initialize(argc, argv); + + std::cout << "starting the exact Riemann solver..." << std::endl; - exact_riemann(); + // initialize the Castro runtime parameters + Castro castro; + castro.doit(); } From d01fa53e94e32ea453e73ef633bf7949088ad9f1 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 12:32:33 -0400 Subject: [PATCH 03/10] this seems to work --- Util/exact_riemann/Make.package | 3 +++ Util/exact_riemann/main.cpp | 48 ++++++++++++--------------------- 2 files changed, 20 insertions(+), 31 deletions(-) diff --git a/Util/exact_riemann/Make.package b/Util/exact_riemann/Make.package index b823fe2404..2cc1befab7 100644 --- a/Util/exact_riemann/Make.package +++ b/Util/exact_riemann/Make.package @@ -6,3 +6,6 @@ CEXE_headers += riemann_sample.H CEXE_headers += riemann_support.H CEXE_sources += extern_parameters.cpp CEXE_headers += extern_parameters.H + +# automatically generated in tmp_build_dir +CEXE_sources += runtime_params.cpp \ No newline at end of file diff --git a/Util/exact_riemann/main.cpp b/Util/exact_riemann/main.cpp index 6e5a5f7150..7f75f6fc36 100644 --- a/Util/exact_riemann/main.cpp +++ b/Util/exact_riemann/main.cpp @@ -14,51 +14,37 @@ #include #include -#include #include -// create a dummy Castro class just to hold the runtime parameters struct +// a global struct to hold the params +#include -class Castro { +int main(int argc, char *argv[]) { - public: + amrex::Initialize(argc, argv); - static params_t params; + std::cout << "starting the exact Riemann solver..." << std::endl; - void doit() { + // initialize the Castro runtime parameters - initialize_cpp_runparams(); + amrex::ParmParse pp("castro"); +#include - // initialize the external runtime parameters in C++ + // initialize the external runtime parameters in C++ - init_prob_parameters(); + init_prob_parameters(); - init_extern_parameters(); + init_extern_parameters(); - // now initialize the C++ Microphysics + // now initialize the C++ Microphysics #ifdef REACTIONS - network_init(); + network_init(); #endif - amrex::Real small_temp = 1.e-200; - amrex::Real small_dens = 1.e-200; - eos_init(small_temp, small_dens); - - exact_riemann(); - - } - -}; - + amrex::Real small_temp = 1.e-200; + amrex::Real small_dens = 1.e-200; + eos_init(small_temp, small_dens); -int main(int argc, char *argv[]) { - - amrex::Initialize(argc, argv); - - std::cout << "starting the exact Riemann solver..." << std::endl; - - // initialize the Castro runtime parameters + exact_riemann(); - Castro castro; - castro.doit(); } From a4a16d352c95772ed11fce64d28289df9a5007ec Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 12:48:45 -0400 Subject: [PATCH 04/10] this works now --- Util/exact_riemann/_prob_params | 2 -- Util/exact_riemann/exact_riemann.H | 4 +++- Util/exact_riemann/inputs.stellarcollapse | 4 ++-- Util/exact_riemann/inputs.test1.helm | 1 + Util/exact_riemann/riemann_rarefaction.H | 8 ++++---- Util/exact_riemann/riemann_sample.H | 6 +++--- Util/exact_riemann/riemann_shock.H | 4 ++-- Util/exact_riemann/riemann_star_state.H | 4 ++-- 8 files changed, 17 insertions(+), 16 deletions(-) diff --git a/Util/exact_riemann/_prob_params b/Util/exact_riemann/_prob_params index 9c43f0d3b5..b8cf88be92 100644 --- a/Util/exact_riemann/_prob_params +++ b/Util/exact_riemann/_prob_params @@ -19,5 +19,3 @@ t real 0.2d0 y npts integer 128 y use_Tinit integer 0 y - -initial_temp_guess real 1.0d5 y diff --git a/Util/exact_riemann/exact_riemann.H b/Util/exact_riemann/exact_riemann.H index 39deafe8d4..52b44db51a 100644 --- a/Util/exact_riemann/exact_riemann.H +++ b/Util/exact_riemann/exact_riemann.H @@ -5,6 +5,8 @@ #include +#include + #include #include #include @@ -88,7 +90,7 @@ exact_riemann() { for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_s[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/inputs.stellarcollapse b/Util/exact_riemann/inputs.stellarcollapse index 92dd346435..412b02c272 100644 --- a/Util/exact_riemann/inputs.stellarcollapse +++ b/Util/exact_riemann/inputs.stellarcollapse @@ -8,8 +8,6 @@ problem.T_r = 5.e10 problem.use_Tinit = 1 -problem.initial_temp_guess = 1.e5 - problem.riemann_max_iter = 20 problem.xmin = 0.0 @@ -20,5 +18,7 @@ problem.t = 1.e-6 problem.npts = 128 +castro.T_guess = 1.e5 + eos.eos_file = LS220_240r_140t_50y_analmu_20120628_SVNr28.h5 diff --git a/Util/exact_riemann/inputs.test1.helm b/Util/exact_riemann/inputs.test1.helm index db1db7d1ce..b36d2f0bab 100644 --- a/Util/exact_riemann/inputs.test1.helm +++ b/Util/exact_riemann/inputs.test1.helm @@ -18,3 +18,4 @@ problem.t = 0.0008 problem.npts = 128 +castro.T_guess = 1.e5 diff --git a/Util/exact_riemann/riemann_rarefaction.H b/Util/exact_riemann/riemann_rarefaction.H index c5d98c6f87..fd2e266342 100644 --- a/Util/exact_riemann/riemann_rarefaction.H +++ b/Util/exact_riemann/riemann_rarefaction.H @@ -185,7 +185,7 @@ rarefaction(const amrex::Real pstar, // initial guess at integration step amrex::Real dp = (pstar - p_s) / static_cast(100); - amrex::Real T = problem::initial_temp_guess; + amrex::Real T = castro::T_guess; // adaptive RK4 loop @@ -253,7 +253,7 @@ rarefaction(const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -311,7 +311,7 @@ rarefaction_to_u(const amrex::Real rho_s, const amrex::Real u_s, const amrex::Re for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -389,7 +389,7 @@ rarefaction_to_u(const amrex::Real rho_s, const amrex::Real u_s, const amrex::Re for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_sample.H b/Util/exact_riemann/riemann_sample.H index 9db07b34df..8f0350fa4c 100644 --- a/Util/exact_riemann/riemann_sample.H +++ b/Util/exact_riemann/riemann_sample.H @@ -25,7 +25,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_l[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -36,7 +36,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_r[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -129,7 +129,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_shock.H b/Util/exact_riemann/riemann_shock.H index c91ac466de..3a74b6e3d3 100644 --- a/Util/exact_riemann/riemann_shock.H +++ b/Util/exact_riemann/riemann_shock.H @@ -26,7 +26,7 @@ W_s_shock(const amrex::Real W_s, const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -113,7 +113,7 @@ shock(const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_star_state.H b/Util/exact_riemann/riemann_star_state.H index 67eeabada7..2448247129 100644 --- a/Util/exact_riemann/riemann_star_state.H +++ b/Util/exact_riemann/riemann_star_state.H @@ -28,7 +28,7 @@ riemann_star_state(const amrex::Real rho_l, const amrex::Real u_l, const amrex:: for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_l[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -42,7 +42,7 @@ riemann_star_state(const amrex::Real rho_l, const amrex::Real u_l, const amrex:: for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_r[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); From b0b48bbd6a0aeb69e6d482933307218b2df06656 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 12:51:41 -0400 Subject: [PATCH 05/10] add missing file --- Util/exact_riemann/struct_params.H | 8 ++++++++ 1 file changed, 8 insertions(+) create mode 100644 Util/exact_riemann/struct_params.H diff --git a/Util/exact_riemann/struct_params.H b/Util/exact_riemann/struct_params.H new file mode 100644 index 0000000000..8ad1776283 --- /dev/null +++ b/Util/exact_riemann/struct_params.H @@ -0,0 +1,8 @@ +#ifndef STRUCT_PARAMS_H +#define STRUCT_PARAMS_H + +#include + +inline params_t params; + +#endif From 9986cc8887736ea9a3c9d542e77874023f79883c Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 12:52:42 -0400 Subject: [PATCH 06/10] fix some namespace issues + add missing headers to runtime_parameters.H (#2891) this will make it easier to use this for non-castro things (like exact_riemann) --- Source/driver/runtime_parameters.H | 17 +++++++++-------- 1 file changed, 9 insertions(+), 8 deletions(-) diff --git a/Source/driver/runtime_parameters.H b/Source/driver/runtime_parameters.H index 0619070811..4881668f2b 100644 --- a/Source/driver/runtime_parameters.H +++ b/Source/driver/runtime_parameters.H @@ -1,6 +1,7 @@ #ifndef RUNTIME_PARAMETERS_H #define RUNTIME_PARAMETERS_H +#include #include #include @@ -25,34 +26,34 @@ initialize_cpp_runparams() { { - ParmParse pp("castro"); + amrex::ParmParse pp("castro"); #include } #ifdef AMREX_PARTICLES { - ParmParse pp("particles"); + amrex::ParmParse pp("particles"); #include } #endif #ifdef DIFFUSION { - ParmParse pp("diffusion"); + amrex::ParmParse pp("diffusion"); #include } #endif #ifdef GRAVITY { - ParmParse pp("gravity"); + amrex::ParmParse pp("gravity"); #include } #endif #ifdef RADIATION { - ParmParse pp("radiation"); + amrex::ParmParse pp("radiation"); #include } #endif @@ -102,14 +103,14 @@ validate_runparams() for (const auto& nm: check_namespaces) { // "castro" - if (ParmParse::hasUnusedInputs(nm)) { + if (amrex::ParmParse::hasUnusedInputs(nm)) { amrex::Print() << "Warning: the following " + nm + ".* parameters are ignored\n"; - auto unused = ParmParse::getUnusedInputs(nm); + auto unused = amrex::ParmParse::getUnusedInputs(nm); for (const auto& p: unused) { amrex::Print() << p << "\n"; } amrex::Print() << std::endl; - if (abort_on_invalid_params) { + if (castro::abort_on_invalid_params) { amrex::Error("Error: invalid parameters"); } } From c47966edd69b2ec99f4a9b4b70bdda166a500a1c Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 14:34:47 -0400 Subject: [PATCH 07/10] move the riemann_util functions into Castro_mol.cpp (#2890) these are only used in the MOL reconstruction and not general Riemann functions --- Source/hydro/Castro_mol.cpp | 144 +++++++++++++++++++++++++++++++ Source/hydro/Make.package | 1 - Source/hydro/riemann_util.cpp | 154 ---------------------------------- 3 files changed, 144 insertions(+), 155 deletions(-) delete mode 100644 Source/hydro/riemann_util.cpp diff --git a/Source/hydro/Castro_mol.cpp b/Source/hydro/Castro_mol.cpp index eb6314769c..fe2004a60c 100644 --- a/Source/hydro/Castro_mol.cpp +++ b/Source/hydro/Castro_mol.cpp @@ -8,6 +8,10 @@ #include #endif +#ifdef HYBRID_MOMENTUM +#include +#endif + #include #include @@ -321,3 +325,143 @@ Castro::mol_diffusive_flux(const Box& bx, }); } + + +void +Castro::compute_flux_from_q(const Box& bx, + Array4 const& qint, + Array4 const& F, + const int idir) { + + // given a primitive state, compute the flux in direction idir + // + + int iu, iv1, iv2; + int im1, im2, im3; + + auto coord = geom.Coord(); + auto mom_check = mom_flux_has_p(idir, idir, coord); + + if (idir == 0) { + iu = QU; + iv1 = QV; + iv2 = QW; + im1 = UMX; + im2 = UMY; + im3 = UMZ; + + } else if (idir == 1) { + iu = QV; + iv1 = QU; + iv2 = QW; + im1 = UMY; + im2 = UMX; + im3 = UMZ; + + } else { + iu = QW; + iv1 = QU; + iv2 = QV; + im1 = UMZ; + im2 = UMX; + im3 = UMY; + } + +#ifdef HYBRID_MOMENTUM + GeometryData geomdata = geom.data(); +#endif + + amrex::ParallelFor(bx, + [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept + { + + Real u_adv = qint(i,j,k,iu); + Real rhoeint = qint(i,j,k,QREINT); + + // Compute fluxes, order as conserved state (not q) + F(i,j,k,URHO) = qint(i,j,k,QRHO)*u_adv; + + F(i,j,k,im1) = F(i,j,k,URHO)*qint(i,j,k,iu); + if (mom_check) { + F(i,j,k,im1) += qint(i,j,k,QPRES); + } + F(i,j,k,im2) = F(i,j,k,URHO)*qint(i,j,k,iv1); + F(i,j,k,im3) = F(i,j,k,URHO)*qint(i,j,k,iv2); + + Real rhoetot = rhoeint + 0.5_rt * qint(i,j,k,QRHO)* + (qint(i,j,k,iu)*qint(i,j,k,iu) + + qint(i,j,k,iv1)*qint(i,j,k,iv1) + + qint(i,j,k,iv2)*qint(i,j,k,iv2)); + + F(i,j,k,UEDEN) = u_adv*(rhoetot + qint(i,j,k,QPRES)); + F(i,j,k,UEINT) = u_adv*rhoeint; + + F(i,j,k,UTEMP) = 0.0; +#ifdef SHOCK_VAR + F(i,j,k,USHK) = 0.0; +#endif + +#ifdef NSE_NET + F(i,j,k,UMUP) = 0.0; + F(i,j,k,UMUN) = 0.0; +#endif + // passively advected quantities + for (int ipassive = 0; ipassive < npassive; ipassive++) { + int n = upassmap(ipassive); + int nqp = qpassmap(ipassive); + + F(i,j,k,n) = F(i,j,k,URHO)*qint(i,j,k,nqp); + } + +#ifdef HYBRID_MOMENTUM + // the hybrid routine uses the Godunov indices, not the full NQ state + GpuArray qgdnv_zone; + qgdnv_zone[GDRHO] = qint(i,j,k,QRHO); + qgdnv_zone[GDU] = qint(i,j,k,QU); + qgdnv_zone[GDV] = qint(i,j,k,QV); + qgdnv_zone[GDW] = qint(i,j,k,QW); + qgdnv_zone[GDPRES] = qint(i,j,k,QPRES); + GpuArray F_zone; + for (int n = 0; n < NUM_STATE; n++) { + F_zone[n] = F(i,j,k,n); + } + compute_hybrid_flux(qgdnv_zone, geomdata, idir, i, j, k, F_zone); + for (int n = 0; n < NUM_STATE; n++) { + F(i,j,k,n) = F_zone[n]; + } +#endif + }); +} + +void +Castro::store_godunov_state(const Box& bx, + Array4 const& qint, +#ifdef RADIATION + Array4 const& lambda, +#endif + Array4 const& qgdnv) { + + // this copies the full interface state (NQ -- one for each primitive + // variable) over to a smaller subset of size NGDNV for use later in the + // hydro advancement. + + amrex::ParallelFor(bx, + [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept + { + + +#ifdef HYBRID_MOMENTUM + qgdnv(i,j,k,GDRHO) = qint(i,j,k,QRHO); +#endif + qgdnv(i,j,k,GDU) = qint(i,j,k,QU); + qgdnv(i,j,k,GDV) = qint(i,j,k,QV); + qgdnv(i,j,k,GDW) = qint(i,j,k,QW); + qgdnv(i,j,k,GDPRES) = qint(i,j,k,QPRES); +#ifdef RADIATION + for (int g = 0; g < NGROUPS; g++) { + qgdnv(i,j,k,GDLAMS+g) = lambda(i,j,k,g); + qgdnv(i,j,k,GDERADS+g) = qint(i,j,k,QRAD+g); + } +#endif + }); +} diff --git a/Source/hydro/Make.package b/Source/hydro/Make.package index b4d67b48fa..6b6ba6c7e8 100644 --- a/Source/hydro/Make.package +++ b/Source/hydro/Make.package @@ -30,7 +30,6 @@ CEXE_sources += riemann.cpp CEXE_headers += HLL_solvers.H CEXE_headers += riemann_solvers.H CEXE_headers += riemann_2shock_solvers.H -CEXE_sources += riemann_util.cpp CEXE_headers += riemann_type.H CEXE_headers += slope.H CEXE_headers += reconstruction.H diff --git a/Source/hydro/riemann_util.cpp b/Source/hydro/riemann_util.cpp deleted file mode 100644 index 2801dad8c2..0000000000 --- a/Source/hydro/riemann_util.cpp +++ /dev/null @@ -1,154 +0,0 @@ -#include -#include - -#ifdef RADIATION -#include -#include -#endif - -#ifdef HYBRID_MOMENTUM -#include -#endif - -#include - -using namespace amrex; - -void -Castro::compute_flux_from_q(const Box& bx, - Array4 const& qint, - Array4 const& F, - const int idir) { - - // given a primitive state, compute the flux in direction idir - // - - int iu, iv1, iv2; - int im1, im2, im3; - - auto coord = geom.Coord(); - auto mom_check = mom_flux_has_p(idir, idir, coord); - - if (idir == 0) { - iu = QU; - iv1 = QV; - iv2 = QW; - im1 = UMX; - im2 = UMY; - im3 = UMZ; - - } else if (idir == 1) { - iu = QV; - iv1 = QU; - iv2 = QW; - im1 = UMY; - im2 = UMX; - im3 = UMZ; - - } else { - iu = QW; - iv1 = QU; - iv2 = QV; - im1 = UMZ; - im2 = UMX; - im3 = UMY; - } - -#ifdef HYBRID_MOMENTUM - GeometryData geomdata = geom.data(); -#endif - - amrex::ParallelFor(bx, - [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept - { - - Real u_adv = qint(i,j,k,iu); - Real rhoeint = qint(i,j,k,QREINT); - - // Compute fluxes, order as conserved state (not q) - F(i,j,k,URHO) = qint(i,j,k,QRHO)*u_adv; - - F(i,j,k,im1) = F(i,j,k,URHO)*qint(i,j,k,iu); - if (mom_check) { - F(i,j,k,im1) += qint(i,j,k,QPRES); - } - F(i,j,k,im2) = F(i,j,k,URHO)*qint(i,j,k,iv1); - F(i,j,k,im3) = F(i,j,k,URHO)*qint(i,j,k,iv2); - - Real rhoetot = rhoeint + 0.5_rt * qint(i,j,k,QRHO)* - (qint(i,j,k,iu)*qint(i,j,k,iu) + - qint(i,j,k,iv1)*qint(i,j,k,iv1) + - qint(i,j,k,iv2)*qint(i,j,k,iv2)); - - F(i,j,k,UEDEN) = u_adv*(rhoetot + qint(i,j,k,QPRES)); - F(i,j,k,UEINT) = u_adv*rhoeint; - - F(i,j,k,UTEMP) = 0.0; -#ifdef SHOCK_VAR - F(i,j,k,USHK) = 0.0; -#endif - -#ifdef NSE_NET - F(i,j,k,UMUP) = 0.0; - F(i,j,k,UMUN) = 0.0; -#endif - // passively advected quantities - for (int ipassive = 0; ipassive < npassive; ipassive++) { - int n = upassmap(ipassive); - int nqp = qpassmap(ipassive); - - F(i,j,k,n) = F(i,j,k,URHO)*qint(i,j,k,nqp); - } - -#ifdef HYBRID_MOMENTUM - // the hybrid routine uses the Godunov indices, not the full NQ state - GpuArray qgdnv_zone; - qgdnv_zone[GDRHO] = qint(i,j,k,QRHO); - qgdnv_zone[GDU] = qint(i,j,k,QU); - qgdnv_zone[GDV] = qint(i,j,k,QV); - qgdnv_zone[GDW] = qint(i,j,k,QW); - qgdnv_zone[GDPRES] = qint(i,j,k,QPRES); - GpuArray F_zone; - for (int n = 0; n < NUM_STATE; n++) { - F_zone[n] = F(i,j,k,n); - } - compute_hybrid_flux(qgdnv_zone, geomdata, idir, i, j, k, F_zone); - for (int n = 0; n < NUM_STATE; n++) { - F(i,j,k,n) = F_zone[n]; - } -#endif - }); -} - -void -Castro::store_godunov_state(const Box& bx, - Array4 const& qint, -#ifdef RADIATION - Array4 const& lambda, -#endif - Array4 const& qgdnv) { - - // this copies the full interface state (NQ -- one for each primitive - // variable) over to a smaller subset of size NGDNV for use later in the - // hydro advancement. - - amrex::ParallelFor(bx, - [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept - { - - -#ifdef HYBRID_MOMENTUM - qgdnv(i,j,k,GDRHO) = qint(i,j,k,QRHO); -#endif - qgdnv(i,j,k,GDU) = qint(i,j,k,QU); - qgdnv(i,j,k,GDV) = qint(i,j,k,QV); - qgdnv(i,j,k,GDW) = qint(i,j,k,QW); - qgdnv(i,j,k,GDPRES) = qint(i,j,k,QPRES); -#ifdef RADIATION - for (int g = 0; g < NGROUPS; g++) { - qgdnv(i,j,k,GDLAMS+g) = lambda(i,j,k,g); - qgdnv(i,j,k,GDERADS+g) = qint(i,j,k,QRAD+g); - } -#endif - }); -} From 8b810445d5b1a6158197f1fc14cee0257bb21126 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Thu, 27 Jun 2024 15:27:02 -0400 Subject: [PATCH 08/10] unused header --- Util/exact_riemann/exact_riemann.H | 1 - 1 file changed, 1 deletion(-) diff --git a/Util/exact_riemann/exact_riemann.H b/Util/exact_riemann/exact_riemann.H index 52b44db51a..e99bf69a4d 100644 --- a/Util/exact_riemann/exact_riemann.H +++ b/Util/exact_riemann/exact_riemann.H @@ -9,7 +9,6 @@ #include #include -#include using namespace amrex::literals; From 97fe9364a30acd09f3a9386c56d2c18be24f36fb Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Mon, 1 Jul 2024 06:59:30 -0400 Subject: [PATCH 09/10] allow the exact_riemann test to work with Castro runtime parameters (#2892) This required adding a option to the parse_castro_params.py to skip the Castro:: scope on the runtime parameter struct. We also now use castro::T_guess as the initial temperature guess. --- Exec/Make.auto_source | 8 +++++++- Source/driver/parse_castro_params.py | 13 +++++++++---- Util/exact_riemann/GNUmakefile | 7 +++++++ Util/exact_riemann/Make.package | 3 +++ Util/exact_riemann/_prob_params | 2 -- Util/exact_riemann/exact_riemann.H | 5 +++-- Util/exact_riemann/inputs.stellarcollapse | 4 ++-- Util/exact_riemann/inputs.test1.helm | 1 + Util/exact_riemann/main.cpp | 11 +++++++++++ Util/exact_riemann/riemann_rarefaction.H | 8 ++++---- Util/exact_riemann/riemann_sample.H | 6 +++--- Util/exact_riemann/riemann_shock.H | 4 ++-- Util/exact_riemann/riemann_star_state.H | 4 ++-- Util/exact_riemann/struct_params.H | 8 ++++++++ 14 files changed, 62 insertions(+), 22 deletions(-) create mode 100644 Util/exact_riemann/struct_params.H diff --git a/Exec/Make.auto_source b/Exec/Make.auto_source index e9e93fc284..dc3dfe194b 100644 --- a/Exec/Make.auto_source +++ b/Exec/Make.auto_source @@ -72,9 +72,15 @@ CPP_PARAMETERS := $(TOP)/Source/driver/_cpp_parameters $(CASTRO_AUTO_SOURCE_DIR)/runtime_params.cpp: $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H +STRUCT_USE_CASTRO_CLASS ?= TRUE +PARSE_ARGS := +ifneq ($(STRUCT_USE_CASTRO_CLASS), TRUE) + PARSE_ARGS += --without-castro-class +endif + $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H: $(CPP_PARAMETERS) @if [ ! -d $(CASTRO_AUTO_SOURCE_DIR) ]; then mkdir -p $(CASTRO_AUTO_SOURCE_DIR); fi - PYTHONPATH=$(MICROPHYSICS_HOME)/util/build_scripts $(TOP)/Source/driver/parse_castro_params.py -o $(CASTRO_AUTO_SOURCE_DIR) $(CPP_PARAMETERS) + PYTHONPATH=$(MICROPHYSICS_HOME)/util/build_scripts $(TOP)/Source/driver/parse_castro_params.py $(PARSE_ARGS) -o $(CASTRO_AUTO_SOURCE_DIR) $(CPP_PARAMETERS) # for debugging test_cxx_params: $(CASTRO_AUTO_SOURCE_DIR)/castro_params.H diff --git a/Source/driver/parse_castro_params.py b/Source/driver/parse_castro_params.py index 89f91a8ea7..4d869d09dd 100755 --- a/Source/driver/parse_castro_params.py +++ b/Source/driver/parse_castro_params.py @@ -129,7 +129,7 @@ def read_param_file(infile): return params -def write_headers_and_source(params, out_directory, struct_name): +def write_headers_and_source(params, out_directory, struct_name, without_castro_class): # output @@ -176,19 +176,23 @@ def write_headers_and_source(params, out_directory, struct_name): cq.write(CWARNING) + class_name = "Castro" + if without_castro_class: + class_name = None + for ifdef in ifdefs: if ifdef is None: for p in [q for q in params_nm if q.ifdef is None]: cq.write(p.get_default_string()) cq.write(p.get_query_string()) - cq.write(p.get_query_struct_string(struct_name=struct_name, class_name="Castro")) + cq.write(p.get_query_struct_string(struct_name=struct_name, class_name=class_name)) cq.write("\n") else: cq.write(f"#ifdef {ifdef}\n") for p in [q for q in params_nm if q.ifdef == ifdef]: cq.write(p.get_default_string()) cq.write(p.get_query_string()) - cq.write(p.get_query_struct_string(struct_name=struct_name, class_name="Castro")) + cq.write(p.get_query_struct_string(struct_name=struct_name, class_name=class_name)) cq.write("\n") cq.write("#endif\n") cq.write("\n") @@ -305,13 +309,14 @@ def main(): help="output directory for the generated files") parser.add_argument("-s", type=str, default="params", help="name for the name struct that will hold the parameters") + parser.add_argument("--without-castro-class", action="store_true", help="don't include Castro:: in the struct namespace") parser.add_argument("input_file", type=str, nargs=1, help="input file containing the list of parameters we will define") args = parser.parse_args() p = read_param_file(args.input_file[0]) - write_headers_and_source(p, args.o, args.s) + write_headers_and_source(p, args.o, args.s, args.without_castro_class) if __name__ == "__main__": main() diff --git a/Util/exact_riemann/GNUmakefile b/Util/exact_riemann/GNUmakefile index 76bc8e8a9f..8be2709e1a 100644 --- a/Util/exact_riemann/GNUmakefile +++ b/Util/exact_riemann/GNUmakefile @@ -23,9 +23,16 @@ EOS_DIR := helmholtz NETWORK_DIR := general_null NETWORK_INPUTS = ignition.net +# for the Castro runtime parameters, we don't want to use Castro:: +STRUCT_USE_CASTRO_CLASS := FALSE + EXTERN_SEARCH += . Bpack := ./Make.package Blocs := . +# we explicitly want runtime_parameters.H so we have access to +# Castro's runtime parameter system +Blocs += $(CASTRO_HOME)/Source/driver + include $(CASTRO_HOME)/Exec/Make.Castro diff --git a/Util/exact_riemann/Make.package b/Util/exact_riemann/Make.package index b823fe2404..2cc1befab7 100644 --- a/Util/exact_riemann/Make.package +++ b/Util/exact_riemann/Make.package @@ -6,3 +6,6 @@ CEXE_headers += riemann_sample.H CEXE_headers += riemann_support.H CEXE_sources += extern_parameters.cpp CEXE_headers += extern_parameters.H + +# automatically generated in tmp_build_dir +CEXE_sources += runtime_params.cpp \ No newline at end of file diff --git a/Util/exact_riemann/_prob_params b/Util/exact_riemann/_prob_params index 9c43f0d3b5..b8cf88be92 100644 --- a/Util/exact_riemann/_prob_params +++ b/Util/exact_riemann/_prob_params @@ -19,5 +19,3 @@ t real 0.2d0 y npts integer 128 y use_Tinit integer 0 y - -initial_temp_guess real 1.0d5 y diff --git a/Util/exact_riemann/exact_riemann.H b/Util/exact_riemann/exact_riemann.H index 39deafe8d4..e99bf69a4d 100644 --- a/Util/exact_riemann/exact_riemann.H +++ b/Util/exact_riemann/exact_riemann.H @@ -5,9 +5,10 @@ #include +#include + #include #include -#include using namespace amrex::literals; @@ -88,7 +89,7 @@ exact_riemann() { for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_s[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/inputs.stellarcollapse b/Util/exact_riemann/inputs.stellarcollapse index 92dd346435..412b02c272 100644 --- a/Util/exact_riemann/inputs.stellarcollapse +++ b/Util/exact_riemann/inputs.stellarcollapse @@ -8,8 +8,6 @@ problem.T_r = 5.e10 problem.use_Tinit = 1 -problem.initial_temp_guess = 1.e5 - problem.riemann_max_iter = 20 problem.xmin = 0.0 @@ -20,5 +18,7 @@ problem.t = 1.e-6 problem.npts = 128 +castro.T_guess = 1.e5 + eos.eos_file = LS220_240r_140t_50y_analmu_20120628_SVNr28.h5 diff --git a/Util/exact_riemann/inputs.test1.helm b/Util/exact_riemann/inputs.test1.helm index db1db7d1ce..b36d2f0bab 100644 --- a/Util/exact_riemann/inputs.test1.helm +++ b/Util/exact_riemann/inputs.test1.helm @@ -18,3 +18,4 @@ problem.t = 0.0008 problem.npts = 128 +castro.T_guess = 1.e5 diff --git a/Util/exact_riemann/main.cpp b/Util/exact_riemann/main.cpp index c2b3ed296d..7f75f6fc36 100644 --- a/Util/exact_riemann/main.cpp +++ b/Util/exact_riemann/main.cpp @@ -13,12 +13,23 @@ #include +#include +#include + +// a global struct to hold the params +#include + int main(int argc, char *argv[]) { amrex::Initialize(argc, argv); std::cout << "starting the exact Riemann solver..." << std::endl; + // initialize the Castro runtime parameters + + amrex::ParmParse pp("castro"); +#include + // initialize the external runtime parameters in C++ init_prob_parameters(); diff --git a/Util/exact_riemann/riemann_rarefaction.H b/Util/exact_riemann/riemann_rarefaction.H index c5d98c6f87..fd2e266342 100644 --- a/Util/exact_riemann/riemann_rarefaction.H +++ b/Util/exact_riemann/riemann_rarefaction.H @@ -185,7 +185,7 @@ rarefaction(const amrex::Real pstar, // initial guess at integration step amrex::Real dp = (pstar - p_s) / static_cast(100); - amrex::Real T = problem::initial_temp_guess; + amrex::Real T = castro::T_guess; // adaptive RK4 loop @@ -253,7 +253,7 @@ rarefaction(const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -311,7 +311,7 @@ rarefaction_to_u(const amrex::Real rho_s, const amrex::Real u_s, const amrex::Re for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -389,7 +389,7 @@ rarefaction_to_u(const amrex::Real rho_s, const amrex::Real u_s, const amrex::Re for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_sample.H b/Util/exact_riemann/riemann_sample.H index 9db07b34df..8f0350fa4c 100644 --- a/Util/exact_riemann/riemann_sample.H +++ b/Util/exact_riemann/riemann_sample.H @@ -25,7 +25,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_l[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -36,7 +36,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_r[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -129,7 +129,7 @@ riemann_sample(const amrex::Real rho_l, const amrex::Real u_l, const amrex::Real for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_shock.H b/Util/exact_riemann/riemann_shock.H index c91ac466de..3a74b6e3d3 100644 --- a/Util/exact_riemann/riemann_shock.H +++ b/Util/exact_riemann/riemann_shock.H @@ -26,7 +26,7 @@ W_s_shock(const amrex::Real W_s, const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -113,7 +113,7 @@ shock(const amrex::Real pstar, for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/riemann_star_state.H b/Util/exact_riemann/riemann_star_state.H index 67eeabada7..2448247129 100644 --- a/Util/exact_riemann/riemann_star_state.H +++ b/Util/exact_riemann/riemann_star_state.H @@ -28,7 +28,7 @@ riemann_star_state(const amrex::Real rho_l, const amrex::Real u_l, const amrex:: for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_l[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); @@ -42,7 +42,7 @@ riemann_star_state(const amrex::Real rho_l, const amrex::Real u_l, const amrex:: for (int n = 0; n < NumSpec; ++n) { eos_state.xn[n] = xn_r[n]; } - eos_state.T = problem::initial_temp_guess; + eos_state.T = castro::T_guess; eos(eos_input_rp, eos_state); diff --git a/Util/exact_riemann/struct_params.H b/Util/exact_riemann/struct_params.H new file mode 100644 index 0000000000..8ad1776283 --- /dev/null +++ b/Util/exact_riemann/struct_params.H @@ -0,0 +1,8 @@ +#ifndef STRUCT_PARAMS_H +#define STRUCT_PARAMS_H + +#include + +inline params_t params; + +#endif From 532cd05a43022b27495c28294c2da59df714e0bb Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Mon, 1 Jul 2024 07:00:01 -0400 Subject: [PATCH 10/10] rename some Riemann runtime parameters (#2893) now they have "riemann" in their name, and also some drop "cg" so that it is clear that they can apply to other solvers --- .../ci-benchmarks/job_info_params.txt | 6 ++-- Exec/science/wdmerger/inputs | 2 +- .../tests/he_double_det/inputs_pakmor | 2 +- .../he_double_det/inputs_pakmor_simp_sdc | 2 +- .../tests/he_double_det/inputs_pakmor_strang | 2 +- .../tests/he_double_det/inputs_scaling | 2 +- .../wdmerger_collision/inputs_2d_collision | 2 +- .../tests/wdmerger_collision_1D/inputs | 2 +- .../wdmerger_retry/inputs_test_wdmerger_retry | 2 +- Source/driver/Castro.cpp | 8 +++--- Source/driver/_cpp_parameters | 15 +++++----- Source/hydro/riemann_2shock_solvers.H | 28 +++++++++---------- 12 files changed, 36 insertions(+), 37 deletions(-) diff --git a/Exec/science/flame_wave/ci-benchmarks/job_info_params.txt b/Exec/science/flame_wave/ci-benchmarks/job_info_params.txt index 127b5bcb64..252b98a482 100644 --- a/Exec/science/flame_wave/ci-benchmarks/job_info_params.txt +++ b/Exec/science/flame_wave/ci-benchmarks/job_info_params.txt @@ -54,9 +54,9 @@ castro.plm_limiter = 2 castro.hybrid_riemann = 0 castro.riemann_solver = 0 - castro.cg_maxiter = 12 - castro.cg_tol = 1e-05 - castro.cg_blend = 2 + castro.riemann_shock_maxiter = 12 + castro.riemann_pstar_tol = 1e-05 + castro.riemann_cg_blend = 2 castro.use_flattening = 1 castro.transverse_use_eos = 0 castro.transverse_reset_density = 1 diff --git a/Exec/science/wdmerger/inputs b/Exec/science/wdmerger/inputs index 5def509a0d..f9b22104b2 100644 --- a/Exec/science/wdmerger/inputs +++ b/Exec/science/wdmerger/inputs @@ -188,7 +188,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor index 2ff2712602..a7028d1544 100644 --- a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor +++ b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor @@ -200,7 +200,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc index db97ad3dac..6316753695 100644 --- a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc +++ b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc @@ -206,7 +206,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang index a57a028d0d..b56d9ba2a5 100644 --- a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang +++ b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang @@ -206,7 +206,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_scaling b/Exec/science/wdmerger/tests/he_double_det/inputs_scaling index 03c95c1636..990c81fd2a 100644 --- a/Exec/science/wdmerger/tests/he_double_det/inputs_scaling +++ b/Exec/science/wdmerger/tests/he_double_det/inputs_scaling @@ -203,7 +203,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/wdmerger_collision/inputs_2d_collision b/Exec/science/wdmerger/tests/wdmerger_collision/inputs_2d_collision index 8a2bca0293..60a7db22d6 100644 --- a/Exec/science/wdmerger/tests/wdmerger_collision/inputs_2d_collision +++ b/Exec/science/wdmerger/tests/wdmerger_collision/inputs_2d_collision @@ -153,7 +153,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/wdmerger_collision_1D/inputs b/Exec/science/wdmerger/tests/wdmerger_collision_1D/inputs index 6dcabf249c..9f40085682 100644 --- a/Exec/science/wdmerger/tests/wdmerger_collision_1D/inputs +++ b/Exec/science/wdmerger/tests/wdmerger_collision_1D/inputs @@ -142,7 +142,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Exec/science/wdmerger/tests/wdmerger_retry/inputs_test_wdmerger_retry b/Exec/science/wdmerger/tests/wdmerger_retry/inputs_test_wdmerger_retry index 0d63ab606c..a6eef171b3 100644 --- a/Exec/science/wdmerger/tests/wdmerger_retry/inputs_test_wdmerger_retry +++ b/Exec/science/wdmerger/tests/wdmerger_retry/inputs_test_wdmerger_retry @@ -133,7 +133,7 @@ castro.riemann_solver = 0 # For the CG Riemann solver, we need to tell the solver not to quit when # the iterations don't converge, but instead to do additional bisection iteration. -castro.cg_blend = 2 +castro.riemann_cg_blend = 2 # Use a lagged predictor estimate of the source terms in the hydro castro.source_term_predictor = 1 diff --git a/Source/driver/Castro.cpp b/Source/driver/Castro.cpp index d79f4d4935..f478f8402d 100644 --- a/Source/driver/Castro.cpp +++ b/Source/driver/Castro.cpp @@ -392,12 +392,12 @@ Castro::read_params () #endif if (riemann_solver == 1) { - if (cg_maxiter > HISTORY_SIZE) { - amrex::Error("error in riemanncg: cg_maxiter > HISTORY_SIZE"); + if (riemann_shock_maxiter > HISTORY_SIZE) { + amrex::Error("riemann_shock_maxiter > HISTORY_SIZE"); } - if (cg_blend == 2 && cg_maxiter < 5) { - amrex::Error("Error: need cg_maxiter >= 5 to do a bisection search on secant iteration failure."); + if (riemann_cg_blend == 2 && riemann_shock_maxiter < 5) { + amrex::Error("Error: need riemann_shock_maxiter >= 5 to do a bisection search on secant iteration failure."); } } #endif diff --git a/Source/driver/_cpp_parameters b/Source/driver/_cpp_parameters index a342478def..f5cdc29bd3 100644 --- a/Source/driver/_cpp_parameters +++ b/Source/driver/_cpp_parameters @@ -114,20 +114,19 @@ hybrid_riemann bool 0 # 2: HLLC riemann_solver int 0 -# for the Colella \& Glaz Riemann solver, the maximum number -# of iterations to take when solving for the star state -cg_maxiter int 12 +# maximum number of iterations to used in the Riemann solver +# when solving for the star state +riemann_shock_maxiter int 12 -# for the Colella \& Glaz Riemann solver, the tolerance to -# demand in finding the star state -cg_tol Real 1.0e-5 +# tolerance to use when finding the star stat +riemann_pstar_tol Real 1.0e-5 # for the Colella \& Glaz Riemann solver, what to do if # we do not converge to a solution for the star state. # 0 = do nothing; print iterations and exit # 1 = revert to the original guess for p-star -# 2 = do a bisection search for another 2 * cg_maxiter iterations. -cg_blend int 2 +# 2 = do a bisection search for another 2 * riemann_shock_maxiter iterations. +riemann_cg_blend int 2 # flatten the reconstructed profiles around shocks to prevent them # from becoming too thin diff --git a/Source/hydro/riemann_2shock_solvers.H b/Source/hydro/riemann_2shock_solvers.H index 208e9878e0..f81b495f29 100644 --- a/Source/hydro/riemann_2shock_solvers.H +++ b/Source/hydro/riemann_2shock_solvers.H @@ -58,7 +58,7 @@ namespace TwoShock { const amrex::Real ur, const amrex::Real pr, const amrex::Real taur, const amrex::Real gamer, const amrex::Real clsqr, const amrex::Real gdot, const amrex::Real gmin, const amrex::Real gmax, - const int lcg_maxiter, const amrex::Real lcg_tol, + const int lriemann_shock_maxiter, const amrex::Real lriemann_pstar_tol, amrex::Real& pstar, amrex::Real& gamstar, bool& converged, amrex::GpuArray& pstar_hist_extra) { @@ -106,7 +106,7 @@ namespace TwoShock { converged = false; amrex::Real pstar_c = 0.0; - for (int iter = 0; iter < PSTAR_BISECT_FACTOR*lcg_maxiter; iter++) { + for (int iter = 0; iter < PSTAR_BISECT_FACTOR * lriemann_shock_maxiter; iter++) { pstar_c = 0.5_rt * (pstar_lo + pstar_hi); pstar_hist_extra[iter] = pstar_c; @@ -125,7 +125,7 @@ namespace TwoShock { amrex::Real f_c = ustar_l - ustar_r; - if ( 0.5_rt * std::abs(pstar_lo - pstar_hi) < lcg_tol * pstar_c ) { + if ( 0.5_rt * std::abs(pstar_lo - pstar_hi) < lriemann_pstar_tol * pstar_c ) { converged = true; break; } @@ -239,7 +239,7 @@ namespace TwoShock { bool converged = false; int iter = 0; - while ((iter < cg_maxiter && !converged) || iter < 2) { + while ((iter < riemann_shock_maxiter && !converged) || iter < 2) { wsqge(ql.p, taul, gamel, gdot, gamstar, gmin, gmax, clsql, pstar, wlsq); @@ -281,7 +281,7 @@ namespace TwoShock { pstar = amrex::max(pstar, small_pres); amrex::Real err = std::abs(pstar - pstar_old); - if (err < cg_tol*pstar) { + if (err < riemann_pstar_tol * pstar) { converged = true; } @@ -300,11 +300,11 @@ namespace TwoShock { if (!converged) { - if (cg_blend == 0) { + if (riemann_cg_blend == 0) { #ifndef AMREX_USE_GPU std::cout << "pstar history: " << std::endl; - for (int iter_l=0; iter_l < cg_maxiter; iter_l++) { + for (int iter_l=0; iter_l < riemann_shock_maxiter; iter_l++) { std::cout << iter_l << " " << pstar_hist[iter_l] << std::endl; } @@ -316,11 +316,11 @@ namespace TwoShock { amrex::Error("ERROR: non-convergence in the Riemann solver"); #endif - } else if (cg_blend == 1) { + } else if (riemann_cg_blend == 1) { pstar = ql.p + ( (qr.p - ql.p) - wr * (qr.un - ql.un) ) * wl / (wl + wr); - } else if (cg_blend == 2) { + } else if (riemann_cg_blend == 2) { // we don't store the history if we are in CUDA, so // we can't do this @@ -328,7 +328,7 @@ namespace TwoShock { // first try to find a reasonable bounds amrex::Real pstarl = 1.e200; amrex::Real pstaru = -1.e200; - for (int n = cg_maxiter-6; n < cg_maxiter; n++) { + for (int n = riemann_shock_maxiter-6; n < riemann_shock_maxiter; n++) { pstarl = amrex::min(pstarl, pstar_hist[n]); pstaru = amrex::max(pstaru, pstar_hist[n]); } @@ -342,17 +342,17 @@ namespace TwoShock { ql.un, ql.p, taul, gamel, clsql, qr.un, qr.p, taur, gamer, clsqr, gdot, gmin, gmax, - cg_maxiter, cg_tol, + riemann_shock_maxiter, riemann_pstar_tol, pstar, gamstar, converged, pstar_hist_extra); if (!converged) { std::cout << "pstar history: " << std::endl; - for (int iter_l = 0; iter_l < cg_maxiter; iter_l++) { + for (int iter_l = 0; iter_l < riemann_shock_maxiter; iter_l++) { std::cout << iter_l << " " << pstar_hist[iter_l] << std::endl; } std::cout << "pstar extra history: " << std::endl; - for (int iter_l = 0; iter_l < PSTAR_BISECT_FACTOR*cg_maxiter; iter_l++) { + for (int iter_l = 0; iter_l < PSTAR_BISECT_FACTOR * riemann_shock_maxiter; iter_l++) { std::cout << iter_l << " " << pstar_hist_extra[iter_l] << std::endl; } @@ -368,7 +368,7 @@ namespace TwoShock { } else { #ifndef AMREX_USE_GPU - amrex::Error("ERROR: unrecognized cg_blend option."); + amrex::Error("ERROR: unrecognized riemann_cg_blend option."); #endif }