From 064eb4e3c594b57a2c5a35d4ac119ed79bb6eefd Mon Sep 17 00:00:00 2001
From: Yohai Meiron <yohai@meiron.org>
Date: Mon, 23 Aug 2021 19:40:55 -0400
Subject: [PATCH] Removed init. of pot_act_ext from main loop; zero filling of
 pot_act_ext now only if external potential active; improved warning on
 minimum timestep; changed CUDAHOME to CUDA_HOME

---
 Makefile     | 20 ++++++++++----------
 phigrape.cpp | 17 ++++++++---------
 2 files changed, 18 insertions(+), 19 deletions(-)

diff --git a/Makefile b/Makefile
index f9d07cb..3ac4642 100644
--- a/Makefile
+++ b/Makefile
@@ -1,18 +1,18 @@
-CUDAHOME     ?= /usr/local/cuda
+CUDA_HOME    ?= /usr/local/cuda
 CPPFLAGS     += -DETICS
 OPTIMIZATION ?= 3
 
-CUDAINC       = -I$(CUDAHOME)/include -I$(CUDAHOME)/samples/common/inc/
-CUDALIB       = -L$(CUDAHOME)/lib64 -lcudart -lcudadevrt -lcuda
+CUDAINC       = -I$(CUDA_HOME)/include -I$(CUDA_HOME)/samples/common/inc/
+CUDALIB       = -L$(CUDA_HOME)/lib64 -lcudart -lcudadevrt -lcuda
 
 default: grapite
-grapite: GRAPEHOME = ../grapite
-grapite: GRAPELIB  = -L$(GRAPEHOME) -lgrapite
-yebisu:  GRAPEHOME = ../yebisu
-yebisu:  GRAPELIB  = -L$(GRAPEHOME) -lyebisug6
-sapporo: GRAPEHOME = ../sapporo2/lib
-sapporo: GRAPELIB  = -L$(GRAPEHOME) -lsapporo
-         GRAPEINC  = -I$(GRAPEHOME)
+grapite: GRAPE_HOME = ../grapite
+grapite: GRAPELIB   = -L$(GRAPE_HOME) -lgrapite
+yebisu:  GRAPE_HOME = ../yebisu
+yebisu:  GRAPELIB   = -L$(GRAPE_HOME) -lyebisug6
+sapporo: GRAPE_HOME = ../sapporo2/lib
+sapporo: GRAPELIB   = -L$(GRAPE_HOME) -lsapporo
+         GRAPEINC   = -I$(GRAPE_HOME)
 
 CXXFLAGS     += -std=c++11 -O$(OPTIMIZATION)
 INC           = $(GRAPEINC) $(CUDAINC)
diff --git a/phigrape.cpp b/phigrape.cpp
index b82cb12..c484045 100644
--- a/phigrape.cpp
+++ b/phigrape.cpp
@@ -71,12 +71,7 @@ public:
     void add_component(External_gravity &component)
     {
         components.push_back(&component);
-    }
-    bool any_active()
-    {
-        for (auto component : components)
-            if (component->is_active) return true;
-        return false;
+        if (component.is_active) any_active = true;
     }
     void operator()(int n, const std::vector<double3> &x, const std::vector<double3> &v, std::vector<double> &pot, std::vector<double3> &acc, std::vector<double3> &jrk)
     {
@@ -92,6 +87,7 @@ public:
         }
         fflush(stdout);
     }
+    bool any_active = false;
 private:
     std::vector<External_gravity*> components;
 };
@@ -497,6 +493,7 @@ int main(int argc, char *argv[])
     std::vector<int> ind_act(N);
     std::vector<double3> x_act_new(N), v_act_new(N), a_act_new(N), adot_act_new(N);
     std::vector<double> t(N, time_cur), pot_act_new(N);
+    std::vector<double> pot_act_ext(N, 0.);
 
     // Functors for the main integration loop
     Active_search active_search(myRank, n_proc, n_loc, N, grapite_active_search_flag);
@@ -555,8 +552,10 @@ int main(int argc, char *argv[])
         }
 
         /* Calculate gravity on active particles due to external forces */
-        std::vector<double> pot_act_ext(N, 0.);
-        calc_ext_grav(n_act, x_act_new, v_act_new, pot_act_ext, a_act_new, adot_act_new);
+        if (calc_ext_grav.any_active) {
+            std::fill_n(begin(pot_act_ext), n_act, 0);
+            calc_ext_grav(n_act, x_act_new, v_act_new, pot_act_ext, a_act_new, adot_act_new);
+        }
 
         /* correct the active particles positions etc... on all the nodes */
         double min_dt = dt_max;
@@ -589,7 +588,7 @@ int main(int argc, char *argv[])
 
             if (config.dt_min_warning && (myRank == 0)) {
                 if (dt_tmp == dt_min) {
-                    printf("!!! Warning1: dt_act = dt_min = %.6E \t ind_act = %07d \n", dt_tmp, ind_act[i]);
+                    printf("!!! Warning1: dt_act = dt_min = %.6E \t ind_act = %07d time_cur=%.16E\n", dt_tmp, ind_act[i], time_cur);
                     fflush(stdout);
                 }
             }