diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index 4309f6f00ef93503cbc92786b1c59b4d44db95aa..cf32a3227a98c11a005b97d2dc13f5a29c623727 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -24,10 +24,10 @@ centos_run:
     - ./run.sh $CI_PROJECT_DIR unused 4 pdata 0 $CI_COMMIT_REF_NAME
     - ./run.sh $CI_PROJECT_DIR unused 5 pdata 0 $CI_COMMIT_REF_NAME
     - cd openfpm_numerics
-    - ./run.sh $CI_PROJECT_DIR unused 1 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 3 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 4 0 0 numerics
+    - ./run.sh $CI_PROJECT_DIR unused 1 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 3 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 4 0 0 numerics $CI_COMMIT_REF_NAME
 
 mac_build:
   variables:
@@ -53,7 +53,7 @@ mac_run:
     - ./run.sh $CI_PROJECT_DIR unused 1 0 0 $CI_COMMIT_REF_NAME
     - ./run.sh $CI_PROJECT_DIR unused 2 0 0 $CI_COMMIT_REF_NAME
     - cd openfpm_numerics
-    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics
+    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics $CI_COMMIT_REF_NAME
 
 ubuntu_build:
   stage: build
@@ -78,8 +78,8 @@ ubuntu_run:
     - ./run.sh $CI_PROJECT_DIR unused 2 pdata 0 $CI_COMMIT_REF_NAME
     - ./run.sh $CI_PROJECT_DIR unused 3 pdata 0 $CI_COMMIT_REF_NAME
     - cd openfpm_numerics
-    - ./run.sh $CI_PROJECT_DIR unused 1 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 3 0 0 numerics
-    - ./run.sh $CI_PROJECT_DIR unused 4 0 0 numerics
+    - ./run.sh $CI_PROJECT_DIR unused 1 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 2 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 3 0 0 numerics $CI_COMMIT_REF_NAME
+    - ./run.sh $CI_PROJECT_DIR unused 4 0 0 numerics $CI_COMMIT_REF_NAME
 
diff --git a/CHANGELOG.md b/CHANGELOG.md
index d5ed2e9d81077b8c78277cc4964c9ba659354d93..1c17bcb16955bef021402290af06bd44677e9799 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -1,6 +1,18 @@
 # Change Log
 All notable changes to this project will be documented in this file.
 
+## [4.0.0] September 2021 (Codename Heisenberg)
+
+- Adding DCPSE, Level-set based numerics (Closest-point)
+
+### Changes
+
+- Particles writer use now the new XML pvtp Paraview format
+- Particles constructor does not accept discordance precision on space example:
+            If I create a Particle set with vector_dist<3,double,.......>
+            I MUST give in the constructor a domain Box<3,double> and a Ghost<3,double>, giving a Ghost<3,float> will produce
+            an error
+
 ## [3.3.0] April 2021 (Codename Vega)
 
 - Adding support for HIP and AMD GPU. (Only particles) 1_gpu_first/7_sph_dlb_gpu/7_sph_dlb_gpu_opt are compatible with HIP
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 8937baf0560bff3f27bbe5343d665bbc50d23ca0..e24e7b719f27c4421ab4ab179aa074cca6227ea2 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -11,7 +11,7 @@ if (POLICY CMP0074)
   cmake_policy(SET CMP0074 NEW)
 endif()
 
-set(openfpm_VERSION 3.3.0)
+set(openfpm_VERSION 4.0.0)
 
 list(APPEND CMAKE_MODULE_PATH ${CMAKE_CURRENT_LIST_DIR}/cmake_modules/)
 
@@ -294,6 +294,12 @@ list(GET VERSION_LIST 2 OPENFPM_VERSION_PATCH)
 
 if (CPACK_RUN_INSTALL_DEPENDENCIES)
 
+        ###### Fix post inst script ######
+
+	
+
+        ######                      ######
+
 	set(CPACK_PACKAGE_DESCRIPTION_SUMMARY "OpenFPM distributed data-structures")
 	set(CPACK_PACKAGE_VENDOR "IBirdSoft")
 	set(CPACK_PACKAGE_DESCRIPTION_FILE "${CMAKE_CURRENT_SOURCE_DIR}/README.txt")
@@ -310,7 +316,7 @@ if (CPACK_RUN_INSTALL_DEPENDENCIES)
 	set(CPACK_DEBIAN_PACKAGE_MAINTAINER Pietro Incardona)
 	set(CPACK_RPM_POST_INSTALL_SCRIPT_FILE ${CMAKE_CURRENT_SOURCE_DIR}/src/scripts/postinst)
 	set(CPACK_POSTFLIGHT_OPENFPM_SCRIPT ${CMAKE_CURRENT_SOURCE_DIR}/src/scripts/postflight)
-	set(CPACK_DEBIAN_PACKAGE_CONTROL_EXTRA "${CMAKE_CURRENT_SOURCE_DIR}/src/script/postinst")
+	set(CPACK_DEBIAN_PACKAGE_CONTROL_EXTRA ${CMAKE_CURRENT_SOURCE_DIR}/src/scripts/postinst)
 	set(CPACK_RESOURCE_FILE_README "${CMAKE_CURRENT_SOURCE_DIR}/README.txt")
 
 	install(FILES $ENV{DEP_PACKING}/openfpm_vars
@@ -361,6 +367,14 @@ if (CPACK_RUN_INSTALL_DEPENDENCIES)
                 DESTINATION dependencies/
                 COMPONENT OpenFPM)
 
+	install(DIRECTORY $ENV{DEP_PACKING}/BLITZ
+                DESTINATION dependencies/
+                COMPONENT OpenFPM)
+
+	install(DIRECTORY $ENV{DEP_PACKING}/ALGOIM
+                DESTINATION dependencies/
+                COMPONENT OpenFPM)
+
 endif()
 
 include(CPack)
diff --git a/example/Numerics/OdeInt/Makefile b/example/Numerics/OdeInt/Makefile
index 7bb01f1ab84fb758dcac072f7a91fa6ce79d3a2d..49784aa93f9eb62a4f0fa3fb77ccffb40999513e 100644
--- a/example/Numerics/OdeInt/Makefile
+++ b/example/Numerics/OdeInt/Makefile
@@ -19,7 +19,7 @@ example_odeint2: $(OBJ2)
 all: example_odeint example_odeint2
 
 run: all
-	mpirun -np 2 ./example_odeint ./example_odeint
+	mpirun  --oversubscribe -np 2 ./example_odeint && mpirun --oversubscribe -np 2 ./example_odeint2
 
 .PHONY: clean all run
 
diff --git a/example/Vector/2_expressions/main.cpp b/example/Vector/2_expressions/main.cpp
index f20b0a2828605a82c2299c3630b689ad5043b6a8..50f04ef53e14f56dd8d47c92bdc6230afd007749 100644
--- a/example/Vector/2_expressions/main.cpp
+++ b/example/Vector/2_expressions/main.cpp
@@ -115,13 +115,13 @@ int main(int argc, char* argv[])
 	openfpm_init(&argc,&argv);
 
 	// Here we define our domain a 2D box with intervals from 0 to 1.0
-	Box<3,float> domain({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> domain({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Here we define the boundary conditions of our problem
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// extended boundary around the domain, and the processor domain
-	Ghost<3,float> g(0.01);
+	Ghost<3,double> g(0.01);
 
 	// Delta t
 	double dt = 0.01;
diff --git a/example/Vector/3_molecular_dynamic/main.cpp b/example/Vector/3_molecular_dynamic/main.cpp
index 8df3ee93947c62e02728e304f82c999e0396347e..749047bb3593ccc28faea6b5669e47ce7311bad6 100644
--- a/example/Vector/3_molecular_dynamic/main.cpp
+++ b/example/Vector/3_molecular_dynamic/main.cpp
@@ -308,13 +308,13 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {10,10,10};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	double dt = 0.0005;
 	double sigma12 = pow(sigma,12);
diff --git a/example/Vector/3_molecular_dynamic/main_expr.cpp b/example/Vector/3_molecular_dynamic/main_expr.cpp
index 090631764a5da0af1e4ab4ae1d95532a7c11badc..352391e519e9ddd4a9bc7fa6b3e2ec478fcadfa6 100644
--- a/example/Vector/3_molecular_dynamic/main_expr.cpp
+++ b/example/Vector/3_molecular_dynamic/main_expr.cpp
@@ -56,9 +56,9 @@ int main(int argc, char* argv[])
 	Vcluster<> & vcl = create_vcluster();
 
 	size_t sz[3] = {10,10,10};
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	ln_force lf(sigma12,sigma6,r_cut*r_cut);
 	ln_potential lp(sigma12,sigma6,r_cut*r_cut,shift);
diff --git a/example/Vector/3_molecular_dynamic/main_expr_paper.cpp b/example/Vector/3_molecular_dynamic/main_expr_paper.cpp
index cfe7b19efa60b40997913e2e0dde9f359d288dac..e166299ae7e268f776e6221c8178160e077e66df 100644
--- a/example/Vector/3_molecular_dynamic/main_expr_paper.cpp
+++ b/example/Vector/3_molecular_dynamic/main_expr_paper.cpp
@@ -53,9 +53,9 @@ int main(int argc, char* argv[])
 	///// Define initialization grid, simulation box, periodicity
 	///// and ghost
 	size_t sz[3] = {10,10,10};
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	///// Lennard_jones potential object used in applyKernel_in
 	ln_force lennard_jones;
diff --git a/example/Vector/3_molecular_dynamic/main_vl.cpp b/example/Vector/3_molecular_dynamic/main_vl.cpp
index 6063ea6cfd5a4603728d66972a7f2b0c94f9267a..1538033fe349173755440b01e9bc090e57fb29c8 100644
--- a/example/Vector/3_molecular_dynamic/main_vl.cpp
+++ b/example/Vector/3_molecular_dynamic/main_vl.cpp
@@ -278,13 +278,13 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {10,10,10};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_gskin);
+	Ghost<3,double> ghost(r_gskin);
 
 	vector_dist<3,double, aggregate<double[3],double[3]> > vd(0,box,bc,ghost);
 
diff --git a/example/Vector/4_reorder/main_comp_ord.cpp b/example/Vector/4_reorder/main_comp_ord.cpp
index e8d4c9f1d66426f358a564a24c4cf794d374d270..81f95d7e9435eb0d243c9f6c75bc3da27dc73ca7 100644
--- a/example/Vector/4_reorder/main_comp_ord.cpp
+++ b/example/Vector/4_reorder/main_comp_ord.cpp
@@ -226,7 +226,7 @@ int main(int argc, char* argv[])
 	//! \cond [vect create] \endcond
 
 	double dt = 0.0001;
-	float r_cut = 0.03;
+	double r_cut = 0.03;
 	double sigma = r_cut/3.0;
 	double sigma12 = pow(sigma,12);
 	double sigma6 = pow(sigma,6);
@@ -241,13 +241,13 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {40,40,40};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	vector_dist<3,double, aggregate<double[3],double[3]> > vd(0,box,bc,ghost);
 
diff --git a/example/Vector/4_reorder/main_data_ord.cpp b/example/Vector/4_reorder/main_data_ord.cpp
index d9b97b4797b0bc8abd0d2298f4ede8d7016d668f..4e5022ed6376761084c2f8d1595813a61c9a28b4 100644
--- a/example/Vector/4_reorder/main_data_ord.cpp
+++ b/example/Vector/4_reorder/main_data_ord.cpp
@@ -41,7 +41,7 @@ int main(int argc, char* argv[])
 	//! \cond [vect create] \endcond
 
 	double dt = 0.0001;
-	float r_cut = 0.03;
+	double r_cut = 0.03;
 	double sigma = r_cut/3.0;
 	double sigma12 = pow(sigma,12);
 	double sigma6 = pow(sigma,6);
@@ -56,13 +56,13 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {40,40,40};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	// Create vector
 	vector_dist<3,double, aggregate<double[3],double[3]> > vd(0,box,bc,ghost);
diff --git a/example/Vector/4_reorder/main_expr.cpp b/example/Vector/4_reorder/main_expr.cpp
index 86c7e31bf41bbd94ef81c469fa90617ed261b0b1..82c0ca9508116928341acd754497cd46d17259c7 100644
--- a/example/Vector/4_reorder/main_expr.cpp
+++ b/example/Vector/4_reorder/main_expr.cpp
@@ -51,9 +51,9 @@ int main(int argc, char* argv[])
 	Vcluster & vcl = create_vcluster();
 
 	size_t sz[3] = {10,10,10};
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
-	Ghost<3,float> ghost(r_cut);
+	Ghost<3,double> ghost(r_cut);
 
 	ln_force lf(sigma12,sigma6);
 	ln_potential lp(sigma12,sigma6);
diff --git a/example/Vector/5_molecular_dynamic_sym/main.cpp b/example/Vector/5_molecular_dynamic_sym/main.cpp
index 1e1cfe2f844450568901f871970e07e58bdc4422..d44029d36c6776cf1a5c489d3ca22082e78a6fd8 100644
--- a/example/Vector/5_molecular_dynamic_sym/main.cpp
+++ b/example/Vector/5_molecular_dynamic_sym/main.cpp
@@ -290,13 +290,13 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {10,10,10};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_gskin);
+	Ghost<3,double> ghost(r_gskin);
 
 	vector_dist<3,double, aggregate<double[3],double[3]> > vd(0,box,bc,ghost);
 
diff --git a/example/Vector/5_molecular_dynamic_sym_crs/main.cpp b/example/Vector/5_molecular_dynamic_sym_crs/main.cpp
index 2985148412526fd4fcb845bc43e053a9b3f8c58a..5b45b06729b69f5d79f9460daa1ba0702c375f5c 100644
--- a/example/Vector/5_molecular_dynamic_sym_crs/main.cpp
+++ b/example/Vector/5_molecular_dynamic_sym_crs/main.cpp
@@ -363,7 +363,7 @@ int main(int argc, char* argv[])
 	size_t sz[3] = {10,10,10};
 
 	// domain
-	Box<3,float> box({0.0,0.0,0.0},{1.0,1.0,1.0});
+	Box<3,double> box({0.0,0.0,0.0},{1.0,1.0,1.0});
 
 	// Boundary conditions
 	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
@@ -371,7 +371,7 @@ int main(int argc, char* argv[])
 	//! \cond [remove_lower_part] \endcond
 
 	// ghost, big enough to contain the interaction radius
-	Ghost<3,float> ghost(r_gskin);
+	Ghost<3,double> ghost(r_gskin);
 	ghost.setLow(0,0.0);
 	ghost.setLow(1,0.0);
 	ghost.setLow(2,0.0);
diff --git a/example/Vector/7_SPH_dlb_gpu/Makefile b/example/Vector/7_SPH_dlb_gpu/Makefile
index 8813b2db0ecd3b2b8764579fea94295107af65ea..b4ebc34312ce5d2227f06fe1a196193e2be27469 100644
--- a/example/Vector/7_SPH_dlb_gpu/Makefile
+++ b/example/Vector/7_SPH_dlb_gpu/Makefile
@@ -8,7 +8,7 @@ CUDA_CC_LINK=
 CC=mpic++
 ifdef HIP
         CUDA_CC=hipcc
-        CUDA_OPTIONS=-D__NVCC__ -D__HIP__ -DCUDART_VERSION=11000 -D__CUDACC__ -D__CUDACC_VER_MAJOR__=11 -D__CUDACC_VER_MINOR__=0 -D__CUDACC_VER_BUILD__=0
+        CUDA_OPTIONS= -D__NVCC__ -D__HIP__ -DCUDART_VERSION=11000 -D__CUDACC__ -D__CUDACC_VER_MAJOR__=11 -D__CUDACC_VER_MINOR__=0 -D__CUDACC_VER_BUILD__=0
         LIBS_SELECT=$(LIBS)
         CC=hipcc
 	CUDA_CC_LINK=hipcc
@@ -44,7 +44,7 @@ sph_dlb_test: OPT += -DTEST_RUN
 sph_dlb_test: sph_dlb
 
 %.o: %.cu
-	$(CUDA_CC) $(CUDA_OPTIONS) -O3 $(OPT) -g -c --std=c++14 -o $@ $< $(INCLUDE_PATH_NVCC)
+	$(CUDA_CC) $(CUDA_OPTIONS) $(OPT) -g -c --std=c++14 -o $@ $< $(INCLUDE_PATH_NVCC)
 
 %.o: %.cpp
 	$(CC) -O3 $(OPT) -g -c --std=c++14 -o $@ $< $(INCLUDE_PATH)
diff --git a/install b/install
index b86dd8d36f4a5b5918f38a86cfb780e41216923f..3d6b4ef6af9e7900eab5a6b74af2780d4c4233c7 100755
--- a/install
+++ b/install
@@ -219,9 +219,9 @@ echo -e "Installing requirements into: $i_dir "
 
 ## call the configure script
 
-### Installing PETSC
+### Installing PETSC and Extended numerics
 
-if [ ! -d "$i_dir/PETSC" -o ! -f "$i_dir/PETSC/include/petsc.h" -o ! -d "$i_dir/EIGEN" ]; then
+if [ ! -d "$i_dir/PETSC" -o ! -f "$i_dir/PETSC/include/petsc.h" -o ! -d "$i_dir/EIGEN" -o ! -d "$i_dir/BLITZ" -o ! -d "$i_dir/ALGOIM" ]; then
     echo -e "\033[1;34m Optional packages  \033[0m"
     echo -e "\033[1mDo you want to install linear algebra packages and extended Level-set functionalities ?(y/n)\033[0m"
     if [ $sq -eq 0 ]; then
diff --git a/script/solve_gdbserver b/script/solve_gdbserver
index a92893f127b4d5dac842896619333271b05fca70..cd7404c8ab5727c3694992efe4d95df81bd72bbb 100755
--- a/script/solve_gdbserver
+++ b/script/solve_gdbserver
@@ -19,7 +19,7 @@ elif [ x"$1" = x"linux"  ]; then
 	elif [ x"$pcman" = x"pacman" ]; then
 		package_name=gdb
 	elif [ x"$pcman" = x"apt-get" ]; then
-		package_name=gdbsever
+		package_name=gdbserver
 	else
 		package_name=gdb-gdbserver
 	fi
diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt
index 3148669a74322677e37a8342ffde69d6b35562d7..eee5aed20acdc3480b6526182051075b789c9e3e 100644
--- a/src/CMakeLists.txt
+++ b/src/CMakeLists.txt
@@ -13,7 +13,8 @@ if(CUDA_FOUND OR CUDA_ON_CPU OR HIP_FOUND)
 	    Decomposition/cuda/decomposition_cuda_tests.cu
 	    Vector/cuda/vector_dist_gpu_unit_tests.cu
 	    Decomposition/cuda/Domain_icells_cart_unit_test.cu
-	    Amr/tests/amr_base_gpu_unit_tests.cu)
+		Amr/tests/amr_base_gpu_unit_tests.cu
+		Grid/tests/grid_dist_id_unit_test.cu)
 endif()
 
 if(CMAKE_CXX_COMPILER_ID STREQUAL "Intel")
@@ -61,8 +62,7 @@ if ( HIP_ENABLE AND HIP_FOUND )
 							  Decomposition/tests/CartDecomposition_unit_test.cpp
 							  Decomposition/tests/shift_vect_converter_tests.cpp
 							  Vector/performance/vector_dist_performance_util.cpp
-							  lib/pdata.cpp
-							  test_multiple_o.cpp
+							  lib/pdata.cpp 
 							  )
 
 
@@ -89,9 +89,10 @@ else()
 							  Decomposition/tests/CartDecomposition_unit_test.cpp
 							  Decomposition/tests/shift_vect_converter_tests.cpp
 							  Vector/performance/vector_dist_performance_util.cpp
-							  lib/pdata.cpp test_multiple_o.cpp)
+							  lib/pdata.cpp)
 
 	add_library(ofpm_pdata STATIC lib/pdata.cpp)
+	add_library(ofpm_pdata_python SHARED lib/pdata_python.cpp)
 
 endif()
 
@@ -111,12 +112,14 @@ if (CMAKE_CXX_COMPILER_ID MATCHES "Clang")
 	target_compile_options(pdata PRIVATE "-Wno-macro-redefined")
 endif()
 if ( CMAKE_COMPILER_IS_GNUCC )
-    target_compile_options(pdata PRIVATE "-Wno-deprecated-declarations")
-	target_compile_options(pdata PRIVATE "-Wno-undefined-var-template")
-	target_compile_options(pdata PRIVATE "-Wno-macro-redefined")
-    #target_compile_options(pdata PRIVATE "-fsanitize=address")
-    #    target_link_options(pdata PRIVATE "-fsanitize=address")
+#target_compile_options(pdata PRIVATE "-Wno-deprecated-declarations")
+#target_compile_options(pdata PRIVATE "-Wno-undefined-var-template")
+#target_compile_options(pdata PRIVATE "-Wno-macro-redefined")
+#target_compile_options(pdata PRIVATE "-fsanitize=address")
+#    target_link_options(pdata PRIVATE "-fsanitize=address")
     if (TEST_COVERAGE)
+######################################## Add mem tests ############################
+
         target_compile_options(pdata PRIVATE $<$<COMPILE_LANGUAGE:CXX>: -fprofile-arcs -ftest-coverage>)
     endif()
 endif()
@@ -212,6 +215,51 @@ add_definitions(-DSCAN_WITH_CUB)
 
 if (TEST_COVERAGE)
     target_link_libraries(pdata -lgcov --coverage)
+
+    get_target_property(PDATA_SOURCES pdata SOURCES)
+    get_target_property(MEM_SOURCES mem SOURCES)
+    get_target_property(DATA_SOURCES mem_map SOURCES)
+    get_target_property(VCLUSTER_SOURCES vcluster_test SOURCES)
+    get_target_property(IO_SOURCES io SOURCES)
+    get_target_property(NUMERIC_SOURCES numerics SOURCES)
+    get_target_property(PDATA_INCLUDES pdata INCLUDE_DIRECTORIES)
+    get_target_property(NUMERIC_INCLUDES numerics INCLUDE_DIRECTORIES)
+    get_target_property(PDATA_LIBS pdata LINK_LIBRARIES)
+    get_target_property(NUMERIC_LIBS numerics LINK_LIBRARIES)
+
+    list(TRANSFORM MEM_SOURCES PREPEND "../openfpm_devices/src/")
+    list(TRANSFORM DATA_SOURCES PREPEND "../openfpm_data/src/")
+    list(TRANSFORM VCLUSTER_SOURCES PREPEND "../openfpm_vcluster/src/")
+    list(TRANSFORM IO_SOURCES PREPEND "../openfpm_io/src/")
+    list(TRANSFORM NUMERIC_SOURCES PREPEND "../openfpm_numerics/src/")
+
+    set(PDATA_FULL_SOURCES ${PDATA_SOURCES} ${MEM_SOURCES} ${DATA_SOURCES} ${VCLUSTER_SOURCES} ${IO_SOURCES} ${NUMERIC_SOURCES})
+    set(PDATA_FULL_SOURCES_CU ${PDATA_SOURCES} ${MEM_SOURCES} ${DATA_SOURCES} ${VCLUSTER_SOURCES} ${IO_SOURCES} ${NUMERIC_SOURCES})
+    list(FILTER PDATA_FULL_SOURCES_CU INCLUDE REGEX ^.*\.cu)
+    set_source_files_properties(${PDATA_FULL_SOURCES_CU} PROPERTIES COMPILE_FLAGS "-D__NVCC__ -DCUDART_VERSION=11000")
+
+    if (CUDA_ON_CPU)
+        set_source_files_properties(${PDATA_FULL_SOURCES_CU} PROPERTIES LANGUAGE CXX)
+    endif()
+
+    add_executable(pdata_full ${PDATA_FULL_SOURCES})
+    
+    target_compile_definitions(pdata_full PUBLIC NO_INIT_AND_MAIN)
+    target_compile_options(pdata_full PRIVATE $<$<COMPILE_LANGUAGE:CXX>: -fprofile-arcs -ftest-coverage>)
+    target_link_options(pdata_full PUBLIC -fprofile-arcs -ftest-coverage)
+
+    target_link_options(pdata_full PUBLIC -fprofile-arcs -ftest-coverage)
+
+    if (ENABLE_ASAN)
+        target_compile_options(pdata_full PUBLIC $<$<COMPILE_LANGUAGE:CUDA>: -Xcompiler "-fsanitize=address -fno-optimize-sibling-calls -fsanitize-address-use-after-scope -fno-omit-frame-pointer -g" >)
+        target_compile_options(pdata_full PRIVATE $<$<COMPILE_LANGUAGE:CXX>: -fsanitize=address -fno-optimize-sibling-calls -fsanitize-address-use-after-scope -fno-omit-frame-pointer -g >)
+    endif()
+
+    target_link_libraries(pdata_full -lgcov --coverage)
+    target_link_libraries(pdata_full ${PDATA_LIBS})
+    target_link_libraries(pdata_full ${NUMERIC_LIBS})
+    target_include_directories(pdata_full PUBLIC ${PDATA_INCLUDES})
+    target_include_directories(pdata_full PUBLIC ${NUMERIC_INCLUDES})
 endif()
 
 
diff --git a/src/Decomposition/common.hpp b/src/Decomposition/common.hpp
index b0994f259a1427c3924842380fd72634572c52ab..4ed7d3cdc6f7ed7ab2fe79e918fc32066f19e6bf 100755
--- a/src/Decomposition/common.hpp
+++ b/src/Decomposition/common.hpp
@@ -126,7 +126,7 @@ struct Box_sub_k
 	//! Where this sub_domain live
 	comb<dim> cmb;
 
-	//! k \see getLocalGhostIBoxE
+	//! k \see getLocalIGhostE
 	long int k;
 
 	Box_sub_k()
@@ -138,21 +138,6 @@ struct Box_sub_k
 
 template<unsigned int dim,typename T> using Box_map = aggregate<Box<dim,T>,long int>;
 
-/*template<unsigned int dim,typename T>
-struct Box_map
-{
-	typedef boost::fusion::vector<Box<dim,T>,long int> type;
-
-	type data;
-
-	static bool noPointers()
-	{
-		return true;
-	}
-
-	static const unsigned int max_prop = 2;
-};*/
-
 //! Case for local ghost box
 template<unsigned int dim, typename T>
 struct lBox_dom
diff --git a/src/Decomposition/ie_loc_ghost.hpp b/src/Decomposition/ie_loc_ghost.hpp
index 231828e1e71503d525a371be2bb14760799f134e..207072bbd165aa88905a3b7ec52a94eeefa329ea 100755
--- a/src/Decomposition/ie_loc_ghost.hpp
+++ b/src/Decomposition/ie_loc_ghost.hpp
@@ -350,9 +350,9 @@ public:
 		return loc_ghost_box;
 	}
 
-	/*! \brief Get the number of local sub-domains
+	/*! \brief Get the number of local ghost boxes
 	 *
-	 * \return the number of local sub-domains
+	 * \return the number of local ghost boxes
 	 *
 	 *
 	 */
@@ -385,7 +385,7 @@ public:
 		return loc_ghost_box.get(id).ibx.size();
 	}
 
-	/*! \brief For the sub-domain i intersected with a surrounding sub-domain enlarged. Produce a internal ghost box from
+	/*! \brief For the sub-domain i intersected with a surrounding sub-domain enlarged j. Produce a internal ghost box from
 	 *        the prospecive of i and an associated external ghost box from the prospective of j.
 	 *        In order to retrieve the information about the external ghost box we have to use getLocalEGhostBox(x,k).
 	 *        where k is the value returned by getLocalIGhostE(i,j) and x is the value returned by
diff --git a/src/Grid/Iterators/grid_dist_id_iterator.hpp b/src/Grid/Iterators/grid_dist_id_iterator.hpp
index 308567b31fff370d39ac7627b29548b65b03c940..fd8dbfa2fa0346b1305cd20d1ebeba5f0b0008b3 100644
--- a/src/Grid/Iterators/grid_dist_id_iterator.hpp
+++ b/src/Grid/Iterators/grid_dist_id_iterator.hpp
@@ -20,6 +20,8 @@
 #include "SparseGridGpu/encap_num.hpp"
 #endif
 
+#include "Grid/cuda/grid_dist_id_kernels.cuh"
+
 template<unsigned int dim>
 struct launch_insert_sparse_lambda_call
 {
@@ -183,6 +185,56 @@ struct launch_insert_sparse
 	}
 };
 
+template<unsigned int dim>
+struct launch_set_dense
+{
+	template<typename grid_type, typename ite_type, typename lambda_f2>
+	__device__ void operator()(grid_type & grid, ite_type itg, lambda_f2 f2)
+	{
+#ifdef __NVCC__
+
+		printf("grid on GPU Dimension %d not implemented, yet\n",(int)dim);
+
+#endif
+	}
+};
+
+template<>
+struct launch_set_dense<2>
+{
+	template<typename grid_type, typename ite_type, typename lambda_f2>
+	__device__ void operator()(grid_type & grid, ite_type itg, lambda_f2 f2)
+	{
+#ifdef __NVCC__
+
+		GRID_ID_2_GLOBAL(itg);
+
+		auto obj = grid.get_o(key);
+
+		f2(obj,keyg.get(0),keyg.get(1));
+
+#endif
+	}
+};
+
+template<>
+struct launch_set_dense<3>
+{
+	template<typename grid_type, typename ite_type, typename lambda_f2>
+	__device__ void operator()(grid_type & grid, ite_type itg, lambda_f2 f2)
+	{
+#ifdef __NVCC__
+
+		GRID_ID_3_GLOBAL(itg);
+
+		auto obj = grid.get_o(key);
+
+		f2(obj,keyg.get(0),keyg.get(1),keyg.get(2));
+
+#endif
+	}
+};
+
 template<bool is_free>
 struct selvg
 {
diff --git a/src/Grid/cuda/grid_dist_id_iterator_gpu.cuh b/src/Grid/cuda/grid_dist_id_iterator_gpu.cuh
index 66d1cbbdad2e316c22a311009a110cd122907de5..b3a322a393767148197f6def300c7254dcaa417d 100644
--- a/src/Grid/cuda/grid_dist_id_iterator_gpu.cuh
+++ b/src/Grid/cuda/grid_dist_id_iterator_gpu.cuh
@@ -77,28 +77,28 @@ class grid_dist_id_iterator_gpu
 	* \return itself
 	*
 	*/
-	grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & operator=(const grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & tmp)
-	{
-		g_c = tmp.g_c;
-		gdb_ext = tmp.gdb_ext;
+//	grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & operator=(const grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & tmp)
+//	{
+//		g_c = tmp.g_c;
+//		gdb_ext = tmp.gdb_ext;
 
-		start = tmp.start;
-		stop = tmp.stop;
-		loc_grids = tmp.loc_grids;
+//		start = tmp.start;
+//		stop = tmp.stop;
+//		loc_grids = tmp.loc_grids;
 
-		return *this;
-	}
+//		return *this;
+//	}
 
 	/*! \brief Copy constructor
 	*
 	* \param tmp iterator to copy
 	*
 	*/
-	grid_dist_id_iterator_gpu(const grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & tmp)
-	:loc_grids(tmp.loc_grids)
-	{
-		this->operator=(tmp);
-	}
+//	grid_dist_id_iterator_gpu(const grid_dist_id_iterator_gpu<Decomposition,deviceGrids> & tmp)
+//	:loc_grids(tmp.loc_grids)
+//	{
+//		this->operator=(tmp);
+//	}
 
 	/*! \brief Constructor of the distributed grid iterator
 	 *
@@ -269,6 +269,7 @@ class grid_dist_id_iterator_gpu
 		}
 	}
 
+
 	/*! \brief Get the starting point of the sub-grid we are iterating
 	 *
 	 * \return the starting point
diff --git a/src/Grid/cuda/grid_dist_id_kernels.cuh b/src/Grid/cuda/grid_dist_id_kernels.cuh
index e01fd5199584642bf9a5fa2d4c98a83647985efc..1095325a0cc6cd35151eb5d2d5aeea439bfb2646 100644
--- a/src/Grid/cuda/grid_dist_id_kernels.cuh
+++ b/src/Grid/cuda/grid_dist_id_kernels.cuh
@@ -8,6 +8,8 @@
 #ifndef GRID_DIST_ID_KERNELS_CUH_
 #define GRID_DIST_ID_KERNELS_CUH_
 
+#include "config.h"
+
 #ifdef CUDA_GPU
 
 template<unsigned int dim>
diff --git a/src/Grid/grid_dist_id.hpp b/src/Grid/grid_dist_id.hpp
index 8d4965bd29b1626140afd1cae55cedc4bce99f64..4b5954dbb05438fcddecef3d7766325306c30e70 100644
--- a/src/Grid/grid_dist_id.hpp
+++ b/src/Grid/grid_dist_id.hpp
@@ -1768,6 +1768,32 @@ public:
 
 #ifdef __NVCC__
 
+	/*! \brief set existing points in the grid
+    *
+	* \param f2 lambda function to set points
+	*/
+	template<typename lambda_t2>
+	void setPoints(lambda_t2 f2)
+	{
+		auto it = getGridIteratorGPU();
+
+		it.template launch<1>(launch_set_dense<dim>(),f2);
+	}
+
+	/*! \brief set point existing in the grid between start and stop
+	*
+	* \param start point
+	* \param stop point
+	* \param f2 lambda function to set points
+	*/
+	template<typename lambda_t2>
+	void setPoints(grid_key_dx<dim> k1, grid_key_dx<dim> k2, lambda_t2 f2)
+	{
+		auto it = getGridIteratorGPU(k1,k2);
+
+		it.template launch<0>(launch_set_dense<dim>(),f2);
+	}
+
 	/*! \brief Insert point in the grid
     *
 	* \param f1 lambda function to insert point
@@ -3318,6 +3344,9 @@ using grid_dist_id_devg = grid_dist_id<dim,St,T,Decomposition,Memory,devg>;
 template<unsigned int dim, typename St, typename T, typename Memory = CudaMemory, typename Decomposition = CartDecomposition<dim,St,CudaMemory,memory_traits_inte> >
 using sgrid_dist_id_gpu = grid_dist_id<dim,St,T,Decomposition,Memory,SparseGridGpu<dim,T>>;
 
+template<unsigned int dim, typename St, typename T, typename Memory = CudaMemory, typename Decomposition = CartDecomposition<dim,St,CudaMemory,memory_traits_inte> >
+using grid_dist_id_gpu = grid_dist_id<dim,St,T,Decomposition,Memory,grid_gpu<dim,T>>;
+
 template<unsigned int dim, typename St, typename T, typename Memory = CudaMemory, typename Decomposition = CartDecomposition<dim,St,CudaMemory,memory_traits_inte> >
 using sgrid_dist_sid_gpu = grid_dist_id<dim,St,T,Decomposition,Memory,SparseGridGpu<dim,T,default_edge<dim>::type::value,default_edge<dim>::tb::value,int>>;
 #endif
diff --git a/src/Grid/grid_dist_id_comm.hpp b/src/Grid/grid_dist_id_comm.hpp
index 5508aa8f1d5fc4ad356f6a01b989f4d72a824a61..d9312ef5ea5bdeefb1c1c7b5a7df7731d9d8854b 100644
--- a/src/Grid/grid_dist_id_comm.hpp
+++ b/src/Grid/grid_dist_id_comm.hpp
@@ -1535,7 +1535,7 @@ public:
 		std::vector<size_t> prp_recv;
 
 		// Create an object of preallocated memory for properties
-		ExtPreAlloc<Memory> & prRecv_prp = *(new ExtPreAlloc<Memory>(tot_recv,g_recv_prp_mem));
+		ExtPreAlloc<Memory> & prRecv_prp = *(new ExtPreAlloc<Memory>(0,g_recv_prp_mem));
 		prRecv_prp.incRef();
 
 		queue_recv_data_put<prp_object>(ig_box,prp_recv,prRecv_prp);
diff --git a/src/Grid/performance/grid_dist_performance.hpp b/src/Grid/performance/grid_dist_performance.hpp
index 368c2b9b21e0297ca58a166ed9d3972539f8af33..13b818ad897cc5abeb2c0f127c5f03515d9ab5d5 100644
--- a/src/Grid/performance/grid_dist_performance.hpp
+++ b/src/Grid/performance/grid_dist_performance.hpp
@@ -356,7 +356,7 @@ BOOST_AUTO_TEST_CASE(grid_iterator_performance_write_report_final)
 
 		if (create_vcluster().getProcessUnitID() == 0)
 		{
-			addUpdtateTime(cg,create_vcluster().size());
+			addUpdtateTime(cg,create_vcluster().size(),"data","grid_performance");
 
 			cg.write("grid_performance.html");
 		}
diff --git a/src/Grid/tests/grid_dist_id_unit_test.cpp b/src/Grid/tests/grid_dist_id_unit_test.cpp
index a13d79a3494da2fc565349fc14c921441bbbf807..d34a1201ad0a8b31a235900edf9c3204933b124e 100644
--- a/src/Grid/tests/grid_dist_id_unit_test.cpp
+++ b/src/Grid/tests/grid_dist_id_unit_test.cpp
@@ -189,7 +189,8 @@ BOOST_AUTO_TEST_CASE( grid_dist_id_domain_grid_unit_converter_test)
 
 			Box<2,size_t> g_box = g_dist.getCellDecomposer().convertDomainSpaceIntoGridUnits(sub,bc);
 
-			vol += g_box.getVolumeKey();
+			if (g_box.isValid() == true)
+			{vol += g_box.getVolumeKey();}
 		}
 
 		v_cl.sum(vol);
diff --git a/src/Grid/tests/grid_dist_id_unit_test.cu b/src/Grid/tests/grid_dist_id_unit_test.cu
new file mode 100644
index 0000000000000000000000000000000000000000..01770490120f6cd8986b253533fb357a3898c35a
--- /dev/null
+++ b/src/Grid/tests/grid_dist_id_unit_test.cu
@@ -0,0 +1,57 @@
+#define BOOST_TEST_DYN_LINK
+#include <boost/test/unit_test.hpp>
+
+#include "Point_test.hpp"
+#include "Grid/grid_dist_id.hpp"
+#include "data_type/aggregate.hpp"
+
+extern void print_test_v(std::string test, size_t sz);
+
+BOOST_AUTO_TEST_SUITE( grid_dist_id_test )
+
+
+BOOST_AUTO_TEST_CASE( grid_dist_id_gpu_test )
+{
+	// Test grid periodic
+
+/*	Box<3,float> domain({-1.0,-1.0,-1.0},{1.0,1.0,1.0});
+
+	Vcluster<> & v_cl = create_vcluster();
+
+	if ( v_cl.getProcessingUnits() > 32 )
+	{return;}
+
+	// grid size
+	size_t sz[3];
+	sz[0] = 32;
+	sz[1] = 32;
+	sz[2] = 32;
+
+	// Ghost
+	Ghost<3,long int> g(1);
+
+	// periodicity
+	periodicity<3> pr = {{PERIODIC,PERIODIC,PERIODIC}};
+
+	// Distributed grid with id decomposition
+    grid_dist_id_gpu<3, float, aggregate<float, float>> g_dist(sz,domain,g,pr);
+    
+	Box<3,size_t> box({1,1,1},{30,30,30});
+    auto it = g_dist.getGridIterator(box.getKP1(),box.getKP2());
+
+    float c = 5.0;
+
+    typedef typename GetSetBlockType<decltype(g_dist)>::type BlockT;
+
+	g_dist.setPoints(box.getKP1(),box.getKP2(),
+			        [c] __device__ (BlockT & data, int i, int j, int k)
+			        {
+			        	data.template get<0>() = c + i*i + j*j + k*k;
+			        }
+                    );
+  */                  
+    
+}
+
+
+BOOST_AUTO_TEST_SUITE_END()
diff --git a/src/Grid/tests/sgrid_dist_id_gpu_unit_tests.cu b/src/Grid/tests/sgrid_dist_id_gpu_unit_tests.cu
index e26743dfb8c54d8b0e91ae85fc3c4fe67d5fa636..432cc219110021621867a3e9d9d69f120f3bbce4 100644
--- a/src/Grid/tests/sgrid_dist_id_gpu_unit_tests.cu
+++ b/src/Grid/tests/sgrid_dist_id_gpu_unit_tests.cu
@@ -300,6 +300,69 @@ BOOST_AUTO_TEST_CASE( sgrid_gpu_test_ghost_get )
 	sgrid_ghost_get(sz2,sz4);
 }
 
+BOOST_AUTO_TEST_CASE( sgrid_gpu_app_point_test_no_box )
+{
+	size_t sz[3] = {75,75,75};
+	periodicity<3> bc = {PERIODIC,PERIODIC,PERIODIC};
+
+	Ghost<3,long int> g(1);
+
+	Box<3,float> domain({0.0,0.0,0.0},{1.0,1.0,1.0});
+
+	sgrid_dist_id_gpu<3,float,aggregate<float,float,float,float>> gdist(sz,domain,g,bc);
+
+	gdist.template setBackgroundValue<0>(666);
+	gdist.template setBackgroundValue<1>(666);
+	gdist.template setBackgroundValue<2>(666);
+	gdist.template setBackgroundValue<3>(666);
+
+	/////// GPU insert + flush
+
+	Box<3,size_t> box({1,1,1},{sz[0],sz[1],sz[2]});
+
+	/////// GPU Run kernel
+
+	float c = 5.0;
+
+	typedef typename GetAddBlockType<decltype(gdist)>::type InsertBlockT;
+
+	CudaMemory cmem;
+	cmem.allocate(sizeof(int));
+
+	*(int *)cmem.getPointer() = 0.0;
+
+	cmem.hostToDevice();
+
+	int * cnt = (int *)cmem.getDevicePointer();
+
+	gdist.addPoints([cnt] __device__ (int i, int j, int k)
+			        {
+						atomicAdd(cnt,1);
+
+						return true;
+			        },
+			        [c] __device__ (InsertBlockT & data, int i, int j, int k)
+			        {
+			        	data.template get<0>() = c + i + j;
+			        	data.template get<1>() = c + 1000 + i + j;
+			        }
+			        );
+
+	gdist.template flush<smax_<0>,smax_<1>>(flush_type::FLUSH_ON_DEVICE);
+	gdist.template ghost_get<0,1>(RUN_ON_DEVICE);
+
+	cmem.deviceToHost();
+
+	int cnt_host = *(int *)cmem.getPointer();
+
+	auto & v_cl = create_vcluster();
+
+	v_cl.sum(cnt_host);
+	v_cl.execute();
+
+	BOOST_REQUIRE_EQUAL(cnt_host,75*75*75);
+}
+
 
 BOOST_AUTO_TEST_CASE( sgrid_gpu_app_point_test )
 {
@@ -352,7 +415,6 @@ BOOST_AUTO_TEST_CASE( sgrid_gpu_app_point_test )
 						{atomicAdd(cnt,1);}
 						else
 						{
-							printf("%d %d %d \n",i,j,k);
 							atomicAdd(cnt_out,1);
 						}
 
diff --git a/src/Vector/cuda/vector_dist_gpu_unit_tests.cu b/src/Vector/cuda/vector_dist_gpu_unit_tests.cu
index c8fcec4e295126ecc66e0f7c06a8fe537e68b3ca..9380e1e72cc5c1ca6b32525735cc374f0d623af7 100644
--- a/src/Vector/cuda/vector_dist_gpu_unit_tests.cu
+++ b/src/Vector/cuda/vector_dist_gpu_unit_tests.cu
@@ -1234,7 +1234,7 @@ void vector_dist_dlb_on_cuda_impl_async(size_t k,double r_cut)
 			++it;
 		}
 
-		//Back to GPU
+		// Back to GPU
 		vd.hostToDevicePos();
 		vd.map(RUN_ON_DEVICE);
 		vd.template Ighost_get<0>(RUN_ON_DEVICE);
diff --git a/src/Vector/performance/cell_list_comp_reorder.hpp b/src/Vector/performance/cell_list_comp_reorder.hpp
index b0b01d3ff614fb7a8b5c0dcc14135a2d05bd4244..bbfae08b35e00051c57953481152e2d10092d256 100644
--- a/src/Vector/performance/cell_list_comp_reorder.hpp
+++ b/src/Vector/performance/cell_list_comp_reorder.hpp
@@ -357,7 +357,9 @@ BOOST_AUTO_TEST_CASE(vector_dist_cl_performance_write_report)
 
 		StandardXMLPerformanceGraph("celllist_performance.xml",file_xml_ref,cg);
 
-		addUpdtateTime(cg,create_vcluster().size());
+		addUpdtateTime(cg,create_vcluster().size(),"pdata","celllist_performance");
+
+		createCommitFile("pdata");
 
 		cg.write("celllist_performance.html");
 	}
diff --git a/src/Vector/performance/cell_list_part_reorder.hpp b/src/Vector/performance/cell_list_part_reorder.hpp
index c31c5c3f896582dcf6d56e72294f2c2104fe013f..ed7a252390cc4dc5deef65bfa5a862ff7af40927 100644
--- a/src/Vector/performance/cell_list_part_reorder.hpp
+++ b/src/Vector/performance/cell_list_part_reorder.hpp
@@ -320,7 +320,7 @@ BOOST_AUTO_TEST_CASE(vector_dist_cl_performance_write_report)
 
 		StandardXMLPerformanceGraph("celllist_partreo_performance.xml",file_xml_ref,cg);
 
-		addUpdtateTime(cg,create_vcluster().size());
+		addUpdtateTime(cg,create_vcluster().size(),"pdata","celllist_part_ord");
 
 		if (create_vcluster().getProcessUnitID() == 0)
 		{cg.write("celllist_part_ord.html");}
diff --git a/src/Vector/performance/vector_dist_gg_map_performance.hpp b/src/Vector/performance/vector_dist_gg_map_performance.hpp
index a75ba1b14c3ffb9feb880404dbea641e9317a1a4..bdd72c375f66e90bf02a663c35072ad429b12c46 100644
--- a/src/Vector/performance/vector_dist_gg_map_performance.hpp
+++ b/src/Vector/performance/vector_dist_gg_map_performance.hpp
@@ -169,7 +169,7 @@ BOOST_AUTO_TEST_CASE(vector_dist_gg_map_performance_write_report)
 
 		StandardXMLPerformanceGraph("particles_map_performance.xml",file_xml_ref,cg);
 
-		addUpdtateTime(cg,create_vcluster().size());
+		addUpdtateTime(cg,create_vcluster().size(),"pdata","particles_map_performance");
 
 		if (create_vcluster().getProcessUnitID() == 0)
 		{cg.write("particles_map_performance.html");}
diff --git a/src/Vector/performance/verlet_performance_tests.hpp b/src/Vector/performance/verlet_performance_tests.hpp
index 97058e3418a41296498e09126104e5b2b42831f1..3c6bf7462fefde745050b78e4fcd2e0c52e572c7 100644
--- a/src/Vector/performance/verlet_performance_tests.hpp
+++ b/src/Vector/performance/verlet_performance_tests.hpp
@@ -217,7 +217,7 @@ BOOST_AUTO_TEST_CASE(vector_dist_verlet_performance_write_report)
 
 		StandardXMLPerformanceGraph("verletlist_performance.xml",file_xml_ref,cg);
 
-		addUpdtateTime(cg,create_vcluster().size());
+		addUpdtateTime(cg,create_vcluster().size(),"pdata","verletlist_performance");
 
 		if (create_vcluster().getProcessUnitID() == 0)
 		{cg.write("verletlist_performance.html");}
diff --git a/src/Vector/vector_dist.hpp b/src/Vector/vector_dist.hpp
index 66ac010b09e1039abcd0c19974cb23f9c1c93a05..5f740f5c78e4fb9223c8ca318be602bfb1f38796 100644
--- a/src/Vector/vector_dist.hpp
+++ b/src/Vector/vector_dist.hpp
@@ -260,7 +260,9 @@ template<unsigned int dim,
          typename prop,
          typename Decomposition = CartDecomposition<dim,St>,
          typename Memory = HeapMemory,
-         template<typename> class layout_base = memory_traits_lin>
+         template<typename> class layout_base = memory_traits_lin,
+		 typename vector_dist_pos = openfpm::vector<Point<dim, St>,Memory,layout_base>,
+		 typename vector_dist_prop = openfpm::vector<prop,Memory,layout_base> >
 class vector_dist : public vector_dist_comm<dim,St,prop,Decomposition,Memory,layout_base>,
 #ifdef CUDA_GPU
 					private vector_dist_ker_list<vector_dist_ker<dim,St,prop,layout_base>>
@@ -281,17 +283,17 @@ private:
 
 	//! Particle position vector, (It has 2 elements) the first has real particles assigned to a processor
 	//! the second element contain unassigned particles
-	openfpm::vector<Point<dim, St>,Memory,layout_base> v_pos;
+	vector_dist_pos v_pos;
 
 	//! Particle properties vector, (It has 2 elements) the first has real particles assigned to a processor
 	//! the second element contain unassigned particles
-	openfpm::vector<prop,Memory,layout_base> v_prp;
+	vector_dist_prop v_prp;
 
 	//! reordered v_pos buffer
-	openfpm::vector<prop,Memory,layout_base> v_prp_out;
+	vector_dist_prop v_prp_out;
 
 	//! reordered v_prp buffer
-	openfpm::vector<Point<dim, St>,Memory,layout_base> v_pos_out;
+	vector_dist_pos v_pos_out;
 
 	//! option used to create this vector
 	size_t opt = 0;
@@ -2528,6 +2530,22 @@ public:
 		g_m -= keys.size();
 	}
 
+	/*! \brief Remove a set of elements from the distributed vector
+	 *
+	 * \warning keys must be sorted
+	 *
+	 * \param keys vector of elements to eliminate
+	 * \param start from where to eliminate
+	 *
+	 */
+	void remove(openfpm::vector<aggregate<int>> & keys, size_t start = 0)
+	{
+		v_pos.remove(keys, start);
+		v_prp.remove(keys, start);
+
+		g_m -= keys.size();
+	}
+
 	/*! \brief Remove one element from the distributed vector
 	 *
 	 * \param key remove one element from the vector
@@ -2687,8 +2705,8 @@ public:
 		if ((opt & 0x0FFF0000) == CSV_WRITER)
 		{
 			// CSVWriter test
-			CSVWriter<openfpm::vector<Point<dim, St>,Memory,layout_base>,
-			          openfpm::vector<prop,Memory,layout_base> > csv_writer;
+			CSVWriter<vector_dist_pos,
+			          vector_dist_prop > csv_writer;
 
 			std::string output = std::to_string(out + "_" + std::to_string(v_cl.getProcessUnitID()) + std::to_string(".csv"));
 
@@ -2703,8 +2721,8 @@ public:
 				ft = file_type::BINARY;
 
 			// VTKWriter for a set of points
-			VTKWriter<boost::mpl::pair<openfpm::vector<Point<dim, St>,Memory,layout_base>,
-									   openfpm::vector<prop,Memory,layout_base>>,
+			VTKWriter<boost::mpl::pair<vector_dist_pos,
+									   vector_dist_prop>,
 			                           VECTOR_POINTS> vtk_writer;
 			vtk_writer.add(v_pos,v_prp,g_m);
 
@@ -2783,8 +2801,8 @@ public:
 		if ((opt & 0x0FFF0000) == CSV_WRITER)
 		{
 			// CSVWriter test
-			CSVWriter<openfpm::vector<Point<dim, St>,Memory,layout_base>,
-					  openfpm::vector<prop,Memory,layout_base> > csv_writer;
+			CSVWriter<vector_dist_pos,
+					  vector_dist_prop > csv_writer;
 
 			std::string output = std::to_string(out + "_" + std::to_string(v_cl.getProcessUnitID()) + "_" + std::to_string(iteration) + std::to_string(".csv"));
 
@@ -2799,8 +2817,8 @@ public:
 				ft = file_type::BINARY;
 
 			// VTKWriter for a set of points
-			VTKWriter<boost::mpl::pair<openfpm::vector<Point<dim, St>,Memory,layout_base>,
-									   openfpm::vector<prop,Memory,layout_base>>, VECTOR_POINTS> vtk_writer;
+			VTKWriter<boost::mpl::pair<vector_dist_pos,
+									   vector_dist_prop>, VECTOR_POINTS> vtk_writer;
 			vtk_writer.add(v_pos,v_prp,g_m);
 
 			std::string output = std::to_string(out + "_" + std::to_string(v_cl.getProcessUnitID()) + "_" + std::to_string(iteration) + std::to_string(".vtp"));
@@ -2872,7 +2890,7 @@ public:
 	 * \return the particle position vector
 	 *
 	 */
-	const openfpm::vector<Point<dim, St>,Memory,layout_base> & getPosVector() const
+	const vector_dist_pos & getPosVector() const
 	{
 		return v_pos;
 	}
@@ -2882,7 +2900,7 @@ public:
 	 * \return the particle position vector
 	 *
 	 */
-	openfpm::vector<Point<dim, St>,Memory,layout_base> & getPosVector()
+	vector_dist_pos & getPosVector()
 	{
 		return v_pos;
 	}
@@ -2892,7 +2910,7 @@ public:
 	 * \return the particle property vector
 	 *
 	 */
-	const openfpm::vector<prop,Memory,layout_base> & getPropVector() const
+	const vector_dist_prop & getPropVector() const
 	{
 		return v_prp;
 	}
@@ -2902,7 +2920,7 @@ public:
 	 * \return the particle property vector
 	 *
 	 */
-	openfpm::vector<prop,Memory,layout_base> & getPropVector()
+	vector_dist_prop & getPropVector()
 	{
 		return v_prp;
 	}
@@ -2912,7 +2930,7 @@ public:
 	 * \return the particle position vector
 	 *
 	 */
-	const openfpm::vector<Point<dim, St>,Memory,layout_base> & getPosVectorSort() const
+	const vector_dist_pos & getPosVectorSort() const
 	{
 		return v_pos_out;
 	}
@@ -2922,7 +2940,7 @@ public:
 	 * \return the particle position vector
 	 *
 	 */
-	openfpm::vector<Point<dim, St>,Memory,layout_base> & getPosVectorSort()
+	vector_dist_pos & getPosVectorSort()
 	{
 		return v_pos_out;
 	}
@@ -2932,7 +2950,7 @@ public:
 	 * \return the particle property vector
 	 *
 	 */
-	const openfpm::vector<prop,Memory,layout_base> & getPropVectorSort() const
+	const vector_dist_prop & getPropVectorSort() const
 	{
 		return v_prp_out;
 	}
@@ -2942,7 +2960,7 @@ public:
 	 * \return the particle property vector
 	 *
 	 */
-	openfpm::vector<prop,Memory,layout_base> & getPropVectorSort()
+	vector_dist_prop & getPropVectorSort()
 	{
 		return v_prp_out;
 	}
diff --git a/src/lib/pdata_python.cpp b/src/lib/pdata_python.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..03710408b3097419999e51ab75695a25218b24f6
--- /dev/null
+++ b/src/lib/pdata_python.cpp
@@ -0,0 +1,7 @@
+
+void init_openfpm_python()
+{
+}
+
+
+
diff --git a/src/scripts/postinst b/src/scripts/postinst
index 763d3ed561a4529800f9222de37392223c86c1b1..f11b3ca9b487dc465991be39a7e9af67b48e813a 100644
--- a/src/scripts/postinst
+++ b/src/scripts/postinst
@@ -4,6 +4,5 @@
 
 sed -i -e 's/projects\/ppm\/rundeck\/openfpm_super_bundles\/$1\/openfpm_dep_$1/usr\/local\/openfpm\/dependencies/g' /usr/local/openfpm/source/openfpm_vars
 sed -i -e 's/projects\/ppm\/rundeck\/openfpm_super_bundles\/$1\/openfpm_dep_$1/usr\/local\/openfpm\/dependencies/g' /usr/local/openfpm/openfpm_pdata/include/example.mk
-echo "export OPAL_PREFIX=/usr/local/openfpm/dependencies/MPI" >> /usr/local/openfpm/source/openfpm_vars
-
 
+echo "export OPAL_PREFIX=/usr/local/openfpm/dependencies/MPI" >> /usr/local/openfpm/source/openfpm_vars
diff --git a/src/test_multiple_o.cpp b/src/test_multiple_o.cpp
deleted file mode 100644
index 5e45e1d284739bb0f6820b12744ba2de72cd08a2..0000000000000000000000000000000000000000
--- a/src/test_multiple_o.cpp
+++ /dev/null
@@ -1,35 +0,0 @@
-/*
- * test_multiple_o.cpp
- *
- *  Created on: Feb 5, 2016
- *      Author: i-bird
- *
- *
- *  It just test that the compilation with multiple translation unit (*.o) does not
- *  produce error, if we have duplicated symbol in the translation unit we will get error
- *
- */
-
-#include "Vector/vector_dist.hpp"
-#include "Grid/grid_dist_id.hpp"
-#include "data_type/aggregate.hpp"
-#include "Decomposition/CartDecomposition.hpp"
-
-void f()
-{
-	// Ghost
-	Ghost<3,float> g(0.01);
-
-	size_t bc[3]={PERIODIC,PERIODIC,PERIODIC};
-
-	Box<3,float> domain({0.0,0.0,0.0},{1.0,1.0,1.0});
-	size_t sz[3];
-	sz[0] = 100;
-	sz[1] = 100;
-	sz[2] = 100;
-
-	vector_dist<3,float, aggregate<float>, CartDecomposition<3,float> > vd(4096,domain,bc,g);
-	grid_dist_id<3, float, aggregate<float[3]>, CartDecomposition<3,float>> g_dist(sz,domain,g);
-}
-
-