Added a Model class supporting checkpoint/restart

Added an EcolabGraph class supporting  forAll method for parallel execution over cells.
Initial implementation of SYCL version of Ecolab model.

include/Makefile

@@ -21,6 +21,7 @@ TIMER=
 # use TkPhotoSurface for CairoItem, overriding platform specific optimisation
 TKPHOTOSURFACE=
 NOGUI=
+DPCPP=
 
 # enables ecolab for the MXE environment 
 MXE=
@@ -103,9 +104,15 @@ CPP=$(CPLUSPLUS) -E
 # disable inclusion of mpi++.h in the MPICH case (don't know what the problem is here)
 FLAGS+=-UHAVE_MPI_CPP
 else
-CC=gcc
+# Engage Intel OneAPI compiler
+ifdef DPCPP
+CPLUSPLUS=icpx
+SYCL+=-fsycl
+else
 CPLUSPLUS=g++
-LINK=$(CPLUSPLUS)
+endif
+CC=gcc
+LINK=$(CPLUSPLUS) $(SYCL)
 CPP=g++ -E
 endif
 
@@ -622,7 +629,7 @@ DCXX=distcc
 endif
 
 .cc.o: 
-	$(DCXX) -c $(FLAGS) $(CXXFLAGS) $(OPT) -o $@ $<
+	$(DCXX) -c $(FLAGS) $(SYCL) $(CXXFLAGS) $(OPT) -o $@ $<
 
 .c.o: 
 	$(CC) -c $(FLAGS) $(OPT) -o $@ $<

include/arrays.h

@@ -1566,7 +1566,7 @@ namespace ecolab
             array_data<T>* oldData=dt;
             bool freeMem = ref()-- == 0;
             dt=alloc(size()); 
-            memcpy(data(),oldData->dt,size()*sizeof(T));
+            memcpy(dt->dt,oldData->dt,size()*sizeof(T));
             if (freeMem) free(oldData);
           }
       }
@@ -2033,7 +2033,7 @@ namespace ecolab
     {
       array<double> r(x.size());
       vdCosh(static_cast<int>(x.size()),x.data(),r.data());
-      return r;
+ r;
     }
 
     inline array<float> sinh(const array<float>& x)

include/cairo_base.h

@@ -10,6 +10,8 @@
 #define CAIRO_BASE_H
 #if defined(CAIRO)
 #include "classdesc.h"
+#undef None
+#undef Success
 #include "arrays.h"
 
 #include <cairo/cairo.h>

include/ecolab.h

@@ -13,6 +13,7 @@
 #define ECOLAB_H
 #include <stdlib.h>
 #include "pythonBuffer.h"
+#include "graphcode.h"
 
 // mpi.h must appear before any standard library stuff
 #ifdef MPI_SUPPORT
@@ -52,6 +53,11 @@ typedef classdesc::string eco_string;
 #include <isa_base.h>
 
 #include "eco_strstream.h"
+#include "random.h"
+
+#ifdef SYCL_LANGUAGE_VERSION
+#include <sycl/sycl.hpp>
+#endif
 
 namespace ecolab
 {
@@ -72,9 +78,37 @@ namespace ecolab
   using classdesc::xdr_pack;
   using classdesc::isa_t;
   using classdesc::is_array;
-}
 
-#include "random.h"
+  template <class M> struct Model
+  {
+    /// checkpoint the model to \a fileName
+    void checkpoint(const char* fileName);
+    /// restart the model from checkpoint saved in \a fileName
+    void restart(const char* fileName);
+  };
+
+#ifdef SYCL_LANGUAGE_VERSION
+  sycl::queue& syclQ();
+#endif
+
+  template <class Cell> struct EcolabGraph: public graphcode::Graph<Cell>
+  {
+    /// apply a functional to all local cells of this processor in parallel
+    /// @param f 
+    template <class F>
+    void forAll(F f) {
+#ifdef SYCL_LANGUAGE_VERSION
+      syclQ().submit([&](sycl::handler& h) {
+        h.parallel_for(this->size(), [=,this](auto i) {
+          f(*(*this)[i]->template as<Cell>());
+        });
+      });
+#else
+      for (auto& i: *this) f(*i->template as<Cell>());
+#endif
+    }
+  };
+}
 
 #ifdef MPI_SUPPORT
 #include <classdescMP.h>
@@ -88,4 +122,5 @@ namespace ecolab
 }
 #endif
 
+#include "ecolab.cd"
 #endif  /* ECOLAB_H */

include/ecolab_epilogue.h

@@ -11,7 +11,28 @@
 
 Put this after including all .cd files
 */
-#include "classdesc_epilogue.h"
 #ifdef GRAPHCODE_H
 #include "graphcode.cd"
 #endif
+#include "classdesc_epilogue.h"
+
+#ifdef ECOLAB_H
+#include "pack_base.h"
+
+namespace ecolab
+{
+  template <class M>
+  void Model<M>::checkpoint(const char* fileName)
+  {
+    pack_t b(fileName,"w");
+    b<<static_cast<M&>(*this);
+  }
+
+  template <class M>
+  void Model<M>::restart(const char* fileName)
+  {
+    pack_t b(fileName,"r");
+    b>>static_cast<M&>(*this);
+  }
+}
+#endif

models/ecolab_model.cc

@@ -26,6 +26,7 @@ using array_ns::pcoord;
 
 #include "ecolab_epilogue.h"
 
+
 using namespace classdesc;
 
 // TODO - move this into main library
@@ -51,7 +52,6 @@ namespace
       r[i]=ROUND(x[i]);
     return r;
   }
-
 }
 
 
@@ -462,7 +462,7 @@ void SpatialModel::setGrid(size_t nx, size_t ny)
 void SpatialModel::generate(unsigned niter)
 {
   if (tstep==0) makeConsistent();
-  for (auto& i: *this) i->as<EcolabCell>()->generate(niter,*this);
+  forAll([=,this](EcolabCell& c) {c.generate(niter,*this);});
   tstep+=niter;
 }
 

models/ecolab_model.h

@@ -8,7 +8,7 @@
 
 #include <netcomplexity.h>
 #include <cairoSurfaceImage.h>
-#include <graphcode.h>
+#include <ecolab.h>
 using classdesc::Object;
 
 #include <vector>
@@ -67,7 +67,7 @@ struct EcolabPoint
   unsigned nsp() const; ///< number of living species in this cell
 };
 
-struct PanmicticModel: public ModelData, public EcolabPoint
+struct PanmicticModel: public ModelData, public EcolabPoint, public ecolab::Model<PanmicticModel>
 {
   ConnectionPlot connectionPlot;
   void updateConnectionPlot() {connectionPlot.update(density,interaction);}
@@ -81,9 +81,10 @@ struct PanmicticModel: public ModelData, public EcolabPoint
 
 struct EcolabCell: public EcolabPoint, public graphcode::Object<EcolabCell> {};
 
-class SpatialModel: public ModelData, public graphcode::Graph<EcolabCell>
+class SpatialModel: public ModelData, public EcolabGraph<EcolabCell>, public ecolab::Model<SpatialModel>
 {
   size_t numX=1, numY=1;
+  CLASSDESC_ACCESS(SpatialModel);
 public:
   size_t makeId(size_t x, size_t y) const {return x%numX + numX*(y%numY);}
   void setGrid(size_t nx, size_t ny);

src/ecolab.cc

@@ -8,19 +8,39 @@
 
 #include "tcl++.h"
 #define TK 1
+#undef None
 #include "cairoSurfaceImage.h"
+#undef None
 #include "plot.h"
 #include "pythonBuffer.h"
 #include "ecolab.h"
 #ifdef MPI_SUPPORT
 #include "graphcode.h"
 #endif
+#include "classdesc.h"
 #include "ecolab_epilogue.h"
 using namespace ecolab;
 using namespace std;
 
 #include <dlfcn.h>
 
+#ifdef SYCL_LANGUAGE_VERSION
+using namespace sycl;
+queue& ecolab::syclQ() {
+  static queue q{default_selector_v};
+  return q;
+}
+void* operator new(size_t s) {
+  cout<<s<<" bytes allocated shared"<<std::endl;
+  if (auto r=malloc_shared(s,syclQ()))
+    return r;
+  throw std::bad_alloc();
+}
+void operator delete(void* p) {return free(p,syclQ());}
+void* operator new[](size_t s) {return operator new(s);}
+void operator delete[](void* p) {return free(p,syclQ());}
+#endif  
+
 namespace ecolab
 {
   bool interpExiting=false;