diff --git a/Makefile b/Makefile
index e8cf23e8e4de6324b1ca97bfac85b18020b4d659..384f6fb5bd9a28619e71fac66f8f123c378ba180 100644
--- a/Makefile
+++ b/Makefile
@@ -3,14 +3,16 @@
 TESTCASE=lj
 PYCMD=python3
 #CC=icc
-CC=mpiicc
+#CC=mpiicpx
+CC=mpiicpc
 NVCC=nvcc
 NVCC_PATH:="$(shell which ${NVCC})"
 LIKWID_INC ?= -I/usr/local/include
 LIKWID_DEFINES ?= -DLIKWID_PERFMON
 LIKWID_LIB ?= -L/usr/local/lib
 LIKWID_FLAGS = -llikwid ${LIKWID_INC} ${LIKWID_DEFINES} ${LIKWID_LIB}
-CFLAGS=-Ofast -xCORE-AVX512 -qopt-zmm-usage=high ${LIKWID_FLAGS}
+CFLAGS=-Ofast -xHost -qopt-zmm-usage=high ${LIKWID_FLAGS}
+#CFLAGS=-Ofast -xCORE-AVX512 -qopt-zmm-usage=high ${LIKWID_FLAGS}
 CUDA_BIN_PATH:="$(shell dirname ${NVCC_PATH})"
 CUDA_PATH:="$(shell dirname ${CUDA_BIN_PATH})"
 OBJ_PATH=obj
@@ -19,7 +21,7 @@ CPU_BIN="$(TESTCASE)_cpu"
 GPU_SRC="$(TESTCASE).cu"
 GPU_BIN="$(TESTCASE)_gpu"
 DEBUG_FLAGS=
-#DEBUG_FLAGS="-DDEBUG"
+#DEBUG_FLAGS=-DDEBUG
 
 all: clean build $(CPU_BIN) $(GPU_BIN)
 	@echo "Everything was done!"
@@ -54,7 +56,7 @@ $(CPU_BIN): $(CPU_SRC) $(OBJ_PATH)/pairs.o $(OBJ_PATH)/regular_6d_stencil.o $(OB
 $(GPU_BIN): $(GPU_SRC) $(OBJ_PATH)/pairs.o $(OBJ_PATH)/regular_6d_stencil.o 
 	$(NVCC) -c -o $(OBJ_PATH)/cuda_runtime.o runtime/devices/cuda.cu ${DEBUG_FLAGS}
 	$(NVCC) -c -o $(OBJ_PATH)/$(GPU_BIN).o $(GPU_SRC) -DDEBUG
-	$(CC) -o $(CFLAGS) $(GPU_BIN) $(OBJ_PATH)/$(GPU_BIN).o $(OBJ_PATH)/cuda_runtime.o $(OBJ_PATH)/pairs.o $(OBJ_PATH)/regular_6d_stencil.o -lcudart -L$(CUDA_PATH)/lib64
+	$(CC) -o $(GPU_BIN) $(CFLAGS) $(OBJ_PATH)/$(GPU_BIN).o $(OBJ_PATH)/cuda_runtime.o $(OBJ_PATH)/pairs.o $(OBJ_PATH)/regular_6d_stencil.o -lcudart -L$(CUDA_PATH)/lib64
 
 clean:
 	@echo "Cleaning..."
diff --git a/runtime/contact_property.hpp b/runtime/contact_property.hpp
index 04950c83ee6b08ed0e7f3121757d88150b91335a..2d1e03d66362b73a365e5bebcae025680e190dfe 100644
--- a/runtime/contact_property.hpp
+++ b/runtime/contact_property.hpp
@@ -32,10 +32,10 @@ public:
     void *getDevicePointer() { return d_ptr; }
     void setPointers(void *h_ptr_, void *d_ptr_) { h_ptr = h_ptr_, d_ptr = d_ptr_; }
     void setSizes(size_t sx_, size_t sy_) { sx = sx_, sy = sy_; }
-    size_t getTotalSize() { return sx * sy * getElemSize(); };
+    size_t getTotalSize() { return sx * sy * getPrimitiveTypeSize(); };
     PropertyType getType() { return type; }
     layout_t getLayout() { return layout; }
-    size_t getElemSize() {
+    size_t getPrimitiveTypeSize() {
         return  (type == Prop_Integer) ? sizeof(int) :
                 (type == Prop_Real) ? sizeof(real_t) :
                 (type == Prop_Vector) ? sizeof(real_t) :
diff --git a/runtime/pairs.cpp b/runtime/pairs.cpp
index 42a4d257cd0808bb1fe3f36eabd6d492b070b43f..bb18af189a7bf0137d4bf15f9fbeb8ddb512801a 100644
--- a/runtime/pairs.cpp
+++ b/runtime/pairs.cpp
@@ -11,7 +11,10 @@
 
 namespace pairs {
 
-void PairsSimulation::initDomain(int *argc, char ***argv, real_t xmin, real_t xmax, real_t ymin, real_t ymax, real_t zmin, real_t zmax) {
+void PairsSimulation::initDomain(
+    int *argc, char ***argv,
+    real_t xmin, real_t xmax, real_t ymin, real_t ymax, real_t zmin, real_t zmax) {
+
     if(dom_part_type == DimRanges) {
         dom_part = new Regular6DStencil(xmin, xmax, ymin, ymax, zmin, zmax);
     } else {
@@ -23,72 +26,114 @@ void PairsSimulation::initDomain(int *argc, char ***argv, real_t xmin, real_t xm
 
 void PairsSimulation::addArray(Array array) {
     int id = array.getId();
-    auto a = std::find_if(arrays.begin(), arrays.end(), [id](Array _a) { return _a.getId() == id; });
+    auto a = std::find_if(
+        arrays.begin(),
+        arrays.end(),
+        [id](Array _a) { return _a.getId() == id; });
+
     PAIRS_ASSERT(a == std::end(arrays));
     arrays.push_back(array);
 }
 
 Array &PairsSimulation::getArray(array_t id) {
-    auto a = std::find_if(arrays.begin(), arrays.end(), [id](Array _a) { return _a.getId() == id; });
+    auto a = std::find_if(
+        arrays.begin(),
+        arrays.end(),
+        [id](Array _a) { return _a.getId() == id; });
+
     PAIRS_ASSERT(a != std::end(arrays));
     return *a;
 }
 
 Array &PairsSimulation::getArrayByName(std::string name) {
-    auto a = std::find_if(arrays.begin(), arrays.end(), [name](Array _a) { return _a.getName() == name; });
+    auto a = std::find_if(
+        arrays.begin(),
+        arrays.end(),
+        [name](Array _a) { return _a.getName() == name; });
+
     PAIRS_ASSERT(a != std::end(arrays));
     return *a;
 }
 
 Array &PairsSimulation::getArrayByHostPointer(const void *h_ptr) {
-    auto a = std::find_if(arrays.begin(), arrays.end(), [h_ptr](Array _a) { return _a.getHostPointer() == h_ptr; });
+    auto a = std::find_if(
+        arrays.begin(),
+        arrays.end(),
+        [h_ptr](Array _a) { return _a.getHostPointer() == h_ptr; });
+
     PAIRS_ASSERT(a != std::end(arrays));
     return *a;
 }
 
 void PairsSimulation::addProperty(Property prop) {
     int id = prop.getId();
-    auto p = std::find_if(properties.begin(), properties.end(), [id](Property _p) { return _p.getId() == id; });
+    auto p = std::find_if(
+        properties.begin(),
+        properties.end(),
+        [id](Property _p) { return _p.getId() == id; });
+
     PAIRS_ASSERT(p == std::end(properties));
     properties.push_back(prop);
 }
 
 Property &PairsSimulation::getProperty(property_t id) {
-    auto p = std::find_if(properties.begin(), properties.end(), [id](Property _p) { return _p.getId() == id; });
+    auto p = std::find_if(
+        properties.begin(),
+        properties.end(),
+        [id](Property _p) { return _p.getId() == id; });
+
     PAIRS_ASSERT(p != std::end(properties));
     return *p;
 }
 
 Property &PairsSimulation::getPropertyByName(std::string name) {
-    auto p = std::find_if(properties.begin(), properties.end(), [name](Property _p) { return _p.getName() == name; });
+    auto p = std::find_if(
+        properties.begin(),
+        properties.end(),
+        [name](Property _p) { return _p.getName() == name; });
+
     PAIRS_ASSERT(p != std::end(properties));
     return *p;
 }
 
 void PairsSimulation::addContactProperty(ContactProperty contact_prop) {
     int id = contact_prop.getId();
-    auto cp = std::find_if(contact_properties.begin(), contact_properties.end(), [id](ContactProperty _cp) { return _cp.getId() == id; });
+    auto cp = std::find_if(
+        contact_properties.begin(),
+        contact_properties.end(),
+        [id](ContactProperty _cp) { return _cp.getId() == id; });
+
     PAIRS_ASSERT(cp == std::end(contact_properties));
     contact_properties.push_back(contact_prop);
 }
 
 ContactProperty &PairsSimulation::getContactProperty(property_t id) {
-    auto cp = std::find_if(contact_properties.begin(), contact_properties.end(), [id](ContactProperty _cp) { return _cp.getId() == id; });
+    auto cp = std::find_if(
+        contact_properties.begin(),
+        contact_properties.end(),
+        [id](ContactProperty _cp) { return _cp.getId() == id; });
+
     PAIRS_ASSERT(cp != std::end(contact_properties));
     return *cp;
 }
 
 ContactProperty &PairsSimulation::getContactPropertyByName(std::string name) {
-    auto cp = std::find_if(contact_properties.begin(), contact_properties.end(), [name](ContactProperty _cp) { return _cp.getName() == name; });
+    auto cp = std::find_if(
+        contact_properties.begin(),
+        contact_properties.end(),
+        [name](ContactProperty _cp) { return _cp.getName() == name; });
+
     PAIRS_ASSERT(cp != std::end(contact_properties));
     return *cp;
 }
 
 void PairsSimulation::addFeatureProperty(FeatureProperty feature_prop) {
     int id = feature_prop.getId();
-    auto fp = std::find_if(feature_properties.begin(),
-                           feature_properties.end(),
-                           [id](FeatureProperty _fp) { return _fp.getId() == id; });
+    auto fp = std::find_if(
+        feature_properties.begin(),
+        feature_properties.end(),
+        [id](FeatureProperty _fp) { return _fp.getId() == id; });
+
     PAIRS_ASSERT(fp == std::end(feature_properties));
     feature_properties.push_back(feature_prop);
 }
@@ -109,7 +154,7 @@ FeatureProperty &PairsSimulation::getFeaturePropertyByName(std::string name) {
     return *fp;
 }
 
-void PairsSimulation::copyArrayToDevice(Array &array, action_t action) {
+void PairsSimulation::copyArrayToDevice(Array &array, action_t action, size_t size) {
     int array_id = array.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
@@ -120,7 +165,10 @@ void PairsSimulation::copyArrayToDevice(Array &array, action_t action) {
                     array.getHostPointer(), array.getDevicePointer(), array.getSize());
             } else {
                 PAIRS_DEBUG("Copying array %s to device\n", array.getName().c_str());
-                pairs::copy_to_device(array.getHostPointer(), array.getDevicePointer(), array.getSize());
+                pairs::copy_to_device(
+                    array.getHostPointer(),
+                    array.getDevicePointer(),
+                    (size > 0) ? size : array.getSize());
             }
         }
     }
@@ -132,7 +180,7 @@ void PairsSimulation::copyArrayToDevice(Array &array, action_t action) {
     array_flags->setDeviceFlag(array_id);
 }
 
-void PairsSimulation::copyArrayToHost(Array &array, action_t action) {
+void PairsSimulation::copyArrayToHost(Array &array, action_t action, size_t size) {
     int array_id = array.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
@@ -143,7 +191,10 @@ void PairsSimulation::copyArrayToHost(Array &array, action_t action) {
                     array.getDevicePointer(), array.getHostPointer(), array.getSize());
             } else {
                 PAIRS_DEBUG("Copying array %s to host\n", array.getName().c_str());
-                pairs::copy_to_host(array.getDevicePointer(), array.getHostPointer(), array.getSize());
+                pairs::copy_to_host(
+                    array.getDevicePointer(),
+                    array.getHostPointer(),
+                    (size > 0) ? size : array.getSize());
             }
 
         }
@@ -156,13 +207,16 @@ void PairsSimulation::copyArrayToHost(Array &array, action_t action) {
     array_flags->setHostFlag(array_id);
 }
 
-void PairsSimulation::copyPropertyToDevice(Property &prop, action_t action) {
+void PairsSimulation::copyPropertyToDevice(Property &prop, action_t action, size_t size) {
     int prop_id = prop.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
         if(action == Ignore || !prop_flags->isDeviceFlagSet(prop_id)) {
             PAIRS_DEBUG("Copying property %s to device\n", prop.getName().c_str());
-            pairs::copy_to_device(prop.getHostPointer(), prop.getDevicePointer(), prop.getTotalSize());
+            pairs::copy_to_device(
+                prop.getHostPointer(),
+                prop.getDevicePointer(),
+                (size > 0) ? size : prop.getTotalSize());
         }
     }
 
@@ -173,13 +227,16 @@ void PairsSimulation::copyPropertyToDevice(Property &prop, action_t action) {
     prop_flags->setDeviceFlag(prop_id);
 }
 
-void PairsSimulation::copyPropertyToHost(Property &prop, action_t action) {
+void PairsSimulation::copyPropertyToHost(Property &prop, action_t action, size_t size) {
     int prop_id = prop.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
         if(action == Ignore || !prop_flags->isHostFlagSet(prop_id)) {
             PAIRS_DEBUG("Copying property %s to host\n", prop.getName().c_str());
-            pairs::copy_to_host(prop.getDevicePointer(), prop.getHostPointer(), prop.getTotalSize());
+            pairs::copy_to_host(
+                prop.getDevicePointer(),
+                prop.getHostPointer(),
+                (size > 0) ? size : prop.getTotalSize());
         }
     }
 
@@ -190,7 +247,9 @@ void PairsSimulation::copyPropertyToHost(Property &prop, action_t action) {
     prop_flags->setHostFlag(prop_id);
 }
 
-void PairsSimulation::copyContactPropertyToDevice(ContactProperty &contact_prop, action_t action) {
+void PairsSimulation::copyContactPropertyToDevice(
+    ContactProperty &contact_prop, action_t action, size_t size) {
+
     int prop_id = contact_prop.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
@@ -199,7 +258,7 @@ void PairsSimulation::copyContactPropertyToDevice(ContactProperty &contact_prop,
             pairs::copy_to_device(
                 contact_prop.getHostPointer(),
                 contact_prop.getDevicePointer(),
-                contact_prop.getTotalSize());
+                (size > 0) ? size : contact_prop.getTotalSize());
 
             contact_prop_flags->setDeviceFlag(prop_id);
         }
@@ -210,7 +269,9 @@ void PairsSimulation::copyContactPropertyToDevice(ContactProperty &contact_prop,
     }
 }
 
-void PairsSimulation::copyContactPropertyToHost(ContactProperty &contact_prop, action_t action) {
+void PairsSimulation::copyContactPropertyToHost(
+    ContactProperty &contact_prop, action_t action, size_t size) {
+
     int prop_id = contact_prop.getId();
 
     if(action == Ignore || action == WriteAfterRead || action == ReadOnly) {
@@ -219,7 +280,7 @@ void PairsSimulation::copyContactPropertyToHost(ContactProperty &contact_prop, a
             pairs::copy_to_host(
                 contact_prop.getDevicePointer(),
                 contact_prop.getHostPointer(),
-                contact_prop.getTotalSize());
+                (size > 0) ? size : contact_prop.getTotalSize());
 
             contact_prop_flags->setHostFlag(prop_id);
         }
@@ -232,7 +293,8 @@ void PairsSimulation::copyContactPropertyToHost(ContactProperty &contact_prop, a
 
 void PairsSimulation::copyFeaturePropertyToDevice(FeatureProperty &feature_prop) {
     PAIRS_DEBUG("Copying static array %s to device\n", feature_prop.getName().c_str());
-    pairs::copy_static_symbol_to_device(feature_prop.getHostPointer(), feature_prop.getDevicePointer(), feature_prop.getArraySize());
+    pairs::copy_static_symbol_to_device(
+        feature_prop.getHostPointer(), feature_prop.getDevicePointer(), feature_prop.getArraySize());
 }
 
 void PairsSimulation::communicateSizes(int dim, const int *send_sizes, int *recv_sizes) {
@@ -258,14 +320,32 @@ void PairsSimulation::communicateData(
     auto nsend_id = getArrayByHostPointer(nsend).getId();
     auto nrecv_id = getArrayByHostPointer(nrecv).getId();
 
-    copyArrayToHost(send_buf_id, ReadOnly);
     copyArrayToHost(send_offsets_id, ReadOnly);
     copyArrayToHost(recv_offsets_id, ReadOnly);
     copyArrayToHost(nsend_id, ReadOnly);
     copyArrayToHost(nrecv_id, ReadOnly);
+
+    /*
+    int nsend_all = 0;
+    int nrecv_all = 0;
+    for(int d = 0; d <= dim; d++) {
+        nsend_all += nsend[d * 2 + 0];
+        nsend_all += nsend[d * 2 + 1];
+        nrecv_all += nrecv[d * 2 + 0];
+        nrecv_all += nrecv[d * 2 + 1];
+    }
+    */
+
+    copyArrayToHost(send_buf_id, ReadOnly);
+    //copyArrayToHost(send_buf_id, Ignore, nsend_all * elem_size * sizeof(real_t));
+
     array_flags->setHostFlag(recv_buf_id);
     array_flags->clearDeviceFlag(recv_buf_id);
-    this->getDomainPartitioner()->communicateData(dim, elem_size, send_buf, send_offsets, nsend, recv_buf, recv_offsets, nrecv);
+    this->getDomainPartitioner()->communicateData(
+        dim, elem_size, send_buf, send_offsets, nsend, recv_buf, recv_offsets, nrecv);
+
+    //copyArrayToHost(recv_buf_id, Ignore);
+    //copyArrayToHost(recv_buf_id, Ignore, nrecv_all * elem_size * sizeof(real_t));
 
     /*
     // Debug messages
diff --git a/runtime/pairs.hpp b/runtime/pairs.hpp
index 6ad84e74adabc3dc84f753689b3f651a3cf686c5..ba987e2ad510c6c67e164e7bae4349dc8120cdba 100644
--- a/runtime/pairs.hpp
+++ b/runtime/pairs.hpp
@@ -53,142 +53,168 @@ public:
         delete timers;
     }
 
-    void initDomain(int *argc, char ***argv, real_t xmin, real_t xmax, real_t ymin, real_t ymax, real_t zmin, real_t zmax);
-    Regular6DStencil *getDomainPartitioner() { return dom_part; }
-    Timers<double> *getTimers() { return timers; }
-
+    // Variables
     template<typename T>
     RuntimeVar<T> addDeviceVariable(T *h_ptr) {
        return RuntimeVar<T>(h_ptr); 
     }
 
-    template<typename T_ptr> void addArray(array_t id, std::string name, T_ptr **h_ptr, std::nullptr_t, size_t size);
-    template<typename T_ptr> void addArray(array_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr, size_t size);
-    template<typename T_ptr> void addStaticArray(array_t id, std::string name, T_ptr *h_ptr, std::nullptr_t, size_t size);
-    template<typename T_ptr> void addStaticArray(array_t id, std::string name, T_ptr *h_ptr, T_ptr *d_ptr, size_t size);
-    void addArray(Array array);
-
-    template<typename T_ptr> void reallocArray(array_t id, T_ptr **h_ptr, std::nullptr_t, size_t size);
-    template<typename T_ptr> void reallocArray(array_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t size);
-
+    // Arrays
     Array &getArray(array_t id);
     Array &getArrayByName(std::string name);
     Array &getArrayByHostPointer(const void *h_ptr);
+    void addArray(Array array);
 
-    template<typename T_ptr> void addProperty(
-        property_t id, std::string name, T_ptr **h_ptr, std::nullptr_t, PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
-    template<typename T_ptr> void addProperty(
-        property_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr, PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
-    void addProperty(Property prop);
+    template<typename T_ptr>
+    void addArray(array_t id, std::string name, T_ptr **h_ptr, std::nullptr_t, size_t size);
+
+    template<typename T_ptr>
+    void addArray(array_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr, size_t size);
+
+    template<typename T_ptr>
+    void addStaticArray(array_t id, std::string name, T_ptr *h_ptr, std::nullptr_t, size_t size);
+
+    template<typename T_ptr>
+    void addStaticArray(array_t id, std::string name, T_ptr *h_ptr, T_ptr *d_ptr, size_t size);
+
+    template<typename T_ptr>
+    void reallocArray(array_t id, T_ptr **h_ptr, std::nullptr_t, size_t size);
 
-    template<typename T_ptr> void reallocProperty(property_t id, T_ptr **h_ptr, std::nullptr_t, size_t sx = 1, size_t sy = 1);
-    template<typename T_ptr> void reallocProperty(property_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t sx = 1, size_t sy = 1);
+    template<typename T_ptr>
+    void reallocArray(array_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t size);
 
+    void copyArrayToDevice(array_t id, action_t action, size_t size = 0) {
+        copyArrayToDevice(getArray(id), action, size);
+    }
+
+    void copyArrayToDevice(Array &array, action_t action, size_t size = 0);
+    void copyArrayToHost(array_t id, action_t action, size_t size = 0) {
+        copyArrayToHost(getArray(id), action, size);
+    }
+
+    void copyArrayToHost(Array &array, action_t action, size_t size = 0);
+
+    // Properties
     Property &getProperty(property_t id);
     Property &getPropertyByName(std::string name);
+    void addProperty(Property prop);
 
-    inline IntProperty &getAsIntegerProperty(Property &prop) { return static_cast<IntProperty&>(prop); }
-    inline FloatProperty &getAsFloatProperty(Property &prop) { return static_cast<FloatProperty&>(prop); }
-    inline VectorProperty &getAsVectorProperty(Property &prop) { return static_cast<VectorProperty&>(prop); }
-    inline MatrixProperty &getAsMatrixProperty(Property &prop) { return static_cast<MatrixProperty&>(prop); }
-    inline QuaternionProperty &getAsQuaternionProperty(Property &prop) { return static_cast<QuaternionProperty&>(prop); }
-
-    inline IntProperty &getIntegerProperty(property_t property) { return static_cast<IntProperty&>(getProperty(property)); }
-    inline FloatProperty &getFloatProperty(property_t property) { return static_cast<FloatProperty&>(getProperty(property)); }
-    inline VectorProperty &getVectorProperty(property_t property) { return static_cast<VectorProperty&>(getProperty(property)); }
-    inline MatrixProperty &getMatrixProperty(property_t property) { return static_cast<MatrixProperty&>(getProperty(property)); }
-    inline QuaternionProperty &getQuaternionProperty(property_t property) { return static_cast<QuaternionProperty&>(getProperty(property)); }
-
-    template<typename T_ptr> void addContactProperty(
-        property_t id, std::string name, T_ptr **h_ptr, std::nullptr_t, PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
-    template<typename T_ptr> void addContactProperty(
-        property_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr, PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
-    void addContactProperty(ContactProperty prop);
+    template<typename T_ptr>
+    void addProperty(
+        property_t id, std::string name, T_ptr **h_ptr, std::nullptr_t,
+        PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
 
-    template<typename T_ptr> void reallocContactProperty(property_t id, T_ptr **h_ptr, std::nullptr_t, size_t sx = 1, size_t sy = 1);
-    template<typename T_ptr> void reallocContactProperty(property_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t sx = 1, size_t sy = 1);
+    template<typename T_ptr> void addProperty(
+        property_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr,
+        PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
 
-    ContactProperty &getContactProperty(property_t id);
-    ContactProperty &getContactPropertyByName(std::string name);
+    template<typename T_ptr>
+    void reallocProperty(property_t id, T_ptr **h_ptr, std::nullptr_t, size_t sx = 1, size_t sy = 1);
 
-    template<typename T_ptr> void addFeatureProperty(property_t id, std::string name, T_ptr *h_ptr, std::nullptr_t, PropertyType type, int nkinds, int array_size);
-    template<typename T_ptr> void addFeatureProperty(property_t id, std::string name, T_ptr *h_ptr, T_ptr *d_ptr, PropertyType type, int nkinds, int array_size);
-    void addFeatureProperty(FeatureProperty feature_prop);
+    template<typename T_ptr>
+    void reallocProperty(property_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t sx = 1, size_t sy = 1);
 
-    FeatureProperty &getFeatureProperty(property_t id);
-    FeatureProperty &getFeaturePropertyByName(std::string name);
+    inline IntProperty &getAsIntegerProperty(Property &prop) {
+        return static_cast<IntProperty&>(prop);
+    }
 
-    void setArrayDeviceFlag(array_t id) { setArrayDeviceFlag(getArray(id)); }
-    void setArrayDeviceFlag(Array &array) { array_flags->setDeviceFlag(array.getId()); }
-    void clearArrayDeviceFlag(array_t id) { clearArrayDeviceFlag(getArray(id)); }
-    void clearArrayDeviceFlag(Array &array) { array_flags->clearDeviceFlag(array.getId()); }
-    void copyArrayToDevice(array_t id, action_t action) { copyArrayToDevice(getArray(id), action); }
-    void copyArrayToDevice(Array &array, action_t action);
-
-    void setArrayHostFlag(array_t id) { setArrayHostFlag(getArray(id)); }
-    void setArrayHostFlag(Array &array) { array_flags->setHostFlag(array.getId()); }
-    void clearArrayHostFlag(array_t id) { clearArrayHostFlag(getArray(id)); }
-    void clearArrayHostFlag(Array &array) { array_flags->clearHostFlag(array.getId()); }
-    void copyArrayToHost(array_t id, action_t action) { copyArrayToHost(getArray(id), action); }
-    void copyArrayToHost(Array &array, action_t action);
-
-    void setPropertyDeviceFlag(property_t id) { setPropertyDeviceFlag(getProperty(id)); }
-    void setPropertyDeviceFlag(Property &prop) { prop_flags->setDeviceFlag(prop.getId()); }
-    void clearPropertyDeviceFlag(property_t id) { clearPropertyDeviceFlag(getProperty(id)); }
-    void clearPropertyDeviceFlag(Property &prop) { prop_flags->clearDeviceFlag(prop.getId()); }
-    void copyPropertyToDevice(property_t id, action_t action) { copyPropertyToDevice(getProperty(id), action); }
-    void copyPropertyToDevice(Property &prop, action_t action);
-
-    void setPropertyHostFlag(property_t id) { setPropertyHostFlag(getProperty(id)); }
-    void setPropertyHostFlag(Property &prop) { prop_flags->setHostFlag(prop.getId()); }
-    void clearPropertyHostFlag(property_t id) { clearPropertyHostFlag(getProperty(id)); }
-    void clearPropertyHostFlag(Property &prop) { prop_flags->clearHostFlag(prop.getId()); }
-    void copyPropertyToHost(property_t id, action_t action) { copyPropertyToHost(getProperty(id), action); }
-    void copyPropertyToHost(Property &prop, action_t action);
-
-    void setContactPropertyDeviceFlag(property_t id) {
-        setContactPropertyDeviceFlag(getContactProperty(id));
+    inline FloatProperty &getAsFloatProperty(Property &prop) {
+        return static_cast<FloatProperty&>(prop);
     }
 
-    void setContactPropertyDeviceFlag(ContactProperty &prop) {
-        contact_prop_flags->setDeviceFlag(prop.getId());
+    inline VectorProperty &getAsVectorProperty(Property &prop) {
+        return static_cast<VectorProperty&>(prop);
     }
 
-    void clearContactPropertyDeviceFlag(property_t id) {
-        clearContactPropertyDeviceFlag(getContactProperty(id));
+    inline MatrixProperty &getAsMatrixProperty(Property &prop) {
+        return static_cast<MatrixProperty&>(prop);
     }
 
-    void clearContactPropertyDeviceFlag(ContactProperty &prop) {
-        contact_prop_flags->clearDeviceFlag(prop.getId());
+    inline QuaternionProperty &getAsQuaternionProperty(Property &prop) {
+        return static_cast<QuaternionProperty&>(prop);
     }
 
-    void copyContactPropertyToDevice(property_t id, action_t action) {
-        copyContactPropertyToDevice(getContactProperty(id), action);
+    inline IntProperty &getIntegerProperty(property_t property) {
+        return static_cast<IntProperty&>(getProperty(property));
     }
 
-    void copyContactPropertyToDevice(ContactProperty &prop, action_t action);
+    inline FloatProperty &getFloatProperty(property_t property) {
+        return static_cast<FloatProperty&>(getProperty(property));
+    }
+
+    inline VectorProperty &getVectorProperty(property_t property) {
+        return static_cast<VectorProperty&>(getProperty(property));
+    }
 
-    void setContactPropertyHostFlag(property_t id) {
-        setContactPropertyHostFlag(getContactProperty(id));
+    inline MatrixProperty &getMatrixProperty(property_t property) {
+        return static_cast<MatrixProperty&>(getProperty(property));
     }
 
-    void setContactPropertyHostFlag(ContactProperty &prop) {
-        contact_prop_flags->setHostFlag(prop.getId());
+    inline QuaternionProperty &getQuaternionProperty(property_t property) {
+        return static_cast<QuaternionProperty&>(getProperty(property));
     }
 
-    void clearContactPropertyHostFlag(property_t id) {
-        clearContactPropertyHostFlag(getContactProperty(id));
+    void copyPropertyToDevice(property_t id, action_t action, size_t size = 0) {
+        copyPropertyToDevice(getProperty(id), action, size);
     }
 
-    void clearContactPropertyHostFlag(ContactProperty &prop) {
-        contact_prop_flags->clearHostFlag(prop.getId());
+    void copyPropertyToDevice(Property &prop, action_t action, size_t size = 0);
+
+    void copyPropertyToHost(property_t id, action_t action, size_t size = 0) {
+        copyPropertyToHost(getProperty(id), action, size);
     }
 
-    void copyContactPropertyToHost(property_t id, action_t action) {
-        copyContactPropertyToHost(getContactProperty(id), action);
+    void copyPropertyToHost(Property &prop, action_t action, size_t size = 0);
+
+    // Contact properties
+    ContactProperty &getContactProperty(property_t id);
+    ContactProperty &getContactPropertyByName(std::string name);
+    void addContactProperty(ContactProperty prop);
+
+    template<typename T_ptr>
+    void addContactProperty(
+        property_t id, std::string name, T_ptr **h_ptr, std::nullptr_t,
+        PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
+
+    template<typename T_ptr>
+    void addContactProperty(
+        property_t id, std::string name, T_ptr **h_ptr, T_ptr **d_ptr,
+        PropertyType type, layout_t layout, size_t sx, size_t sy = 1);
+
+    template<typename T_ptr>
+    void reallocContactProperty(
+        property_t id, T_ptr **h_ptr, std::nullptr_t, size_t sx = 1, size_t sy = 1);
+
+    template<typename T_ptr>
+    void reallocContactProperty(
+        property_t id, T_ptr **h_ptr, T_ptr **d_ptr, size_t sx = 1, size_t sy = 1);
+
+    void copyContactPropertyToDevice(property_t id, action_t action, size_t size = 0) {
+        copyContactPropertyToDevice(getContactProperty(id), action, size);
     }
 
-    void copyContactPropertyToHost(ContactProperty &prop, action_t action);
+    void copyContactPropertyToDevice(ContactProperty &prop, action_t action, size_t size = 0);
+
+    void copyContactPropertyToHost(property_t id, action_t action, size_t size = 0) {
+        copyContactPropertyToHost(getContactProperty(id), action, size);
+    }
+
+    void copyContactPropertyToHost(ContactProperty &prop, action_t action, size_t size = 0);
+
+    // Feature properties
+    FeatureProperty &getFeatureProperty(property_t id);
+    FeatureProperty &getFeaturePropertyByName(std::string name);
+    void addFeatureProperty(FeatureProperty feature_prop);
+
+    template<typename T_ptr>
+    void addFeatureProperty(
+        property_t id, std::string name, T_ptr *h_ptr, std::nullptr_t,
+        PropertyType type, int nkinds, int array_size);
+
+    template<typename T_ptr>
+    void addFeatureProperty(
+        property_t id, std::string name, T_ptr *h_ptr, T_ptr *d_ptr,
+        PropertyType type, int nkinds, int array_size);
 
     void copyFeaturePropertyToDevice(property_t id) {
         copyFeaturePropertyToDevice(getFeatureProperty(id));
@@ -196,6 +222,12 @@ public:
 
     void copyFeaturePropertyToDevice(FeatureProperty &feature_prop);
 
+    // Communication
+    void initDomain(
+        int *argc, char ***argv,
+        real_t xmin, real_t xmax, real_t ymin, real_t ymax, real_t zmin, real_t zmax);
+
+    Regular6DStencil *getDomainPartitioner() { return dom_part; }
     void communicateSizes(int dim, const int *send_sizes, int *recv_sizes);
     void communicateData(
         int dim, int elem_size,
@@ -203,7 +235,12 @@ public:
         real_t *recv_buf, const int *recv_offsets, const int *nrecv);
 
     void fillCommunicationArrays(int neighbor_ranks[], int pbc[], real_t subdom[]);
+
+    // Device functions
     void sync() { device_synchronize(); }
+
+    // Timers
+    Timers<double> *getTimers() { return timers; }
     void printTimers() {
         if(this->getDomainPartitioner()->getRank() == 0) {
             this->getTimers()->print();
@@ -310,7 +347,7 @@ void PairsSimulation::reallocProperty(property_t id, T_ptr **h_ptr, std::nullptr
 			  [id](Property _p) { return _p.getId() == id; });
     PAIRS_ASSERT(p != std::end(properties));
 
-    size_t size = sx * sy * p->getElemSize();
+    size_t size = sx * sy * p->getPrimitiveTypeSize();
     PAIRS_ASSERT(size > 0);
 
     size_t old_size = p->getTotalSize();
@@ -329,7 +366,7 @@ void PairsSimulation::reallocProperty(property_t id, T_ptr **h_ptr, T_ptr **d_pt
 			  [id](Property _p) { return _p.getId() == id; });
     PAIRS_ASSERT(p != std::end(properties));
 
-    size_t size = sx * sy * p->getElemSize();
+    size_t size = sx * sy * p->getPrimitiveTypeSize();
     PAIRS_ASSERT(size > 0);
 
     size_t old_size = p->getTotalSize();
@@ -380,7 +417,7 @@ void PairsSimulation::reallocContactProperty(property_t id, T_ptr **h_ptr, std::
 			   [id](ContactProperty _cp) { return _cp.getId() == id; });
     PAIRS_ASSERT(cp != std::end(contact_properties));
 
-    size_t size = sx * sy * cp->getElemSize();
+    size_t size = sx * sy * cp->getPrimitiveTypeSize();
     PAIRS_ASSERT(size > 0);
 
     size_t old_size = cp->getTotalSize();
@@ -399,7 +436,7 @@ void PairsSimulation::reallocContactProperty(property_t id, T_ptr **h_ptr, T_ptr
 			   [id](ContactProperty _cp) { return _cp.getId() == id; });
     PAIRS_ASSERT(cp != std::end(contact_properties));
 
-    size_t size = sx * sy * cp->getElemSize();
+    size_t size = sx * sy * cp->getPrimitiveTypeSize();
     PAIRS_ASSERT(size > 0);
 
     size_t old_size = cp->getTotalSize();
diff --git a/runtime/property.hpp b/runtime/property.hpp
index 58b55949160cf1647f5b74b3474cdf204f8a84fe..741594d1745edf923f63bf8eae422aa11218f41a 100644
--- a/runtime/property.hpp
+++ b/runtime/property.hpp
@@ -32,10 +32,10 @@ public:
     void *getDevicePointer() { return d_ptr; }
     void setPointers(void *h_ptr_, void *d_ptr_) { h_ptr = h_ptr_, d_ptr = d_ptr_; }
     void setSizes(size_t sx_, size_t sy_) { sx = sx_, sy = sy_; }
-    size_t getTotalSize() { return sx * sy * getElemSize(); };
+    size_t getTotalSize() { return sx * sy * getPrimitiveTypeSize(); };
     PropertyType getType() { return type; }
     layout_t getLayout() { return layout; }
-    size_t getElemSize() {
+    size_t getPrimitiveTypeSize() {
         return  (type == Prop_Integer) ? sizeof(int) :
                 (type == Prop_Real) ? sizeof(real_t) :
                 (type == Prop_Vector) ? sizeof(real_t) :
diff --git a/src/pairs/code_gen/cgen.py b/src/pairs/code_gen/cgen.py
index 80716a8383a26aaebc8ee58e769b8c7adc1be2f9..e5c0300d3232f5edff01a14583f9e6bc5edb32a8 100644
--- a/src/pairs/code_gen/cgen.py
+++ b/src/pairs/code_gen/cgen.py
@@ -478,14 +478,16 @@ class CGen:
             prop_name = ast_node.contact_prop().name()
             action = Actions.c_keyword(ast_node.action())
             ctx_suffix = "Device" if ast_node.context() == Contexts.Device else "Host"
-            self.print(f"pairs->copyContactPropertyTo{ctx_suffix}({prop_id}, {action}); // {prop_name}")
+            size = self.generate_expression(ast_node.contact_prop().copy_size())
+            self.print(f"pairs->copyContactPropertyTo{ctx_suffix}({prop_id}, {action}, {size}); // {prop_name}")
 
         if isinstance(ast_node, CopyProperty):
             prop_id = ast_node.prop().id()
             prop_name = ast_node.prop().name()
             action = Actions.c_keyword(ast_node.action())
             ctx_suffix = "Device" if ast_node.context() == Contexts.Device else "Host"
-            self.print(f"pairs->copyPropertyTo{ctx_suffix}({prop_id}, {action}); // {prop_name}")
+            size = self.generate_expression(ast_node.prop().copy_size())
+            self.print(f"pairs->copyPropertyTo{ctx_suffix}({prop_id}, {action}, {size}); // {prop_name}")
 
         if isinstance(ast_node, CopyVar):
             var_name = ast_node.variable().name()
diff --git a/src/pairs/ir/arrays.py b/src/pairs/ir/arrays.py
index a6f5da5a0016f686213de510f18edccc10d46d46..15c55503e9229def0dd6258d7185b5fbc4eb9e01 100644
--- a/src/pairs/ir/arrays.py
+++ b/src/pairs/ir/arrays.py
@@ -48,8 +48,8 @@ class Array(ASTNode):
         super().__init__(sim)
         self.arr_id = Array.last_array_id
         self.arr_name = a_name
-        self.arr_sizes = [Lit.cvt(sim, a_sizes)] if not isinstance(a_sizes, list) \
-                         else [Lit.cvt(sim, s) for s in a_sizes]
+        self.arr_sizes = [Lit.cvt(sim, a_sizes)] if not isinstance(a_sizes, list) else \
+                         [Lit.cvt(sim, s) for s in a_sizes]
         self.arr_type = a_type
         self.arr_layout = a_layout
         self.arr_sync = a_sync
@@ -93,7 +93,8 @@ class Array(ASTNode):
         self.arr_strides[dim] = stride
 
     def strides(self):
-        return [self.arr_strides[i] if i in self.arr_strides else self.arr_sizes[i] for i in range(self.arr_ndims)]
+        return [self.arr_strides[i] if i in self.arr_strides else self.arr_sizes[i] \
+                for i in range(self.arr_ndims)]
 
     def alloc_size(self):
         return reduce((lambda x, y: x * y), [s for s in self.arr_sizes])
diff --git a/src/pairs/ir/properties.py b/src/pairs/ir/properties.py
index 85024f1d5448ad092f31aa77d87aec77520ed9aa..f841199ea91822434fe29415e799b3811682d978 100644
--- a/src/pairs/ir/properties.py
+++ b/src/pairs/ir/properties.py
@@ -6,6 +6,7 @@ from pairs.ir.layouts import Layouts
 from pairs.ir.lit import Lit
 from pairs.ir.operator_class import OperatorClass
 from pairs.ir.scalars import ScalarOp
+from pairs.ir.sizeof import Sizeof
 from pairs.ir.types import Types
 
 
@@ -85,9 +86,14 @@ class Property(ASTNode):
     def ndims(self):
         return 1 if Types.is_scalar(self.prop_type) else 2
 
+    def copy_size(self):
+        return ScalarOp.inline((self.sim.nlocal + self.sim.nghost) *
+                               Types.number_of_elements(self.sim, self.prop_type) *
+                               Sizeof(self.sim, self.type()))
+
     def sizes(self):
-        return [self.sim.particle_capacity] if Types.is_scalar(self.prop_type) \
-               else [Types.number_of_elements(self.sim, self.prop_type), self.sim.particle_capacity]
+        return [self.sim.particle_capacity] if Types.is_scalar(self.prop_type) else \
+               [Types.number_of_elements(self.sim, self.prop_type), self.sim.particle_capacity]
 
 
 class PropertyAccess(ASTTerm):
@@ -249,10 +255,16 @@ class ContactProperty(ASTNode):
     def ndims(self):
         return 1 if Types.is_scalar(self.contact_prop_type) else 2
 
+    def copy_size(self):
+        return ScalarOp.inline(self.sim.nlocal *
+                               self.sim.neighbor_capacity *
+                               Types.number_of_elements(self.sim, self.prop_type) *
+                               Sizeof(self.sim, self.type()))
+
     def sizes(self):
         neighbor_list_sizes = [self.sim.particle_capacity, self.sim.neighbor_capacity]
-        return neighbor_list_sizes if Types.is_scalar(self.contact_prop_type) \
-               else [Types.number_of_elements(self.sim, self.contact_prop_type)] + neighbor_list_sizes
+        return neighbor_list_sizes if Types.is_scalar(self.contact_prop_type) else \
+               [Types.number_of_elements(self.sim, self.contact_prop_type)] + neighbor_list_sizes
 
 
 class ContactPropertyAccess(ASTTerm):