From 1215a4225b2caa62a904f57a3d03212f868a60c2 Mon Sep 17 00:00:00 2001 From: Andy Nonaka Date: Thu, 5 Sep 2024 11:42:28 -0700 Subject: [PATCH] cleanup to both SUNDIALS tutorials --- .../Exec/inputs_sundials_mri | 2 +- .../Reaction-Diffusion/Source/main.cpp | 35 +++++++------------ .../SUNDIALS/{ => Single-Rate}/CMakeLists.txt | 0 .../{ => Single-Rate}/Exec/GNUmakefile | 6 ++-- .../{ => Single-Rate}/Exec/README_sundials | 0 .../SUNDIALS/{ => Single-Rate}/Exec/inputs | 5 +++ .../{ => Single-Rate}/Source/Make.package | 0 .../{ => Single-Rate}/Source/main.cpp | 9 ++--- .../{ => Single-Rate}/Source/myfunc.H | 0 9 files changed, 25 insertions(+), 32 deletions(-) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/CMakeLists.txt (100%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Exec/GNUmakefile (88%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Exec/README_sundials (100%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Exec/inputs (88%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Source/Make.package (100%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Source/main.cpp (96%) rename ExampleCodes/SUNDIALS/{ => Single-Rate}/Source/myfunc.H (100%) diff --git a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri index e380326c..a5079043 100644 --- a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri +++ b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri @@ -29,7 +29,7 @@ use_MRI = true fast_dt_ratio = 0.02 -# Use adaptive time stepping and set integrator relative and absolute tolerances +# Use adaptive time stepping (multi-stage integrators only!) and set integrator relative and absolute tolerances # adapt_dt = true # reltol = 1.0e-4 # abstol = 1.0e-9 diff --git a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp index 32089415..0f78e831 100644 --- a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp +++ b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp @@ -48,6 +48,10 @@ void main_main () Real reltol = 1.0e-4; Real abstol = 1.0e-9; + // Reaction and Diffusion Coefficients + Real reaction_coef = 0.0; + Real diffusion_coef = 1.0; + // inputs parameters { // ParmParse is way of reading inputs from the inputs file @@ -85,6 +89,8 @@ void main_main () pp.query("reltol",reltol); pp.query("abstol",abstol); + pp.query("reaction_coef",reaction_coef); + pp.query("diffusion_coef",diffusion_coef); } // ********************************** @@ -132,7 +138,7 @@ void main_main () // How Boxes are distrubuted among MPI processes DistributionMapping dm(ba); - // we allocate two phi multifabs; one will store the old state, the other the new. + // allocate phi MultiFab MultiFab phi(ba, dm, Ncomp, Nghost); // time = starting time in the simulation @@ -171,8 +177,7 @@ void main_main () WriteSingleLevelPlotfile(pltfile, phi, {"phi"}, geom, time, 0); } - auto rhs_fast_function = [&](MultiFab& S_rhs, MultiFab& S_data, - const Real /* time */) { + auto rhs_fast_function = [&](MultiFab& S_rhs, MultiFab& S_data, const Real /* time */) { // fill periodic ghost cells S_data.FillBoundary(geom.periodicity()); @@ -181,12 +186,6 @@ void main_main () auto& phi_data = S_data; auto& phi_rhs = S_rhs; - //Reaction and Diffusion Coefficients - Real reaction_coef = 0.0; - Real diffusion_coef = 1.0; - ParmParse pp; - pp.query("diffusion_coef",diffusion_coef); - for ( MFIter mfi(phi_data); mfi.isValid(); ++mfi ) { const Box& bx = mfi.validbox(); @@ -200,14 +199,14 @@ void main_main () phi_rhs_array(i,j,k) = diffusion_coef*( (phi_array(i+1,j,k) - 2.*phi_array(i,j,k) + phi_array(i-1,j,k)) / (dx[0]*dx[0]) +(phi_array(i,j+1,k) - 2.*phi_array(i,j,k) + phi_array(i,j-1,k)) / (dx[1]*dx[1]) #if (AMREX_SPACEDIM == 3) - +(phi_array(i,j,k+1) - 2.*phi_array(i,j,k) + phi_array(i,j,k-1)) / (dx[2]*dx[2]) ); + +(phi_array(i,j,k+1) - 2.*phi_array(i,j,k) + phi_array(i,j,k-1)) / (dx[2]*dx[2]) #endif + ); }); } }; - auto rhs_function = [&](MultiFab& S_rhs, MultiFab& S_data, - const Real /* time */) { + auto rhs_function = [&](MultiFab& S_rhs, MultiFab& S_data, const Real /* time */) { // fill periodic ghost cells S_data.FillBoundary(geom.periodicity()); @@ -216,13 +215,6 @@ void main_main () auto& phi_data = S_data; auto& phi_rhs = S_rhs; - //Reaction and Diffusion Coefficients - Real reaction_coef = 0.0; - Real diffusion_coef = 1.0; - ParmParse pp; - pp.query("reaction_coef",reaction_coef); - pp.query("diffusion_coef",diffusion_coef); - for ( MFIter mfi(phi_data); mfi.isValid(); ++mfi ) { const Box& bx = mfi.validbox(); @@ -239,15 +231,14 @@ void main_main () phi_rhs_array(i,j,k) = diffusion_coef*( (phi_array(i+1,j,k) - 2.*phi_array(i,j,k) + phi_array(i-1,j,k)) / (dx[0]*dx[0]) +(phi_array(i,j+1,k) - 2.*phi_array(i,j,k) + phi_array(i,j-1,k)) / (dx[1]*dx[1]) #if (AMREX_SPACEDIM == 3) - +(phi_array(i,j,k+1) - 2.*phi_array(i,j,k) + phi_array(i,j,k-1)) / (dx[2]*dx[2]) ) + +(phi_array(i,j,k+1) - 2.*phi_array(i,j,k) + phi_array(i,j,k-1)) / (dx[2]*dx[2]) #endif - -1.*reaction_coef*phi_array(i,j,k); + ) - 1.*reaction_coef*phi_array(i,j,k); } }); } }; - TimeIntegrator integrator(phi, time); integrator.set_rhs(rhs_function); if (use_MRI) { diff --git a/ExampleCodes/SUNDIALS/CMakeLists.txt b/ExampleCodes/SUNDIALS/Single-Rate/CMakeLists.txt similarity index 100% rename from ExampleCodes/SUNDIALS/CMakeLists.txt rename to ExampleCodes/SUNDIALS/Single-Rate/CMakeLists.txt diff --git a/ExampleCodes/SUNDIALS/Exec/GNUmakefile b/ExampleCodes/SUNDIALS/Single-Rate/Exec/GNUmakefile similarity index 88% rename from ExampleCodes/SUNDIALS/Exec/GNUmakefile rename to ExampleCodes/SUNDIALS/Single-Rate/Exec/GNUmakefile index 448b77a0..7947b0b7 100644 --- a/ExampleCodes/SUNDIALS/Exec/GNUmakefile +++ b/ExampleCodes/SUNDIALS/Single-Rate/Exec/GNUmakefile @@ -1,5 +1,5 @@ # AMREX_HOME defines the directory in which we will find all the AMReX code. -AMREX_HOME ?= ../../../../amrex +AMREX_HOME ?= ../../../../../amrex DEBUG = FALSE USE_MPI = TRUE @@ -17,9 +17,9 @@ INCLUDE_LOCATIONS += ../Source ifeq ($(USE_SUNDIALS),TRUE) ifeq ($(USE_CUDA),TRUE) - SUNDIALS_ROOT ?= $(TOP)../../../../sundials/instdir_cuda + SUNDIALS_ROOT ?= $(TOP)../../../../../sundials/instdir_cuda else - SUNDIALS_ROOT ?= $(TOP)../../../../sundials/instdir + SUNDIALS_ROOT ?= $(TOP)../../../../../sundials/instdir endif ifeq ($(NERSC_HOST),perlmutter) SUNDIALS_LIB_DIR ?= $(SUNDIALS_ROOT)/lib64 diff --git a/ExampleCodes/SUNDIALS/Exec/README_sundials b/ExampleCodes/SUNDIALS/Single-Rate/Exec/README_sundials similarity index 100% rename from ExampleCodes/SUNDIALS/Exec/README_sundials rename to ExampleCodes/SUNDIALS/Single-Rate/Exec/README_sundials diff --git a/ExampleCodes/SUNDIALS/Exec/inputs b/ExampleCodes/SUNDIALS/Single-Rate/Exec/inputs similarity index 88% rename from ExampleCodes/SUNDIALS/Exec/inputs rename to ExampleCodes/SUNDIALS/Single-Rate/Exec/inputs index be6b0821..2d2d32d0 100644 --- a/ExampleCodes/SUNDIALS/Exec/inputs +++ b/ExampleCodes/SUNDIALS/Single-Rate/Exec/inputs @@ -13,6 +13,11 @@ dt = 1.e-4 #plot_int = 20 #dt = 2.5e-5 +# Use adaptive time stepping (multi-stage integrators only!) and set integrator relative and absolute tolerances +# adapt_dt = true +# reltol = 1.0e-4 +# abstol = 1.0e-9 + # INTEGRATION ## integration.type can take on the following values: ## 0 or "ForwardEuler" => Native AMReX Forward Euler integrator diff --git a/ExampleCodes/SUNDIALS/Source/Make.package b/ExampleCodes/SUNDIALS/Single-Rate/Source/Make.package similarity index 100% rename from ExampleCodes/SUNDIALS/Source/Make.package rename to ExampleCodes/SUNDIALS/Single-Rate/Source/Make.package diff --git a/ExampleCodes/SUNDIALS/Source/main.cpp b/ExampleCodes/SUNDIALS/Single-Rate/Source/main.cpp similarity index 96% rename from ExampleCodes/SUNDIALS/Source/main.cpp rename to ExampleCodes/SUNDIALS/Single-Rate/Source/main.cpp index 0db29934..ce6b221c 100644 --- a/ExampleCodes/SUNDIALS/Source/main.cpp +++ b/ExampleCodes/SUNDIALS/Single-Rate/Source/main.cpp @@ -38,7 +38,7 @@ void main_main () Real dt; // use adaptive time step (dt used to set output times) - int adapt_dt = 0; + bool adapt_dt = false; // adaptive time step relative and absolute tolerances Real reltol = 1.0e-4; @@ -129,9 +129,6 @@ void main_main () // time = starting time in the simulation Real time = 0.0; - Print() << "dt = " << dt << std::endl; - Print() << "Diffusive CFL time step = " << dx[0]*dx[0]/(2.*AMREX_SPACEDIM) << std::endl; - // ********************************** // INITIALIZE DATA @@ -165,8 +162,7 @@ void main_main () WriteSingleLevelPlotfile(pltfile, phi, {"phi"}, geom, time, 0); } - auto rhs_function = [&](MultiFab& S_rhs, MultiFab& S_data, - const Real /* time */) { + auto rhs_function = [&](MultiFab& S_rhs, MultiFab& S_data, const Real /* time */) { // fill periodic ghost cells S_data.FillBoundary(geom.periodicity()); @@ -174,6 +170,7 @@ void main_main () // loop over boxes auto& phi_data = S_data; auto& phi_rhs = S_rhs; + for ( MFIter mfi(phi_data); mfi.isValid(); ++mfi ) { const Box& bx = mfi.validbox(); diff --git a/ExampleCodes/SUNDIALS/Source/myfunc.H b/ExampleCodes/SUNDIALS/Single-Rate/Source/myfunc.H similarity index 100% rename from ExampleCodes/SUNDIALS/Source/myfunc.H rename to ExampleCodes/SUNDIALS/Single-Rate/Source/myfunc.H