diff --git a/ExampleCodes/SUNDIALS/Exec/inputs_forward_euler b/ExampleCodes/SUNDIALS/Exec/inputs_forward_euler deleted file mode 100644 index 39b09830..00000000 --- a/ExampleCodes/SUNDIALS/Exec/inputs_forward_euler +++ /dev/null @@ -1,14 +0,0 @@ -n_cell = 32 -max_grid_size = 16 - -nsteps = 1000 -plot_int = 100 - -dt = 1.e-5 - -# INTEGRATION -## integration.type can take on the following values: -## 0 or "ForwardEuler" => Native AMReX Forward Euler integrator -## 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type -## 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.strategy -integration.type = ForwardEuler diff --git a/ExampleCodes/SUNDIALS/Exec/inputs_rk3 b/ExampleCodes/SUNDIALS/Exec/inputs_rk3 deleted file mode 100644 index 6eb60137..00000000 --- a/ExampleCodes/SUNDIALS/Exec/inputs_rk3 +++ /dev/null @@ -1,24 +0,0 @@ -n_cell = 32 -max_grid_size = 16 - -nsteps = 1000 -plot_int = 100 - -dt = 1.e-5 - -# INTEGRATION -## integration.type can take on the following values: -## 0 or "ForwardEuler" => Native AMReX Forward Euler integrator -## 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type -## 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.strategy -integration.type = RungeKutta - -## Native AMReX Explicit Runge-Kutta parameters -# -## integration.rk.type can take the following values: -### 0 = User-specified Butcher Tableau -### 1 = Forward Euler -### 2 = Trapezoid Method -### 3 = SSPRK3 Method -### 4 = RK4 Method -integration.rk.type = 3 diff --git a/ExampleCodes/SUNDIALS/Exec/inputs_sundials b/ExampleCodes/SUNDIALS/Exec/inputs_sundials deleted file mode 100644 index bbdf2c69..00000000 --- a/ExampleCodes/SUNDIALS/Exec/inputs_sundials +++ /dev/null @@ -1,44 +0,0 @@ -n_cell = 32 -max_grid_size = 16 - -nsteps = 1000 -plot_int = 100 - -dt = 1.e-5 - -# Use adaptive time stepping 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 -# 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type -# 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.type -# -# If using the SUNDIALS Submodule, then compile with USE_SUNDIALS=TRUE or -# AMReX_SUNDIALS=ON -integration.type = SUNDIALS - -# Set the SUNDIALS method type: -# ERK = Explicit Runge-Kutta method -# DIRK = Diagonally Implicit Runge-Kutta method -# IMEX-RK = Implicit-Explicit Additive Runge-Kutta method -# EX-MRI = Explicit Multirate Infatesimal method -# IM-MRI = Implicit Multirate Infatesimal method -# IMEX-MRI = Implicit-Explicit Multirate Infatesimal method -# -# Optionally select a specific SUNDIALS method by name, see the SUNDIALS -# documentation for the supported method names - -# Use the SUNDIALS default method for the chosen type (fixed or adaptive step sizes) -integration.sundials.type = ERK - -# Use forward Euler (fixed step sizes only) -# integration.sundials.type = ERK -# integration.sundials.method = ARKODE_FORWARD_EULER_1_1 - -# Use backward Euler (fixed step sizes only) -# integration.sundials.type = DIRK -# integration.sundials.method = ARKODE_BACKWARD_EULER_1_1 diff --git a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri index e380326c..92d09735 100644 --- a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri +++ b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Exec/inputs_sundials_mri @@ -12,13 +12,12 @@ max_grid_size = 16 nsteps = 12 # 60 # 600 plot_int = 2 # 10 # 100 -dt = 8.5e-3 # 1.7e-3 #1.7e-4 - # With integration.sundials.type = ERK, dt = 1.6e4 is stable, dt = 1.7e-4 is unstable # With integration.sundials.type = EX-MRI, dt = 1.7e-4 with fast_dt_ratio = 0.1 is stable # With integration.sundials.type = EX-MRI, dt = 1.7e-3 with fast_dt_ratio = 0.1 is stable # With integration.sundials.type = EX-MRI, dt = 8.5e-3 with fast_dt_ratio = 0.1 FAILS # With integration.sundials.type = EX-MRI, dt = 8.5e-3 with fast_dt_ratio = 0.02 is stable +dt = 8.5e-3 # To replicate heat equation diffusion_coef = 1.0 @@ -29,19 +28,18 @@ 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 # INTEGRATION -# integration.type can take on the following values: -# 0 or "ForwardEuler" => Native AMReX Forward Euler integrator -# 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type -# 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.type -# -# If using the SUNDIALS Submodule, then compile with USE_SUNDIALS=TRUE or -# AMReX_SUNDIALS=ON +## integration.type can take on the following values: +## 0 or "ForwardEuler" => Native AMReX Forward Euler integrator +## 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type +## 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.type +#integration.type = ForwardEuler +#integration.type = RungeKutta integration.type = SUNDIALS # Set the SUNDIALS method type: 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/Single-Rate/Exec/inputs b/ExampleCodes/SUNDIALS/Single-Rate/Exec/inputs new file mode 100644 index 00000000..02631eb7 --- /dev/null +++ b/ExampleCodes/SUNDIALS/Single-Rate/Exec/inputs @@ -0,0 +1,53 @@ +n_cell = 32 +max_grid_size = 16 + +nsteps = 5 +plot_int = 5 +dt = 1.e-4 + +#nsteps = 10 +#plot_int = 10 +#dt = 5.e-5 + +#nsteps = 20 +#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 +## 1 or "RungeKutta" => Native AMReX Explicit Runge Kutta controlled by integration.rk.type +## 2 or "SUNDIALS" => SUNDIALS backend controlled by integration.sundials.type +#integration.type = ForwardEuler +#integration.type = RungeKutta +integration.type = SUNDIALS + +## Native AMReX Explicit Runge-Kutta parameters +# +## integration.rk.type can take the following values: +### 0 = User-specified Butcher Tableau +### 1 = Forward Euler +### 2 = Trapezoid Method +### 3 = SSPRK3 Method +### 4 = RK4 Method +integration.rk.type = 1 + +# Set the SUNDIALS method type: +# ERK = Explicit Runge-Kutta method +# DIRK = Diagonally Implicit Runge-Kutta method +# +# Optionally select a specific SUNDIALS method by name, see the SUNDIALS +# documentation for the supported method names + +# Use forward Euler (fixed step sizes only) +integration.sundials.type = ERK +integration.sundials.method = ARKODE_FORWARD_EULER_1_1 + +# Use backward Euler (fixed step sizes only) +#integration.sundials.type = DIRK +#integration.sundials.method = ARKODE_BACKWARD_EULER_1_1 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 88% rename from ExampleCodes/SUNDIALS/Source/main.cpp rename to ExampleCodes/SUNDIALS/Single-Rate/Source/main.cpp index ba5b5b7b..ce6b221c 100644 --- a/ExampleCodes/SUNDIALS/Source/main.cpp +++ b/ExampleCodes/SUNDIALS/Single-Rate/Source/main.cpp @@ -123,7 +123,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 @@ -162,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()); @@ -171,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(); @@ -181,7 +181,7 @@ void main_main () // fill the right-hand-side for phi amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) { - phi_rhs_array(i,j,k) = ( (phi_array(i+1,j,k) - 2.*phi_array(i,j,k) + phi_array(i-1,j,k)) / (dx[0]*dx[0]) + phi_rhs_array(i,j,k) = ( (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]) @@ -200,16 +200,23 @@ void main_main () integrator.set_time_step(dt); } + Real evolution_start_time = ParallelDescriptor::second(); + for (int step = 1; step <= nsteps; ++step) { // Set time to evolve to time += dt; + Real step_start_time = ParallelDescriptor::second(); + // Advance to output time integrator.evolve(phi, time); + Real step_stop_time = ParallelDescriptor::second() - step_start_time; + ParallelDescriptor::ReduceRealMax(step_stop_time); + // Tell the I/O Processor to write out which step we're doing - amrex::Print() << "Advanced step " << step << "\n"; + amrex::Print() << "Advanced step " << step << " in " << step_stop_time << " seconds; dt = " << dt << " time = " << time << "\n"; // Write a plotfile of the current data (plot_int was defined in the inputs file) if (plot_int > 0 && step%plot_int == 0) @@ -218,4 +225,8 @@ void main_main () WriteSingleLevelPlotfile(pltfile, phi, {"phi"}, geom, time, step); } } + + Real evolution_stop_time = ParallelDescriptor::second() - evolution_start_time; + ParallelDescriptor::ReduceRealMax(evolution_stop_time); + amrex::Print() << "Total evolution time = " << evolution_stop_time << " seconds\n"; } 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