diff --git a/src/DCPSE/DCPSE_op/DCPSE_op.hpp b/src/DCPSE/DCPSE_op/DCPSE_op.hpp
index 2ebf7845a2fed58ecde4353aae8b1a1734b93c08..9420beaf283fae447571975d6f494eba687d18ab 100644
--- a/src/DCPSE/DCPSE_op/DCPSE_op.hpp
+++ b/src/DCPSE/DCPSE_op/DCPSE_op.hpp
@@ -726,6 +726,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -810,6 +830,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -877,6 +917,26 @@ public:
return vector_dist_expression_op<operand_type, dcpse_type, VECT_DCPSE>(arg, *(dcpse_type *) dcpse);
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1419,6 +1479,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1501,6 +1581,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1583,6 +1683,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1665,6 +1785,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1747,6 +1887,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1827,6 +1987,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1882,6 +2062,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1937,6 +2137,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -1992,6 +2212,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -2047,6 +2287,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -2103,6 +2363,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -2159,6 +2439,26 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -2214,6 +2514,96 @@ public:
}
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
+ /*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
+ *
+ *
+ * \param parts particle set
+ */
+ template<typename particles_type>
+ void update(particles_type &particles) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->initializeUpdate(particles);
+
+ }
+};
+
+
+template<template<unsigned int, typename, typename...> class Dcpse_type = Dcpse>
+class Derivative_G_T {
+
+ void *dcpse;
+
+public:
+
+ template<typename particles_type>
+ Derivative_G_T(particles_type &parts, unsigned int ord, typename particles_type::stype rCut,
+ const Point<particles_type::dims, unsigned int> &p,double oversampling_factor = dcpse_oversampling_factor,
+ support_options opt = support_options::RADIUS) {
+ dcpse = new Dcpse_type<particles_type::dims, particles_type>(parts, p, ord, rCut, oversampling_factor, opt);
+ }
+
+ template<typename operand_type>
+ vector_dist_expression_op<operand_type, Dcpse_type<operand_type::vtype::dims, typename operand_type::vtype>, VECT_DCPSE>
+ operator()(operand_type arg) {
+ typedef Dcpse_type<operand_type::vtype::dims, typename operand_type::vtype> dcpse_type;
+
+ return vector_dist_expression_op<operand_type, dcpse_type, VECT_DCPSE>(arg, *(dcpse_type *) dcpse);
+ }
+
+ template<typename particles_type>
+ void checkMomenta(particles_type &particles) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->checkMomenta(particles);
+
+ }
+
+ template<unsigned int prp, typename particles_type>
+ void DrawKernel(particles_type &particles, int k) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->template DrawKernel<prp>(particles, k);
+
+ }
+
+ /*! \brief Method for Saving the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be saved.
+ */
+ template<typename particles_type>
+ void save(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->save(file);
+ }
+ /*! \brief Method for Loading the DCPSE Operator.
+ *
+ * \param parts particle set
+ * \param file name for data to be loaded from.
+ */
+ template<typename particles_type>
+ void load(particles_type &particles, const std::string &file) {
+ auto dcpse_temp = (Dcpse_type<particles_type::dims, particles_type> *) dcpse;
+ dcpse_temp->load(file);
+ }
/*! \brief Method for Updating the DCPSE Operator by recomputing DCPSE Kernels.
*
*
@@ -2249,6 +2639,8 @@ typedef Derivative_yyy_T<Dcpse> Derivative_yyy;
typedef Derivative_xxxx_T<Dcpse> Derivative_xxxx;
typedef Derivative_yyyy_T<Dcpse> Derivative_yyyy;
typedef Derivative_xxyy_T<Dcpse> Derivative_xxyy;
+typedef Derivative_G_T<Dcpse> Derivative_G;
+
#if defined(__NVCC__)
typedef Derivative_x_T<Dcpse_gpu> Derivative_x_gpu;
@@ -2269,6 +2661,7 @@ typedef Derivative_xxx_T<Dcpse_gpu> Derivative_xxx_gpu;
typedef Derivative_xxy_T<Dcpse_gpu> Derivative_xxy_gpu;
typedef Derivative_yyx_T<Dcpse_gpu> Derivative_yyx_gpu;
typedef Derivative_yyy_T<Dcpse_gpu> Derivative_yyy_gpu;
+typedef Derivative_G_T<Dcpse_gpu> Derivative_G_gpu;
#endif
diff --git a/src/DCPSE/DCPSE_op/tests/DCPSE_op_test_base_tests.cpp b/src/DCPSE/DCPSE_op/tests/DCPSE_op_test_base_tests.cpp
index 372206c9d1a80354f00f9bbd57d229dcd97a80e6..47286818d5b7ff349f2f37825338a339e3a5f86b 100644
--- a/src/DCPSE/DCPSE_op/tests/DCPSE_op_test_base_tests.cpp
+++ b/src/DCPSE/DCPSE_op/tests/DCPSE_op_test_base_tests.cpp
@@ -61,7 +61,7 @@ BOOST_AUTO_TEST_CASE(dcpse_op_tests) {
domain.getLastPos()[1] = y;//+gaussian(rng);
// Here fill the function value
domain.template getLastProp<0>() = sin(domain.getLastPos()[0]) + sin(domain.getLastPos()[1]);
- domain.template getLastProp<2>() = cos(domain.getLastPos()[0]) + cos(domain.getLastPos()[1]);
+ domain.template getLastProp<2>() = 2*cos(domain.getLastPos()[0]) + cos(domain.getLastPos()[1]);
++counter;
++it;
}
@@ -69,7 +69,6 @@ BOOST_AUTO_TEST_CASE(dcpse_op_tests) {
domain.map();
domain.ghost_get<0>();
-
Derivative_x Dx(domain, 2, rCut);
Derivative_y Dy(domain, 2, rCut);
Gradient Grad(domain, 2, rCut);
@@ -77,7 +76,7 @@ BOOST_AUTO_TEST_CASE(dcpse_op_tests) {
auto v = getV<1>(domain);
auto P = getV<0>(domain);
- v = Dx(P) + Dy(P);
+ v = 2*Dx(P) + Dy(P);
auto it2 = domain.getDomainIterator();
double worst = 0.0;
@@ -97,6 +96,145 @@ BOOST_AUTO_TEST_CASE(dcpse_op_tests) {
}
+ BOOST_AUTO_TEST_CASE(dcpse_op_save_load) {
+ size_t edgeSemiSize = 40;
+ const size_t sz[2] = {2 * edgeSemiSize, 2 * edgeSemiSize};
+ Box<2, double> box({0, 0}, {2 * M_PI, 2 * M_PI});
+ size_t bc[2] = {NON_PERIODIC, NON_PERIODIC};
+ double spacing[2];
+ spacing[0] = 2 * M_PI / (sz[0] - 1);
+ spacing[1] = 2 * M_PI / (sz[1] - 1);
+ Ghost<2, double> ghost(spacing[0] * 3.9);
+ double rCut = 3.9 * spacing[0];
+ BOOST_TEST_MESSAGE("Init vector_dist...");
+ double sigma2 = spacing[0] * spacing[1] / (2 * 4);
+
+ vector_dist<2, double, aggregate<double, double, double, VectorS<2, double>, VectorS<2, double>>> domain(0, box,
+ bc,
+ ghost);
+
+ //Init_DCPSE(domain)
+ BOOST_TEST_MESSAGE("Init domain...");
+
+ auto it = domain.getGridIterator(sz);
+ size_t pointId = 0;
+ size_t counter = 0;
+ double minNormOne = 999;
+ while (it.isNext()) {
+ domain.add();
+ auto key = it.get();
+ mem_id k0 = key.get(0);
+ double x = k0 * spacing[0];
+ domain.getLastPos()[0] = x;//+ gaussian(rng);
+ mem_id k1 = key.get(1);
+ double y = k1 * spacing[1];
+ domain.getLastPos()[1] = y;//+gaussian(rng);
+ // Here fill the function value
+ domain.template getLastProp<0>() = sin(domain.getLastPos()[0]) + sin(domain.getLastPos()[1]);
+ domain.template getLastProp<2>() = 2*cos(domain.getLastPos()[0]) + cos(domain.getLastPos()[1]);
+ ++counter;
+ ++it;
+ }
+ BOOST_TEST_MESSAGE("Sync domain across processors...");
+
+ domain.map();
+ domain.ghost_get<0>();
+ Derivative_x Dx(domain, 2, rCut);
+ Derivative_y Dy(domain, 2, rCut);
+ auto v = getV<1>(domain);
+ auto v2 = getV<3>(domain);
+ auto P = getV<0>(domain);
+ v2 = 2*Dx(P) + Dy(P);
+ Dx.save(domain,"DX_test");
+ Dy.save(domain,"DY_test");
+ Derivative_x DxLoaded(domain, 2, rCut,1,support_options::LOAD);
+ Derivative_y DyLoaded(domain, 2, rCut,1,support_options::LOAD);
+ DxLoaded.load(domain,"DX_test");
+ DyLoaded.load(domain,"DY_test");
+ v= 2*DxLoaded(P)+DyLoaded(P);
+ auto it2 = domain.getDomainIterator();
+ double worst = 0.0;
+ while (it2.isNext()) {
+ auto p = it2.get();
+
+ if (fabs(domain.getProp<1>(p) - domain.getProp<2>(p)) > worst) {
+ worst = fabs(domain.getProp<1>(p) - domain.getProp<2>(p));
+ }
+
+ ++it2;
+ }
+ domain.deleteGhost();
+ //std::cout<<worst;
+ BOOST_REQUIRE(worst < 0.03);
+ }
+
+ BOOST_AUTO_TEST_CASE(dcpse_op_save_load2) {
+ size_t edgeSemiSize = 40;
+ const size_t sz[2] = {2 * edgeSemiSize, 2 * edgeSemiSize};
+ Box<2, double> box({0, 0}, {2 * M_PI, 2 * M_PI});
+ size_t bc[2] = {NON_PERIODIC, NON_PERIODIC};
+ double spacing[2];
+ spacing[0] = 2 * M_PI / (sz[0] - 1);
+ spacing[1] = 2 * M_PI / (sz[1] - 1);
+ Ghost<2, double> ghost(spacing[0] * 3.9);
+ double rCut = 3.9 * spacing[0];
+ BOOST_TEST_MESSAGE("Init vector_dist...");
+ double sigma2 = spacing[0] * spacing[1] / (2 * 4);
+
+ vector_dist<2, double, aggregate<double, double, double, VectorS<2, double>, VectorS<2, double>>> domain(0, box,
+ bc,
+ ghost);
+
+ //Init_DCPSE(domain)
+ BOOST_TEST_MESSAGE("Init domain...");
+
+ auto it = domain.getGridIterator(sz);
+ size_t pointId = 0;
+ size_t counter = 0;
+ double minNormOne = 999;
+ while (it.isNext()) {
+ domain.add();
+ auto key = it.get();
+ mem_id k0 = key.get(0);
+ double x = k0 * spacing[0];
+ domain.getLastPos()[0] = x;//+ gaussian(rng);
+ mem_id k1 = key.get(1);
+ double y = k1 * spacing[1];
+ domain.getLastPos()[1] = y;//+gaussian(rng);
+ // Here fill the function value
+ domain.template getLastProp<0>() = sin(domain.getLastPos()[0]) + sin(domain.getLastPos()[1]);
+ domain.template getLastProp<2>() = 2*cos(domain.getLastPos()[0]) + cos(domain.getLastPos()[1]);
+ ++counter;
+ ++it;
+ }
+ BOOST_TEST_MESSAGE("Sync domain across processors...");
+
+ domain.map();
+ domain.ghost_get<0>();
+ auto v = getV<1>(domain);
+ auto v2 = getV<3>(domain);
+ auto P = getV<0>(domain);
+ Derivative_x DxLoaded(domain, 2, rCut,1,support_options::LOAD);
+ Derivative_y DyLoaded(domain, 2, rCut,1,support_options::LOAD);
+ DxLoaded.load(domain,"DX_test");
+ DyLoaded.load(domain,"DY_test");
+ v= 2*DxLoaded(P)+DyLoaded(P);
+ auto it2 = domain.getDomainIterator();
+ double worst = 0.0;
+ while (it2.isNext()) {
+ auto p = it2.get();
+
+ if (fabs(domain.getProp<1>(p) - domain.getProp<2>(p)) > worst) {
+ worst = fabs(domain.getProp<1>(p) - domain.getProp<2>(p));
+ }
+
+ ++it2;
+ }
+ domain.deleteGhost();
+ //std::cout<<worst;
+ BOOST_REQUIRE(worst < 0.03);
+ }
+
BOOST_AUTO_TEST_CASE(dcpse_op_tests_fa) {
size_t edgeSemiSize = 40;
const size_t sz[2] = {2 * edgeSemiSize, 2 * edgeSemiSize};
@@ -157,7 +295,7 @@ BOOST_AUTO_TEST_CASE(dcpse_op_tests) {
}
//std::cout<<"Worst:"<<worst<<std::endl;
domain.deleteGhost();
- //domain.write("test");
+ //domain.write_frame("test",0,0.024,BINARY);
BOOST_REQUIRE(worst < 0.03);
}
diff --git a/src/DCPSE/Dcpse.hpp b/src/DCPSE/Dcpse.hpp
index a973d032fa5b954c38c414737ce4156e48302ff5..f2db680be66e3ac45ee789cce62120c0ffbf4bc3 100644
--- a/src/DCPSE/Dcpse.hpp
+++ b/src/DCPSE/Dcpse.hpp
@@ -75,7 +75,7 @@ private:
vector_type & particlesFrom;
vector_type2 & particlesTo;
- double rCut,supportSizeFactor,nSpacing;
+ double rCut,supportSizeFactor=1,nSpacing;
unsigned int convergenceOrder,nCount;
bool isSurfaceDerivative=false;
@@ -334,7 +334,89 @@ public:
++epsItInvPow;
}
}
+ /*! \brief Save the DCPSE computations
+ *
+ */
+ void save(const std::string &file){
+ auto & v_cl=create_vcluster();
+ size_t req = 0;
+
+ Packer<decltype(localSupports),HeapMemory>::packRequest(localSupports,req);
+ Packer<decltype(localEps),HeapMemory>::packRequest(localEps,req);
+ Packer<decltype(localEpsInvPow),HeapMemory>::packRequest(localEpsInvPow,req);
+ Packer<decltype(calcKernels),HeapMemory>::packRequest(calcKernels,req);
+ Packer<decltype(kerOffsets),HeapMemory>::packRequest(kerOffsets,req);
+
+ // allocate the memory
+ HeapMemory pmem;
+ //pmem.allocate(req);
+ ExtPreAlloc<HeapMemory> mem(req,pmem);
+
+ //Packing
+ Pack_stat sts;
+ Packer<decltype(localSupports),HeapMemory>::pack(mem,localSupports,sts);
+ Packer<decltype(localEps),HeapMemory>::pack(mem,localEps,sts);
+ Packer<decltype(localEpsInvPow),HeapMemory>::pack(mem,localEpsInvPow,sts);
+ Packer<decltype(calcKernels),HeapMemory>::pack(mem,calcKernels,sts);
+ Packer<decltype(kerOffsets),HeapMemory>::pack(mem,kerOffsets,sts);
+
+ // Save into a binary file
+ std::ofstream dump (file+"_"+std::to_string(v_cl.rank()), std::ios::out | std::ios::binary);
+ if (dump.is_open() == false)
+ { std::cerr << __FILE__ << ":" << __LINE__ <<" Unable to write since dump is open at rank "<<v_cl.rank()<<std::endl;
+ return;
+ }
+ dump.write ((const char *)pmem.getPointer(), pmem.size());
+ return;
+ }
+ /*! \brief Load the DCPSE computations
+ *
+ *
+ */
+ void load(const std::string & file)
+ {
+ auto & v_cl=create_vcluster();
+ std::ifstream fs (file+"_"+std::to_string(v_cl.rank()), std::ios::in | std::ios::binary | std::ios::ate );
+ if (fs.is_open() == false)
+ {
+ std::cerr << __FILE__ << ":" << __LINE__ << " error, opening file: " << file << std::endl;
+ return;
+ }
+
+ // take the size of the file
+ size_t sz = fs.tellg();
+ fs.close();
+
+ // reopen the file without ios::ate to read
+ std::ifstream input (file+"_"+std::to_string(v_cl.rank()), std::ios::in | std::ios::binary );
+ if (input.is_open() == false)
+ {//some message here maybe
+ return;}
+
+ // Create the HeapMemory and the ExtPreAlloc memory
+ size_t req = 0;
+ req += sz;
+ HeapMemory pmem;
+ ExtPreAlloc<HeapMemory> mem(req,pmem);
+
+ mem.allocate(pmem.size());
+
+ // read
+ input.read((char *)pmem.getPointer(), sz);
+
+ //close the file
+ input.close();
+
+ //Unpacking
+ Unpack_stat ps;
+ Unpacker<decltype(localSupports),HeapMemory>::unpack(mem,localSupports,ps);
+ Unpacker<decltype(localEps),HeapMemory>::unpack(mem,localEps,ps);
+ Unpacker<decltype(localEpsInvPow),HeapMemory>::unpack(mem,localEpsInvPow,ps);
+ Unpacker<decltype(calcKernels),HeapMemory>::unpack(mem,calcKernels,ps);
+ Unpacker<decltype(kerOffsets),HeapMemory>::unpack(mem,kerOffsets,ps);
+ return;
+ }
void checkMomenta(vector_type &particles)
@@ -805,6 +887,12 @@ private:
this->rCut=rCut;
this->supportSizeFactor=supportSizeFactor;
this->convergenceOrder=convergenceOrder;
+ auto & v_cl=create_vcluster();
+ if(this->opt==LOAD){
+ if(v_cl.rank()==0)
+ {std::cout<<"Warning: Creating empty DC-PSE operator! Please use update or load to get kernels."<<std::endl;}
+ return;
+ }
SupportBuilder<vector_type,vector_type2>
supportBuilder(particlesFrom,particlesTo, differentialSignature, rCut, differentialOrder == 0);
unsigned int requiredSupportSize = monomialBasis.size() * supportSizeFactor;
@@ -875,7 +963,6 @@ private:
++Counter;
}
- auto & v_cl=create_vcluster();
v_cl.sum(avgSpacingGlobal);
v_cl.sum(Counter);
v_cl.execute();
diff --git a/src/DCPSE/DcpseInterpolation.hpp b/src/DCPSE/DcpseInterpolation.hpp
index f7095b41eae26c1e39695478d5d0c40819c202b5..50c01ae21b26681aa2257338c9c83231a10ec2bb 100644
--- a/src/DCPSE/DcpseInterpolation.hpp
+++ b/src/DCPSE/DcpseInterpolation.hpp
@@ -4,6 +4,7 @@
#ifndef OPENFPM_PDATA_DCPSEINTERPOLATION_HPP
#define OPENFPM_PDATA_DCPSEINTERPOLATION_HPP
+#include "DCPSE/Dcpse.hpp"
/*! \brief Class for Creating the DCPSE Operator For the function approximation objects and computes DCPSE Kernels.
*
diff --git a/src/DCPSE/Support.hpp b/src/DCPSE/Support.hpp
index dbfe0d43bd6f3bd6ffa0c3117a3e07b122063822..49282ab8202f66893216b0079913e156cc00c446 100644
--- a/src/DCPSE/Support.hpp
+++ b/src/DCPSE/Support.hpp
@@ -59,6 +59,34 @@ public:
return keys;
}
+ static bool pack()
+ {
+ return true;
+ }
+
+ static bool packRequest()
+ {
+ return true;
+ }
+
+ template<int ... prp> inline void packRequest(size_t & req) const
+ {
+ req += sizeof(size_t);
+ keys.packRequest(req);
+ }
+
+ template<int ... prp> inline void pack(ExtPreAlloc<HeapMemory> & mem, Pack_stat & sts) const
+ {
+ Packer<size_t,HeapMemory>::pack(mem,referencePointKey,sts);
+ keys.template pack<prp ...>(mem,sts);
+ }
+
+ template<unsigned int ... prp, typename MemType> inline void unpack(ExtPreAlloc<MemType> & mem, Unpack_stat & ps)
+ {
+ Unpacker<size_t,MemType>::unpack(mem,referencePointKey,ps);
+ keys.template unpack<prp ...>(mem,ps);
+ }
+
};
diff --git a/src/DCPSE/SupportBuilder.hpp b/src/DCPSE/SupportBuilder.hpp
index 3ab8584ea308931651a1cb2798173035a7a64dc5..ee3ec25479aff7addda56ffb6379a508fd6e3a26 100644
--- a/src/DCPSE/SupportBuilder.hpp
+++ b/src/DCPSE/SupportBuilder.hpp
@@ -17,7 +17,8 @@
enum support_options
{
N_PARTICLES,
- RADIUS
+ RADIUS,
+ LOAD
};
diff --git a/src/Operators/Vector/cuda/vector_dist_operators_cuda.cuh b/src/Operators/Vector/cuda/vector_dist_operators_cuda.cuh
index 0060b134ac7edd6c30961064d781f28a1c199e73..1c15185d6e0f6f2103d5f2469fd8c656366c22f6 100644
--- a/src/Operators/Vector/cuda/vector_dist_operators_cuda.cuh
+++ b/src/Operators/Vector/cuda/vector_dist_operators_cuda.cuh
@@ -26,11 +26,11 @@ struct SubsetSelector_impl<true>
{
template<typename particle_type,typename subset_type>
static void check(particle_type &particles,subset_type &particle_subset){
-
- if(particles.getMapCtr()!=particle_subset.getUpdateCtr())
+ //This getMapCtr needs to be created or fixed for cuda!
+ /* if(particles.getMapCtr()!=particle_subset.getUpdateCtr())
{
std::cerr<<__FILE__<<":"<<__LINE__<<" Error: You forgot a subset update after map."<<std::endl;
- }
+ }*/
}
};
#endif
diff --git a/src/Operators/Vector/vector_dist_operators.hpp b/src/Operators/Vector/vector_dist_operators.hpp
index fed9c142e0d8fd924945697dbb0901d452a737d1..0eb79d6af216fa33036970269397106b55288916 100644
--- a/src/Operators/Vector/vector_dist_operators.hpp
+++ b/src/Operators/Vector/vector_dist_operators.hpp
@@ -2039,9 +2039,9 @@ operator+(const vector_dist_expression<prp1,v1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int prp1 , typename v1>
+template<typename T, unsigned int prp1, typename v1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_SUM>
-operator+(const vector_dist_expression<prp1,v1> & va, float d)
+operator+(const vector_dist_expression<prp1,v1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_SUM> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2107,9 +2107,9 @@ operator+(const vector_dist_expression_op<exp1,exp2,op1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<typename exp1 , typename exp2, unsigned int op1>
+template<typename T, typename exp1 , typename exp2, unsigned int op1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_SUM>
-operator+(const vector_dist_expression_op<exp1,exp2,op1> & va, float d)
+operator+(const vector_dist_expression_op<exp1,exp2,op1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_SUM> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2243,9 +2243,10 @@ operator-(const vector_dist_expression<prp1,v1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int prp1, typename v1>
+//template<unsigned int prp1, typename v1>
+template<typename T, unsigned int prp1,typename v1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_SUB>
-operator-(const vector_dist_expression<prp1,v1> & va, float d)
+operator-(const vector_dist_expression<prp1,v1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_SUB> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2277,9 +2278,9 @@ operator-(double d, const vector_dist_expression<prp1,v1> & vb)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int prp1, typename v1>
+template<typename T, unsigned int prp1,typename v1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<prp1,v1>,VECT_SUB>
-operator-(float d, const vector_dist_expression<prp1,v1> & vb)
+operator-(T d, const vector_dist_expression<prp1,v1> & vb)
{
vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<prp1,v1>,VECT_SUB> exp_sum(vector_dist_expression<0,float>(d),vb);
@@ -2311,9 +2312,9 @@ operator*(double d, const vector_dist_expression<p2,v2> & vb)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int p2, typename v2>
+template<typename T, unsigned int p2,typename v2, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<p2,v2>,VECT_MUL>
-operator*(float d, const vector_dist_expression<p2,v2> & vb)
+operator*(T d, const vector_dist_expression<p2,v2> & vb)
{
vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<p2,v2>,VECT_MUL> exp_sum(vector_dist_expression<0,float>(d),vb);
@@ -2345,9 +2346,9 @@ operator*(const vector_dist_expression<p2,v2> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int p2, typename v2>
+template<typename T, unsigned int p2,typename v2, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<p2,v2>,vector_dist_expression<0,float>,VECT_MUL>
-operator*(const vector_dist_expression<p2,v2> & va, float d)
+operator*(const vector_dist_expression<p2,v2> & va, T d)
{
vector_dist_expression_op<vector_dist_expression<p2,v2>,vector_dist_expression<0,float>,VECT_MUL> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2447,9 +2448,9 @@ operator*(const vector_dist_expression_op<exp1,exp2,op1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<typename exp1 , typename exp2, unsigned int op1>
+template<typename T, typename exp1 , typename exp2, unsigned int op1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_MUL>
-operator*(const vector_dist_expression_op<exp1,exp2,op1> & va, float d)
+operator*(const vector_dist_expression_op<exp1,exp2,op1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_MUL> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2481,9 +2482,9 @@ operator*(double d, const vector_dist_expression_op<exp1,exp2,op1> & vb)
* \return an object that encapsulate the expression
*
*/
-template<typename exp1 , typename exp2, unsigned int op1>
+template<typename T, typename exp1 , typename exp2, unsigned int op1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression_op<exp1,exp2,op1>,VECT_MUL>
-operator*(float d, const vector_dist_expression_op<exp1,exp2,op1> & vb)
+operator*(T d, const vector_dist_expression_op<exp1,exp2,op1> & vb)
{
vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression_op<exp1,exp2,op1>,VECT_MUL> exp_sum(vector_dist_expression<0,float>(d),vb);
@@ -2515,9 +2516,9 @@ operator/(const vector_dist_expression_op<exp1,exp2,op1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<typename exp1, typename exp2, unsigned int op1>
+template<typename T, typename exp1 , typename exp2, unsigned int op1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_DIV>
-operator/(const vector_dist_expression_op<exp1,exp2,op1> & va, float d)
+operator/(const vector_dist_expression_op<exp1,exp2,op1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_DIV> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2549,9 +2550,9 @@ operator/(double d, const vector_dist_expression_op<exp1,exp2,op1> & va)
* \return an object that encapsulate the expression
*
*/
-template<typename exp1, typename exp2, unsigned int op1>
+template<typename T, typename exp1 , typename exp2, unsigned int op1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_DIV>
-operator/(float d, const vector_dist_expression_op<exp1,exp2,op1> & va)
+operator/(T d, const vector_dist_expression_op<exp1,exp2,op1> & va)
{
vector_dist_expression_op<vector_dist_expression_op<exp1,exp2,op1>,vector_dist_expression<0,float>,VECT_DIV> exp_sum(vector_dist_expression<0,float>(d),va);
@@ -2583,9 +2584,9 @@ operator/(const vector_dist_expression<prp1,v1> & va, double d)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int prp1, typename v1>
+template<typename T, unsigned int prp1,typename v1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_DIV>
-operator/(const vector_dist_expression<prp1,v1> & va, float d)
+operator/(const vector_dist_expression<prp1,v1> & va, T d)
{
vector_dist_expression_op<vector_dist_expression<prp1,v1>,vector_dist_expression<0,float>,VECT_DIV> exp_sum(va,vector_dist_expression<0,float>(d));
@@ -2617,9 +2618,9 @@ operator/(double d, const vector_dist_expression<prp1,v1> & va)
* \return an object that encapsulate the expression
*
*/
-template<unsigned int prp1, typename v1>
+template<typename T, unsigned int prp1,typename v1, typename sfinae = typename std::enable_if<std::is_same<T,float>::value>::type >
inline vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<prp1,v1>,VECT_DIV>
-operator/(float d, const vector_dist_expression<prp1,v1> & va)
+operator/(T d, const vector_dist_expression<prp1,v1> & va)
{
vector_dist_expression_op<vector_dist_expression<0,float>,vector_dist_expression<prp1,v1>,VECT_DIV> exp_sum(vector_dist_expression<0,float>(d),va);