From: wolf Date: Tue, 23 Mar 2004 19:55:56 +0000 (+0000) Subject: Add parallel vectors. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=9ec4ed79a320aa0278b301fbe467ab63c93e9006;p=dealii-svn.git Add parallel vectors. git-svn-id: https://svn.dealii.org/trunk@8851 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/petsc_parallel_vector.h b/deal.II/lac/include/lac/petsc_parallel_vector.h new file mode 100644 index 0000000000..9faf31ef19 --- /dev/null +++ b/deal.II/lac/include/lac/petsc_parallel_vector.h @@ -0,0 +1,177 @@ +//---------------------------- petsc_vector.h --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- petsc_vector.h --------------------------- +#ifndef __deal2__petsc_vector_h +#define __deal2__petsc_vector_h + +#include +#include +#include + +#include + +#ifdef DEAL_II_USE_PETSC + +#include + + +namespace PETScWrappers +{ +/** + * Implementation of a parallel vector class based on PETSC and using MPI + * communication to synchronise distributed operations. All the functionality + * is actually in the base class, except for the calls to generate a + * parallel vector. This is possible since PETSc only works on an abstract + * vector type and internally distributes to functions that do the actual work + * depending on the actual vector type (much like using virtual + * functions). Only the functions creating a vector of specific type differ, + * and are implemented in this particular class. + * + * @author Wolfgang Bangerth, 2004 + */ + class ParallelVector : public VectorBase + { + public: + /** + * Default constructor. Initialize the + * vector as empty. + */ + ParallelVector (); + + /** + * Constructor. Set dimension to + * @p{n} and initialize all + * elements with zero. + * + * The constructor is made explicit to + * avoid accidents like this: + * @p{v=0;}. Presumably, the user wants + * to set every element of the vector to + * zero, but instead, what happens is + * this call: @p{v=Vector(0);}, + * i.e. the vector is replaced by one of + * length zero. + */ + explicit ParallelVector (const unsigned int n, + const unsigned int local_size, + const MPI_Comm &communicator); + + /** + * Copy-constructor from deal.II + * vectors. Sets the dimension to that + * of the given vector, and copies all + * elements. + */ + template + explicit ParallelVector (const ::Vector &v, + const unsigned int local_size, + const MPI_Comm &communicator); + + /** + * Copy-constructor the + * values from a PETSc wrapper vector + * class. + */ + explicit ParallelVector (const VectorBase &v, + const unsigned int local_size, + const MPI_Comm &communicator); + + /** + * Set all components of the vector to + * the given number @p{s}. Simply pass + * this down to the base class, but we + * still need to declare this function + * to make the example given in the + * discussion about making the + * constructor explicit work. + */ + ParallelVector & operator = (const PetscScalar s); + + /** + * Copy the values of a deal.II vector + * (as opposed to those of the PETSc + * vector wrapper class) into this + * object. + */ + template + ParallelVector & operator = (const ::Vector &v); + + protected: + /** + * Create a vector of length @p{n}. For + * this class, we create a parallel + * vector. @arg n denotes the total + * size of the vector to be + * created. @arg local_size denotes how + * many of these elements shall be + * stored locally. The last argument is + * ignored for sequential vectors. + */ + virtual void create_vector (const unsigned int n, + const unsigned int local_size); + + private: + /** + * Copy of the communicator object to + * be used for this parallel vector. + */ + MPI_Comm communicator; + }; + + + +// ------------------ template and inline functions ------------- + + + template + ParallelVector::ParallelVector (const ::Vector &v, + const unsigned int local_size, + const MPI_Comm &communicator) + : + communicator (communicator) + { + ParallelVector::create_vector (v.size(), local_size); + + VectorBase::operator = (v); + } + + + + inline + ParallelVector & + ParallelVector::operator = (const PetscScalar s) + { + VectorBase::operator = (s); + + return *this; + } + + + template + inline + ParallelVector & + ParallelVector::operator = (const ::Vector &v) + { + VectorBase::operator = (v); + + return *this; + } + +} + + +#endif // DEAL_II_USE_PETSC + +/*---------------------------- petsc_vector.h ---------------------------*/ + +#endif +/*---------------------------- petsc_vector.h ---------------------------*/ diff --git a/deal.II/lac/source/petsc_parallel_vector.cc b/deal.II/lac/source/petsc_parallel_vector.cc new file mode 100644 index 0000000000..7cc8a18451 --- /dev/null +++ b/deal.II/lac/source/petsc_parallel_vector.cc @@ -0,0 +1,79 @@ +//--------------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2004 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//--------------------------------------------------------------------------- + + +#include + +#include + +#ifdef DEAL_II_USE_PETSC + + +namespace PETScWrappers +{ + + + ParallelVector::ParallelVector () + { + // this is an invalid empty vector, so we + // can just as well create a sequential + // one to avoid all the overhead incurred + // by parallelism + const int n = 0; + const int ierr + = VecCreateSeq (PETSC_COMM_SELF, n, &vector); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + } + + + + ParallelVector::ParallelVector (const unsigned int n, + const unsigned int local_size, + const MPI_Comm &communicator) + : + communicator (communicator) + { + ParallelVector::create_vector (n, local_size); + } + + + + ParallelVector::ParallelVector (const VectorBase &v, + const unsigned int local_size, + const MPI_Comm &communicator) + : + communicator (communicator) + { + ParallelVector::create_vector (v.size(), local_size); + + VectorBase::operator = (v); + } + + + void + ParallelVector::create_vector (const unsigned int n, + const unsigned int local_size) + { + Assert (local_size < n, ExcIndexRange (local_size, 0, n)); + + const int ierr + = VecCreateMPI (PETSC_COMM_SELF, local_size, n, &vector); + AssertThrow (ierr == 0, ExcPETScError(ierr)); + } +} + +#else +// On gcc2.95 on Alpha OSF1, the native assembler does not like empty +// files, so provide some dummy code +namespace { void dummy () {} } +#endif // DEAL_II_USE_PETSC