From ffa94a0e83e283390ad3f8faa447d222a3539ce4 Mon Sep 17 00:00:00 2001
From: prkkumar <prabhatams@gmail.com>
Date: Wed, 17 Jul 2024 09:45:32 -0700
Subject: [PATCH] fast rhs fun

---
 .../Reaction-Diffusion/Source/main.cpp        | 38 +++++++++++++++++++
 1 file changed, 38 insertions(+)

diff --git a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp
index 4eab09d4..4add8ebd 100644
--- a/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp
+++ b/ExampleCodes/SUNDIALS/Reaction-Diffusion/Source/main.cpp
@@ -202,9 +202,43 @@ void main_main ()
         }
     };
 
+    auto rhs_fast_function = [&](MultiFab& S_rhs, const MultiFab& S_data,
+                            const Real /* time */) {
+
+        // loop over boxes
+        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();
+
+            const Array4<const Real>& phi_array = phi_data.array(mfi);
+            const Array4<Real>& phi_rhs_array = phi_rhs.array(mfi);
+
+            // 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) = 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]) );
+#endif
+            });
+        }
+    };
+
     TimeIntegrator<MultiFab> integrator(phi, time);
     integrator.set_pre_rhs_action(pre_rhs_function);
     integrator.set_rhs(rhs_function);
+    integrator.set_fast_rhs(rhs_fast_function);
     if (adapt_dt) {
         integrator.set_adaptive_step();
         integrator.set_tolerances(reltol, abstol);
@@ -212,6 +246,10 @@ void main_main ()
         integrator.set_time_step(dt);
     }
 
+    // This sets the ratio of slow timestep size to fast timestep size as an integer,
+    // or equivalently, the number of fast timesteps per slow timestep.
+    integrator.set_slow_fast_timestep_ratio(10);
+
     for (int step = 1; step <= nsteps; ++step)
     {
         // Set time to evolve to