Fixing conflicts

CHANGELOG.md

 # Change Log
 All notable changes to this project will be documented in this file.
 
-## [1.0.0] 13 September 2017
+## [1.1.0] February 2018
+
+### Added
+
+- Interface for Multi-vector dynamic load balancing
+- Increaded performance for grid ghost get
+- Introduced forms to increase the performance of the grid iterator in case of stencil code (see example 5_GrayScott)
+- EMatrix wrapped eigen matrices compatibles with vector_dist_id
+- General tuning for high dimension vector_dist_id (up to 50 dimensions) + PS_CMA_ES (Particle-Swarm Covariant Matrix Adaptation Evolution Strategy) example in Numerics
+- Added Discrete element Method example (8_DEM)
+- Added serial_to_parallel example VCluster (2_serial_to_parallel). The example it show how to port a serial example into openfpm gradually swtiching from
+  a serial section to a parallel section
+- Introduced map(LOCAL) for fast communication in case we have small movement
+
+### Fixed
+
+- Installation/detection of PETSC
+- CRITICAL-BUG scalar product in combination with vector product is broken (it return 0)
+- Fixing 2D IO in binary for vector
+- Fixing 1D grid writer in ASCII mode
+- Fixing Intel compilation of Linear algebra
+
+## [1.0.0] 13 September 2017 (Codename: Vortex)
 
 ### Added
 - Introduced getDomainIterator for Cell-list

build.sh

@@ -81,12 +81,8 @@ then
 
  source /etc/profile
  echo "$PATH"
- module load eigen/3.2.0
- module load suitesparse/4.2.1-gnu-multimkl
- module load boost/1.60.0
- module load gcc/5.3.0
- module load openmpi/1.10.2-gnu
- module unload bullxmpi
+ module load gcc/7.1.0
+ module load openmpi/3.0.0-gnu7.1
  
  export LD_LIBRARY_PATH="$LD_LIBRARY_PATH:/home/incard/PARMETIS/lib:/home/incard/METIS/lib:/home/incard/HDF5/lib"
 
@@ -105,20 +101,27 @@ else
  echo "Compiling general"
 
  source ~/.bashrc
+ 
+ installation_dir=""
+ if [ x"$2" == x"sbalzarini-mac-15" ]; then
+  installation_dir="--prefix=/Users/jenkins/openfpm_install"
+ else
+  installation_dir="--prefix=/home/jenkins/openfpm_install"
+ fi
 
  mkdir $HOME/$branch
  if [ x"$4" == x"full" ]; then
-  ./install -i $HOME/$branch  -s -c "--prefix=/Users/jenkins/openfpm_install"
+  ./install -i $HOME/$branch  -s -c "$installation_dir"
   mv $HOME/openfpm_vars $HOME/openfpm_vars_$branch
   source $HOME/openfpm_vars_$branch
  elif [ x"$3" == x"numerics" ]; then
   branch=$(git ls-remote --heads origin | grep $(git rev-parse HEAD) | cut -d / -f 3)
-  ./install -i $HOME/$branch  -m -s -c "--prefix=/home/jenkins/openfpm_install"
+  ./install -i $HOME/$branch  -m -s -c "$installation_dir"
   mv $HOME/openfpm_vars $HOME/openfpm_vars_$branch
   source $HOME/openfpm_vars_$branch
   make $3
  else
-  ./install -i $HOME/$branch -m -s -c "--prefix=/Users/jenkins/openfpm_install --no-recursion"
+  ./install -i $HOME/$branch -m -s -c "$installation_dir --no-recursion"
   mv $HOME/openfpm_vars $HOME/openfpm_vars_$branch
   source $HOME/openfpm_vars_$branch
   make $3

build_sec_OS.sh

@@ -18,6 +18,7 @@ if [ "$2" == "windows10" ]; then
 
     branch=$(git ls-remote --heads origin | grep $(git rev-parse HEAD) | cut -d / -f 3)    
 
+    git pull origin master
     ./install -i "/home/jenkins/$branch" < input_install
 
 fi

configure.ac

@@ -64,8 +64,7 @@ INCLUDES_PATH=" "
 echo "$base" > install_dir
 
 # Needed for build library
-AC_PROG_RANLIB
-AM_PROG_AR
+LT_INIT
 
 # Checks for programs.
 AC_PROG_CXX
@@ -139,14 +138,6 @@ fi
 AX_LIB_HILBERT([],[echo "Cannot detect libhilbert, use the --with-libhilbert option if it is not installed in the default location"
                     exit 210])
 
-##########
-
-## Check for PETSC
-
-AX_LIB_PETSC()
-
-## Check for quadmath
-
 
 have_quad_lib=no
 have_quad_head=no
@@ -245,6 +236,7 @@ AC_COMPILE_IFELSE([AC_LANG_PROGRAM([])],
 )
 AC_LANG_POP([C++])
 CXXFLAGS="$my_save_cflags"
+
 AC_SUBST([AM_CXXFLAGS])
 
 
@@ -377,6 +369,12 @@ AX_LAPACK([],[])
 
 AX_SUITESPARSE([],[])
 
+##########
+
+## Check for PETSC
+
+AX_LIB_PETSC()
+
 ###### Checking for EIGEN
 
 AX_EIGEN([],[])

dep_dir

+master

example/Debugging/1_classes/main.cpp

@@ -17,7 +17,6 @@
 #define SE_CLASS3
 #define THROW_ON_ERROR
 #include "Memleak_check.hpp"
-#include "data_type/scalar.hpp"
 #include "Grid/grid_dist_id.hpp"
 #include "Decomposition/CartDecomposition.hpp"
 #include "Point_test.hpp"
@@ -72,7 +71,7 @@ int main(int argc, char* argv[])
 	// * g: ghost extension
 	//
 	//
-	grid_dist_id<3, float, scalar<float[3]>, CartDecomposition<3,float>> * g_dist = new grid_dist_id<3, float, scalar<float[3]>, CartDecomposition<3,float>>(sz,domain,g);
+	grid_dist_id<3, float, aggregate<float[3]>, CartDecomposition<3,float>> * g_dist = new grid_dist_id<3, float, aggregate<float[3]> >(sz,domain,g);
 
 	//
 	// ### WIKI 6 ###

example/Grid/0_simple/main.cpp

@@ -9,6 +9,7 @@
  * \subpage Grid_2_solve_eq
  * \subpage Grid_3_gs
  * \subpage Grid_3_gs_3D
+ * \subpage Grid_3_gs_3D_vector
  *
  */
 

example/Grid/3_gray_scott/main.cpp

@@ -321,7 +321,7 @@ int main(int argc, char* argv[])
 		// visualization
 		if (i % 100 == 0)
 		{
-			Old.write("output",count);
+			Old.write_frame("output",count);
 			count++;
 		}
 	}

example/Grid/3_gray_scott_3d/Makefile

@@ -7,7 +7,7 @@ LDIR =
 OBJ = main.o
 
 %.o: %.cpp
-	$(CC) -O3 -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+	$(CC) -O3 -g -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
 
 gray_scott: $(OBJ)
 	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)

example/Grid/3_gray_scott_3d/main.cpp

@@ -6,6 +6,8 @@
  *
  * \page Grid_3_gs_3D Gray Scott in 3D
  *
+ * [TOC]
+ *
  * # Solving a gray scott-system in 3D # {#e3_gs_gray_scott}
  *
  * This example is just an extension of the 2D Gray scott example.
@@ -17,9 +19,25 @@
  * <img src="http://ppmcore.mpi-cbg.de/web/images/examples/gray_scott_3d/gs_alpha.png"/>
  * \endhtmlonly
  *
+ * More or less this example is the adaptation of the previous example to 3D
+ * with the improvement of using stencil iterator.
+ *
+ * ## Stencil iterator {#e3_gs_grat_scott_si}
+ *
+ * Stencil iterator require that you define a stencil,
+ *
+ * \snippet Grid/3_gray_scott_3d/main.cpp stencil def
+ *
+ * once is defined it is
+ * possible get and use a stencil iterator
+ *
+ * \snippet Grid/3_gray_scott_3d/main.cpp stencil get and use
+ *
+ * The rest of the example remain the same with the exception
+ * that the code has been extended in 3D.
+ *
  * \see \ref Grid_2_solve_eq
  *
- * \snippet Grid/3_gray_scott/main.cpp constants
  * 
  */
 
@@ -107,29 +125,13 @@ int main(int argc, char* argv[])
         size_t timeSteps = 5000;
 
 	// K and F (Physical constant in the equation)
-        double K = 0.014;
-        double F = 0.053;
-
-	//! \cond [init lib] \endcond
-
-	/*!
-	 * \page Grid_3_gs_3D Gray Scott in 3D
-	 *
-	 * Here we create 2 distributed grid in 2D Old and New. In particular because we want that
-	 * the second grid is distributed across processors in the same way we pass the decomposition
-	 * of the Old grid to the New one in the constructor with **Old.getDecomposition()**. Doing this,
-	 * we force the two grid to have the same decomposition.
-	 *
-	 * \snippet Grid/3_gray_scott/main.cpp init grid
-	 *
-	 */
-
-	//! \cond [init grid] \endcond
+        double K = 0.053;
+        double F = 0.014;
 
 	grid_dist_id<3, double, aggregate<double,double>> Old(sz,domain,g,bc);
 
 	// New grid with the decomposition of the old grid
-        grid_dist_id<3, double, aggregate<double,double>> New(Old.getDecomposition(),sz,g);
+    grid_dist_id<3, double, aggregate<double,double>> New(Old.getDecomposition(),sz,g);
 
 	
 	// spacing of the grid on x and y
@@ -149,46 +151,71 @@ int main(int argc, char* argv[])
 	timer tot_sim;
 	tot_sim.start();
 
+	//! \cond [stencil def] \endcond
+
+	static grid_key_dx<3> star_stencil_3D[7] = {{0,0,0},
+                                         	    {0,0,-1},
+												{0,0,1},
+												{0,-1,0},
+												{0,1,0},
+												{-1,0,0},
+												{1,0,0}};
+
+	//! \cond [stencil def] \endcond
+
 	for (size_t i = 0; i < timeSteps; ++i)
 	{
 		if (i % 300 == 0)
-			std::cout << "STEP: " << i << std::endl;
+		{std::cout << "STEP: " << i << std::endl;}
+
+		//! \cond [stencil get and use] \endcond
 
-		auto it = Old.getDomainIterator();
+		auto it = Old.getDomainIteratorStencil(star_stencil_3D);
 
 		while (it.isNext())
 		{
-			auto key = it.get();
+			// center point
+			auto Cp = it.getStencil<0>();
+
+			// plus,minus X,Y,Z
+			auto mx = it.getStencil<1>();
+			auto px = it.getStencil<2>();
+			auto my = it.getStencil<3>();
+			auto py = it.getStencil<4>();
+			auto mz = it.getStencil<5>();
+			auto pz = it.getStencil<6>();
 
 			// update based on Eq 2
-			New.get<U>(key) = Old.get<U>(key) + uFactor * (
-										Old.get<U>(key.move(x,1)) +
-										Old.get<U>(key.move(x,-1)) +
-										Old.get<U>(key.move(y,1)) +
-										Old.get<U>(key.move(y,-1)) +
-										Old.get<U>(key.move(z,1)) +
-										Old.get<U>(key.move(z,-1)) -
-										6.0*Old.get<U>(key)) +
-										- deltaT * Old.get<U>(key) * Old.get<V>(key) * Old.get<V>(key) +
-										- deltaT * F * (Old.get<U>(key) - 1.0);
+			New.get<U>(Cp) = Old.get<U>(Cp) + uFactor * (
+										Old.get<U>(mz) +
+										Old.get<U>(pz) +
+										Old.get<U>(my) +
+										Old.get<U>(py) +
+										Old.get<U>(mx) +
+										Old.get<U>(px) -
+										6.0*Old.get<U>(Cp)) +
+										- deltaT * Old.get<U>(Cp) * Old.get<V>(Cp) * Old.get<V>(Cp) +
+										- deltaT * F * (Old.get<U>(Cp) - 1.0);
 
 
 			// update based on Eq 2
-			New.get<V>(key) = Old.get<V>(key) + vFactor * (
-										Old.get<V>(key.move(x,1)) +
-										Old.get<V>(key.move(x,-1)) +
-										Old.get<V>(key.move(y,1)) +
-										Old.get<V>(key.move(y,-1)) +
-										Old.get<V>(key.move(z,1)) +
-                                                                                Old.get<V>(key.move(z,-1)) -
-										6*Old.get<V>(key)) +
-										deltaT * Old.get<U>(key) * Old.get<V>(key) * Old.get<V>(key) +
-										- deltaT * (F+K) * Old.get<V>(key);
+			New.get<V>(Cp) = Old.get<V>(Cp) + vFactor * (
+										Old.get<V>(mz) +
+										Old.get<V>(pz) +
+										Old.get<V>(my) +
+										Old.get<V>(py) +
+										Old.get<V>(mx) +
+                                        Old.get<V>(px) -
+										6*Old.get<V>(Cp)) +
+										deltaT * Old.get<U>(Cp) * Old.get<V>(Cp) * Old.get<V>(Cp) +
+										- deltaT * (F+K) * Old.get<V>(Cp);
 
 			// Next point in the grid
 			++it;
 		}
 
+		//! \cond [stencil get and use] \endcond
+
 		// Here we copy New into the old grid in preparation of the new step
 		// It would be better to alternate, but using this we can show the usage
 		// of the function copy. To note that copy work only on two grid of the same
@@ -199,11 +226,11 @@ int main(int argc, char* argv[])
 		// After copy we synchronize again the ghost part U and V
 		Old.ghost_get<U,V>();
 
-		// Every 30 time step we output the configuration for
+		// Every 500 time step we output the configuration for
 		// visualization
-		if (i % 60 == 0)
+		if (i % 500 == 0)
 		{
-			Old.write_frame("output",count,VTK_WRITER | FORMAT_BINARY);
+			Old.save("output_" + std::to_string(count));
 			count++;
 		}
 	}
@@ -229,4 +256,13 @@ int main(int argc, char* argv[])
 	openfpm_finalize();
 
 	//! \cond [finalize] \endcond
+
+	/*!
+	 * \page Grid_3_gs_3D Gray Scott in 3D
+	 *
+	 * # Full code # {#code}
+	 *
+	 * \include Grid/3_gray_scott_3d/main.cpp
+	 *
+	 */
 }

example/Grid/3_gray_scott_3d_vectorization/Makefile

+include ../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o update_new.o
+
+%.o: %.f90
+	mpif90 -ffree-line-length-none -fno-range-check -fno-second-underscore  -fimplicit-none  -mavx -O3 -c -g -o $@ $<
+
+%.o: %.cpp
+	$(CC) -O3 -mavx  -g -c --std=c++11 -Wno-ignored-attributes  -o  $@ $< $(INCLUDE_PATH) -I/where/is/vc/installation/include
+
+gray_scott: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS) -L/where/is/vc/installation/lib  # -lVc (Add -lVc if you use VCDevel)
+
+all: gray_scott
+
+run: all
+	mpirun -np 4 ./gray_scott
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core gray_scott
+

example/Grid/3_gray_scott_3d_vectorization/config.cfg

+[pack]
+files = main.cpp update_new.f90 Makefile

example/Grid/3_gray_scott_3d_vectorization/update_new.f90

+subroutine update_new ( &
+                              lo, hi, &
+                               u, ulo, uhi, &
+                               v, vlo, vhi, &
+                               u_new, nulo, nuhi, &
+                               v_new, nvlo, nvhi, &
+                               dt, uFactor, vFactor,F,Kf) bind(C, name="update_new")
+
+  implicit none
+
+  integer, intent(in) :: lo(3), hi(3)
+  integer, intent(in) :: ulo(3), uhi(3)
+  integer, intent(in) :: vlo(3), vhi(3)
+  integer, intent(in) :: nulo(3), nuhi(3), nvlo(3), nvhi(3)
+  real*8, intent(in)    :: u  (ulo(1):uhi(1),ulo(2):uhi(2),ulo(3):uhi(3))
+  real*8, intent(in)    :: v  (vlo(1):vhi(1),vlo(2):vhi(2),vlo(3):vhi(3))
+  real*8, intent(inout) :: u_new( nulo(1): nuhi(1), nulo(2): nuhi(2), nulo(3): nuhi(3))
+  real*8, intent(inout) :: v_new( nvlo(1): nvhi(1), nvlo(2): nvhi(2), nvlo(3): nvhi(3))
+  real*8, intent(in) :: dt, F, Kf, uFactor, vFactor
+
+  ! local variables
+  integer i,j,k
+
+  ! x-fluxes
+  do       k = lo(3), hi(3)
+     do    j = lo(2), hi(2)
+        do i = lo(1), hi(1)
+           u_new(i,j,k) = u(i,j,k) + uFactor * ( u(i+1,j,k) + u(i-1,j,k) + &
+                                    u(i,j+1,k) + u(i,j-1,k) + &
+                                    u(i,j,k-1) + u(i,j,k+1) - &
+                                    6.0*u(i,j,k) ) - &
+                                    dt * u(i,j,k)*v(i,j,k)*v(i,j,k) - &
+                                    dt * F * (u(i,j,k) - 1.0)
+
+
+           v_new(i,j,k) = v(i,j,k) + vFactor * ( v(i+1,j,k) + v(i-1,j,k) + &
+                                    v(i,j+1,k) + v(i,j-1,k) + &
+                                    v(i,j,k-1) + v(i,j,k+1) - &
+                                    6.0*v(i,j,k) ) + &
+                                    dt * u(i,j,k)*v(i,j,k)*v(i,j,k) - &
+                                    dt * (F+Kf) * v(i,j,k)
+        end do
+     end do
+  end do
+
+
+end subroutine update_new
+

example/Numerics/PS-CMA-ES/Makefile

+include ../../example.mk
+
+CC=mpic++
+
+LDIR =
+
+OBJ = main.o
+
+%.o: %.cpp
+	$(CC) -O3 -g -c --std=c++11 -o $@ $< $(INCLUDE_PATH)
+
+ps_cma_es: $(OBJ)
+	$(CC) -o $@ $^ $(CFLAGS) $(LIBS_PATH) $(LIBS)
+
+all: ps_cma_es
+
+run: all
+	mpirun -np 2 ./ps_cma_es
+
+.PHONY: clean all run
+
+clean:
+	rm -f *.o *~ core ps_cma_es
+

example/Numerics/PS-CMA-ES/config.cfg

+[pack]
+files = main.cpp f15_cec_const.hpp  f15_cec_fun.hpp Makefile

example/Numerics/PS-CMA-ES/f15_cec_fun.hpp

+/*
+ * f15_cec_fun.hpp
+ *
+ *  Created on: Jan 14, 2018
+ *      Author: i-bird
+ */
+
+#ifndef EXAMPLE_NUMERICS_PS_CMA_ES_F15_CEC_FUN_HPP_
+#define EXAMPLE_NUMERICS_PS_CMA_ES_F15_CEC_FUN_HPP_
+
+#include "f15_cec_const.hpp"
+#include <limits>
+#include <math.h>
+
+template<unsigned int dim>
+void Job15(int funcnr,Eigen::VectorXd & vars,double & res)
+{
+      // local used vars
+      double sum,sum1,sum2,prod,e1,e2;
+      int i,j,k;
+      // weierstrass vars
+      int Kmax = 20;
+      const double a_c = 0.5;
+      const double b_c = 3.0;
+
+
+      if (funcnr < 2)
+      {
+        // rastrigin
+         sum = 10.0 * dim;
+         for (size_t i = 0 ; i < dim ; i++)
+         {
+            sum += vars(i)*vars(i);
+            sum -= 10.0*cos(2*M_PI*vars[i]);
+         }
+
+         res = sum;
+      }
+      else if (funcnr < 4)
+      {
+        // weierstrass
+         sum1 = 0.0;
+         sum2 = 0.0;
+         double a_k = 1.0;
+         double b_k = 1.0;
+         for (size_t i = 0 ; i < dim ; i++)
+         {
+            a_k = 1.0;
+            b_k = 1.0;
+            for (size_t j = 0 ; j <= Kmax ; j++, a_k *= a_c,b_k *= b_c)
+            {
+                sum1 = sum1 + a_k * cos((M_PI)*2.0 * b_k * (vars(i)+0.5));
+            }
+         }
+         a_k = 1.0;
+         b_k = 1.0;
+         for (size_t j = 0 ; j <= Kmax ; j++, a_k *= a_c, b_k *= b_c)
+         {
+            sum2 = sum2 + a_k * cos((M_PI)*2.0 * b_k * (0.5));
+         }
+         res = sum1 - sum2*dim;
+      }
+      else if (funcnr < 6)
+      {
+         // griewank
+         prod = 1;
+         sum = 0.0;
+         for (size_t i = 1 ; i <= dim ; i++)
+         {
+        	 sum= sum + (vars(i-1)*vars(i-1))/4000.0;
+        	 prod=prod * cos(vars(i-1)/(sqrt(double(i))));
+         }
+         res = sum-prod+1;
+      }
+      else if (funcnr < 8)
+      {
+        // ackley
+        e1 = 0.0;
+        e2 = 0.0;
+        for (size_t i = 0 ; i < dim ; i++)
+        {
+            e1 = e1 + vars(i)*vars(i);
+            e2 = e2 + cos(2.0*M_PI*vars(i));
+        }
+        res = exp(1.0) + 20.0 - 20*exp(-0.2*sqrt(e1/dim));
+        res = res - exp(e2/dim);
+      }
+      else if (funcnr <= 10)
+      {
+        // sphere
+        sum = vars.transpose() * vars;
+        res = sum;
+      }