PETSC introduction

configure.ac

@@ -25,6 +25,7 @@ m4_ifdef([AX_LAPACK],,[m4_include([m4/ax_lapack.m4])])
 m4_ifdef([AX_SUITESPARSE],,[m4_include([m4/ax_suitesparse.m4])])
 m4_ifdef([AX_EIGEN],,[m4_include([m4/ax_eigen.m4])])
 m4_ifdef([AX_LIB_HDF5],,[m4_include([m4/ax_lib_hdf5.m4])]])
+m4_ifdef([AX_PETSC_LIB],,[m4/ax_petsc_lib.m4])
 
 case $host_os in
    *cygwin*)
@@ -107,6 +108,10 @@ IMMDX_LIB_METIS([],[echo "Cannot detect metis, use the --with-metis option if it
 IMMDX_LIB_PARMETIS([],[echo "Cannot detect parmetis, use the --with-parmetis option if it is not installed in the default location"
                     exit 203])
 
+## Check for PETSC
+
+AX_LIB_PETSC()
+
 #########
 
 ## Check for HDF5

src/Makefile.am

 
-LINKLIBS =  $(SUITESPARSE_LIBS) $(LAPACK_LIBS) $(BLAS_LIBS)  $(METIS_LIB) $(PARMETIS_LIB) $(DEFAULT_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_PROGRAM_OPTIONS_LIB) $(BOOST_IOSTREAMS_LIB) $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(LIBQUADMATH)
+LINKLIBS =  $(SUITESPARSE_LIBS) $(LAPACK_LIBS) $(BLAS_LIBS)  $(METIS_LIB) $(PARMETIS_LIB) $(DEFAULT_LIB) $(PTHREAD_LIBS) $(OPT_LIBS) $(BOOST_LDFLAGS) $(BOOST_PROGRAM_OPTIONS_LIB) $(BOOST_IOSTREAMS_LIB) $(HDF5_LDFLAGS)  $(HDF5_LIBS) $(LIBQUADMATH) $(PETSC_LIB)
 
 noinst_PROGRAMS = numerics
 numerics_SOURCES = main.cpp ../../openfpm_vcluster/src/VCluster.cpp ../../openfpm_devices/src/memory/HeapMemory.cpp ../../openfpm_devices/src/memory/PtrMemory.cpp ../../openfpm_devices/src/Memleak_check.cpp
-numerics_CXXFLAGS = $(AM_CXXFLAGS) $(HDF5_CPPFLAGS) $(INCLUDES_PATH) $(BOOST_CPPFLAGS) $(SUITESPARSE_INCLUDE) $(METIS_INCLUDE) $(PARMETIS_INCLUDE)  $(EIGEN_INCLUDE) -Wno-deprecated-declarations -Wno-unused-local-typedefs
+numerics_CXXFLAGS = $(AM_CXXFLAGS) $(HDF5_CPPFLAGS) $(INCLUDES_PATH) $(BOOST_CPPFLAGS) $(SUITESPARSE_INCLUDE) $(METIS_INCLUDE) $(PARMETIS_INCLUDE)  $(EIGEN_INCLUDE) $(PETSC_INCLUDE) -Wno-deprecated-declarations -Wno-unused-local-typedefs
 numerics_CFLAGS = $(CUDA_CFLAGS)
 numerics_LDADD = $(LINKLIBS) -lparmetis -lmetis
 nobase_include_HEADERS = PSE/Kernels.hpp PSE/Kernels_test_util.hpp

src/Matrix/SparseMatrix.hpp

@@ -50,7 +50,7 @@ template<typename T, unsigned int impl> struct triplet
  *
  * \tparam T Type of the sparse Matrix store on each row,colums
  * \tparam id_t type of id
- * \tparam impl implementation
+ * \tparam Mi implementation
  *
  */
 template<typename T,typename id_t ,typename Mi = Eigen::SparseMatrix<T,0,id_t>>

src/Matrix/SparseMatrix_Eigen.hpp

@@ -139,18 +139,6 @@ public:
 	{
 	}
 
-	/*! \brief Create a Matrix from a set of triplet
-	 *
-	 * \param N1 number of row
-	 * \param N2 number of colums
-	 * \param trpl triplet set
-	 *
-	 */
-	SparseMatrix(size_t N1, size_t N2, const openfpm::vector<Eigen::Triplet<T,id_t>> & trpl)
-	:mat(N1,N2)
-	{
-		this->trpl = trpl;
-	}
 
 	/*! \brief Create an empty Matrix
 	 *

src/Vector/Vector.hpp

@@ -34,6 +34,6 @@ class Vector
 };
 
 #include "Vector_eigen.hpp"
-
+#include "Vector_petsc.hpp"
 
 #endif /* OPENFPM_NUMERICS_SRC_VECTOR_VECTOR_HPP_ */

src/Vector/Vector_eigen.hpp

@@ -89,7 +89,7 @@ class Vector<T,Eigen::Matrix<T, Eigen::Dynamic, 1>>
 	mutable openfpm::vector<size_t> sz;
 
 
-	/*! \brief Here we collect the full matrix on master
+	/*! \brief Here we collect the full vector on master
 	 *
 	 */
 	void collect() const

src/Vector/Vector_unit_tests.hpp

@@ -130,6 +130,96 @@ BOOST_AUTO_TEST_CASE(vector_eigen_parallel)
 }
 
 
+BOOST_AUTO_TEST_CASE(vector_petsc_parallel)
+{
+	Vcluster & vcl = create_vcluster();
+
+	if (vcl.getProcessingUnits() != 3)
+		return;
+
+	// 3 Processors 9x9 Matrix to invert
+
+	Vector<double,Vec> v(9);
+
+	if (vcl.getProcessUnitID() == 0)
+	{
+		// v row1
+		v.insert(0,0);
+		v.insert(1,1);
+		v.insert(2,2);
+	}
+	else if (vcl.getProcessUnitID() == 1)
+	{
+		v.insert(3,3);
+		v.insert(4,4);
+		v.insert(5,5);
+	}
+	else if (vcl.getProcessUnitID() == 2)
+	{
+		v.insert(6,6);
+		v.insert(7,7);
+		v.insert(8,8);
+	}
+
+	Vector<double,Vec> v3;
+	v3 = v;
+
+	if (vcl.getProcessUnitID() == 0)
+	{
+		// v row1
+		v(0) = 8;
+		v(1) = 7;
+		v(2) = 6;
+	}
+	else if (vcl.getProcessUnitID() == 1)
+	{
+		v(3) = 5;
+		v(4) = 4;
+		v(5) = 3;
+	}
+	else if (vcl.getProcessUnitID() == 2)
+	{
+		v(6) = 2;
+		v(7) = 1;
+		v(8) = 0;
+	}
+
+	// Master has the full vector
+
+	if (vcl.getProcessUnitID() == 0)
+	{
+		BOOST_REQUIRE_EQUAL(v(0),8);
+		BOOST_REQUIRE_EQUAL(v(1),7);
+		BOOST_REQUIRE_EQUAL(v(2),6);
+
+		BOOST_REQUIRE_EQUAL(v3(0),0);
+		BOOST_REQUIRE_EQUAL(v3(1),1);
+		BOOST_REQUIRE_EQUAL(v3(2),2);
+	}
+	else if (vcl.getProcessUnitID() == 1)
+	{
+		BOOST_REQUIRE_EQUAL(v(3),5);
+		BOOST_REQUIRE_EQUAL(v(4),4);
+		BOOST_REQUIRE_EQUAL(v(5),3);
+
+		BOOST_REQUIRE_EQUAL(v3(0),3);
+		BOOST_REQUIRE_EQUAL(v3(1),4);
+		BOOST_REQUIRE_EQUAL(v3(2),5);
+	}
+	else if (vcl.getProcessUnitID() == 2)
+	{
+		BOOST_REQUIRE_EQUAL(v(6),2);
+		BOOST_REQUIRE_EQUAL(v(7),1);
+		BOOST_REQUIRE_EQUAL(v(8),0);
+
+		BOOST_REQUIRE_EQUAL(v3(0),6);
+		BOOST_REQUIRE_EQUAL(v3(1),7);
+		BOOST_REQUIRE_EQUAL(v3(2),9);
+	}
+
+	Vec & v2 = v.getVec();
+}
+
 BOOST_AUTO_TEST_SUITE_END()