From 9a682f9d233d2f4634f50b38b8db02c32414fdf5 Mon Sep 17 00:00:00 2001 From: Thomas Witkowski Date: Tue, 21 Aug 2012 08:53:21 +0000 Subject: [PATCH] Modified some files to make use of const FiniteElemSpace pointers, where this is possible. --- AMDiS/src/BoundaryCondition.h | 44 ++++------ AMDiS/src/ComponentTraverseInfo.h | 10 +-- AMDiS/src/DirichletBC.cc | 4 +- AMDiS/src/DirichletBC.h | 4 +- AMDiS/src/MatrixVector.cc | 4 +- AMDiS/src/MatrixVector.h | 2 +- AMDiS/src/PeriodicBC.cc | 2 +- AMDiS/src/PeriodicBC.h | 2 +- AMDiS/src/ProblemStat.cc | 6 +- AMDiS/src/ProblemStat.h | 16 ++-- AMDiS/src/RobinBC.cc | 12 +-- AMDiS/src/RobinBC.h | 20 ++--- AMDiS/src/SystemVector.cc | 14 ++- AMDiS/src/SystemVector.h | 85 ++++++++----------- .../src/parallel/ParallelCoarseSpaceMatVec.cc | 2 +- AMDiS/src/parallel/PetscSolver.cc | 36 -------- AMDiS/src/parallel/PetscSolver.h | 6 -- AMDiS/src/parallel/PetscSolverFeti.cc | 4 +- AMDiS/src/parallel/PetscSolverGlobalMatrix.cc | 19 +---- AMDiS/src/parallel/PetscSolverGlobalMatrix.h | 4 - 20 files changed, 110 insertions(+), 186 deletions(-) diff --git a/AMDiS/src/BoundaryCondition.h b/AMDiS/src/BoundaryCondition.h index e74d980e..d177ef9b 100644 --- a/AMDiS/src/BoundaryCondition.h +++ b/AMDiS/src/BoundaryCondition.h @@ -80,22 +80,18 @@ namespace AMDiS { /// Destructor. virtual ~BoundaryCondition() {} - /** \brief - * Adds the local boundary condition for elInfo to object. - * The dofIndices and localBound as well as nBasFcts are determined by - * the calling BoundaryManager. - */ + /// Adds the local boundary condition for elInfo to object. + /// The dofIndices and localBound as well as nBasFcts are determined by + // the calling BoundaryManager. virtual void fillBoundaryCondition(DOFMatrix *matrix, ElInfo *elInfo, const DegreeOfFreedom *dofIndices, const BoundaryType *localBound, int nBasFcts) {} - /** \brief - * Adds the local boundary condition for elInfo to vector. - * The dofIndices and localBound as well as nBasFcts are determined by - * the calling BoundaryManager. - */ + /// Adds the local boundary condition for elInfo to vector. + /// The dofIndices and localBound as well as nBasFcts are determined by + /// the calling BoundaryManager. virtual void fillBoundaryCondition(DOFVectorBase *vector, ElInfo *elInfo, const DegreeOfFreedom *dofIndices, @@ -110,10 +106,8 @@ namespace AMDiS { return 0.0; } - /** \brief - * Returns whether the condition must be treated as Dirichlet condition - * while assemblage. - */ + /// Returns whether the condition must be treated as Dirichlet condition + /// while assemblage. virtual bool isDirichlet() { return false; @@ -125,25 +119,21 @@ namespace AMDiS { return false; } - /** \brief - * In some situations it may be required to set Dirichlet boundary conditions, - * but not to apply them to the matrix. This is for example the case, if the - * boundary condition is set to a couple matrix. Then, the boundary conditions - * must be applied to the couple matrix, but they are set to all matrices in this - * row (to ensure that there are no other element entries in the Dirichlet boundary - * condition rows). - */ + /// In some situations it may be required to set Dirichlet boundary + /// conditions, but not to apply them to the matrix. This is for example the + /// case, if the boundary condition is set to a couple matrix. Then, the + /// boundary conditions must be applied to the couple matrix, but they are + /// set to all matrices in this row (to ensure that there are no other + /// element entries in the Dirichlet boundary condition rows). virtual bool applyBoundaryCondition() { return true; } protected: - /** \brief - * Speciefies for which parts of the boundary the condition holds. - * This id corresponds to the boundary numbers spcified in the - * macro file. - */ + /// Speciefies for which parts of the boundary the condition holds. + /// This id corresponds to the boundary numbers spcified in the + /// macro file. BoundaryType boundaryType; /// FiniteElemSpace for this BoundaryCondition. diff --git a/AMDiS/src/ComponentTraverseInfo.h b/AMDiS/src/ComponentTraverseInfo.h index f8f38eba..5e55fccf 100644 --- a/AMDiS/src/ComponentTraverseInfo.h +++ b/AMDiS/src/ComponentTraverseInfo.h @@ -38,7 +38,7 @@ namespace AMDiS { status(0) {} - void setFeSpace(FiniteElemSpace *row, FiniteElemSpace *col = NULL) + void setFeSpace(const FiniteElemSpace *row, const FiniteElemSpace *col = NULL) { rowFeSpace = row; colFeSpace = col; @@ -66,12 +66,12 @@ namespace AMDiS { return auxFeSpaces.size(); } - FiniteElemSpace *getRowFeSpace() + const FiniteElemSpace *getRowFeSpace() { return rowFeSpace; } - FiniteElemSpace *getColFeSpace() + const FiniteElemSpace *getColFeSpace() { return colFeSpace; } @@ -94,9 +94,9 @@ namespace AMDiS { } protected: - FiniteElemSpace *rowFeSpace; + const FiniteElemSpace *rowFeSpace; - FiniteElemSpace *colFeSpace; + const FiniteElemSpace *colFeSpace; std::set auxFeSpaces; diff --git a/AMDiS/src/DirichletBC.cc b/AMDiS/src/DirichletBC.cc index f53ca3d1..face0ddd 100644 --- a/AMDiS/src/DirichletBC.cc +++ b/AMDiS/src/DirichletBC.cc @@ -20,8 +20,8 @@ namespace AMDiS { DirichletBC::DirichletBC(BoundaryType type, AbstractFunction > *fct, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace, + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace, bool apply) : BoundaryCondition(type, rowFeSpace, colFeSpace), f(fct), diff --git a/AMDiS/src/DirichletBC.h b/AMDiS/src/DirichletBC.h index 312bde44..2babebca 100644 --- a/AMDiS/src/DirichletBC.h +++ b/AMDiS/src/DirichletBC.h @@ -45,8 +45,8 @@ namespace AMDiS { /// Constructor. DirichletBC(BoundaryType type, AbstractFunction > *fct, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL, + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL, bool apply = true); /// Constructor. diff --git a/AMDiS/src/MatrixVector.cc b/AMDiS/src/MatrixVector.cc index 361fc0cd..7909901e 100644 --- a/AMDiS/src/MatrixVector.cc +++ b/AMDiS/src/MatrixVector.cc @@ -16,9 +16,9 @@ namespace AMDiS { - vector getFeSpaces(Matrix &mat) + vector getComponentFeSpaces(Matrix &mat) { - FUNCNAME("getFeSpaces()"); + FUNCNAME("getComponentFeSpaces()"); int nComponents = mat.getNumRows(); vector result(nComponents); diff --git a/AMDiS/src/MatrixVector.h b/AMDiS/src/MatrixVector.h index 9edba99b..b8d06ce8 100644 --- a/AMDiS/src/MatrixVector.h +++ b/AMDiS/src/MatrixVector.h @@ -521,7 +521,7 @@ namespace AMDiS { /// Returns a vector with the FE spaces of each row in the matrix. Thus, the /// resulting vector may contain the same FE space multiple times. - vector getFeSpaces(Matrix &mat); + vector getComponentFeSpaces(Matrix &mat); } #endif // AMDIS_MATRIXVECTOR_H diff --git a/AMDiS/src/PeriodicBC.cc b/AMDiS/src/PeriodicBC.cc index 26ef5d8e..1fb18b6a 100644 --- a/AMDiS/src/PeriodicBC.cc +++ b/AMDiS/src/PeriodicBC.cc @@ -97,7 +97,7 @@ namespace AMDiS { } - PeriodicBC::PeriodicBC(BoundaryType type, FiniteElemSpace *rowSpace) + PeriodicBC::PeriodicBC(BoundaryType type, const FiniteElemSpace *rowSpace) : BoundaryCondition(type, rowSpace, NULL), masterMatrix(NULL) { diff --git a/AMDiS/src/PeriodicBC.h b/AMDiS/src/PeriodicBC.h index 5f83a690..6ab6785d 100644 --- a/AMDiS/src/PeriodicBC.h +++ b/AMDiS/src/PeriodicBC.h @@ -88,7 +88,7 @@ namespace AMDiS { { public: /// Constructor. - PeriodicBC(BoundaryType type, FiniteElemSpace *rowFeSpace); + PeriodicBC(BoundaryType type, const FiniteElemSpace *rowFeSpace); ~PeriodicBC(); diff --git a/AMDiS/src/ProblemStat.cc b/AMDiS/src/ProblemStat.cc index 8e178769..27fddd4d 100644 --- a/AMDiS/src/ProblemStat.cc +++ b/AMDiS/src/ProblemStat.cc @@ -1586,9 +1586,7 @@ namespace AMDiS { FUNCNAME("ProblemStat::addPeriodicBC()"); boundaryConditionSet = true; - - FiniteElemSpace *feSpace = componentSpaces[row]; - PeriodicBC *periodic = new PeriodicBC(type, feSpace); + PeriodicBC *periodic = new PeriodicBC(type, componentSpaces[row]); if (systemMatrix && (*systemMatrix)[row][col]) (*systemMatrix)[row][col]->getBoundaryManager()->addBoundaryCondition(periodic); @@ -1628,7 +1626,7 @@ namespace AMDiS { } - void ProblemStatSeq::assembleOnOneMesh(FiniteElemSpace *feSpace, + void ProblemStatSeq::assembleOnOneMesh(const FiniteElemSpace *feSpace, Flag assembleFlag, DOFMatrix *matrix, DOFVector *vector) diff --git a/AMDiS/src/ProblemStat.h b/AMDiS/src/ProblemStat.h index faa4f8cb..bcd9cae8 100644 --- a/AMDiS/src/ProblemStat.h +++ b/AMDiS/src/ProblemStat.h @@ -256,7 +256,7 @@ namespace AMDiS { /// This function assembles a DOFMatrix and a DOFVector for the case, /// the meshes from row and col FE-space are equal. - void assembleOnOneMesh(FiniteElemSpace *feSpace, + void assembleOnOneMesh(const FiniteElemSpace *feSpace, Flag assembleFlag, DOFMatrix *matrix, DOFVector *vector); @@ -321,7 +321,7 @@ namespace AMDiS { } /// Returns \ref feSpace_. - inline FiniteElemSpace* getFeSpace(int comp = 0) + inline const FiniteElemSpace* getFeSpace(int comp = 0) { FUNCNAME("ProblemStatSeq::getFeSpace()"); TEST_EXIT(comp < static_cast(componentSpaces.size()) && comp >= 0) @@ -330,13 +330,13 @@ namespace AMDiS { } /// Returns \ref feSpaces. - inline vector getFeSpaces() + inline vector getFeSpaces() { return feSpaces; } /// Returns \ref componentSpaces; - inline vector getComponentFeSpaces() + inline vector getComponentFeSpaces() { return componentSpaces; } @@ -428,12 +428,12 @@ namespace AMDiS { } /// Sets the FE space for the given component. - inline void setFeSpace(FiniteElemSpace *feSpace, int comp = 0) + inline void setFeSpace(const FiniteElemSpace *feSpace, int comp = 0) { feSpaces[comp] = feSpace; } - void setComponentSpace(int comp, FiniteElemSpace *feSpace) + void setComponentSpace(int comp, const FiniteElemSpace *feSpace) { if (static_cast(componentSpaces.size()) < nComponents) componentSpaces.resize(nComponents); @@ -565,13 +565,13 @@ namespace AMDiS { int nMeshes; /// FE spaces of this problem. - vector feSpaces; + vector feSpaces; /// Meshes of this problem. vector meshes; /// Pointer to the fe spaces for the different problem components - vector componentSpaces; + vector componentSpaces; /// Pointer to the meshes for the different problem components vector componentMeshes; diff --git a/AMDiS/src/RobinBC.cc b/AMDiS/src/RobinBC.cc index fa0f6d94..b16423be 100644 --- a/AMDiS/src/RobinBC.cc +++ b/AMDiS/src/RobinBC.cc @@ -24,8 +24,8 @@ namespace AMDiS { RobinBC::RobinBC(BoundaryType type, AbstractFunction > *j, AbstractFunction > *alpha, - FiniteElemSpace *rowFeSpace_, - FiniteElemSpace *colFeSpace_) + const FiniteElemSpace *rowFeSpace_, + const FiniteElemSpace *colFeSpace_) : BoundaryCondition(type, rowFeSpace_, colFeSpace_), neumannOperators(NULL), robinOperators(NULL) @@ -75,8 +75,8 @@ namespace AMDiS { RobinBC::RobinBC(BoundaryType type, DOFVectorBase *j, DOFVectorBase *alpha, - FiniteElemSpace *rowFeSpace_, - FiniteElemSpace *colFeSpace_) + const FiniteElemSpace *rowFeSpace_, + const FiniteElemSpace *colFeSpace_) : BoundaryCondition(type, rowFeSpace_, colFeSpace_), neumannOperators(NULL), robinOperators(NULL) @@ -125,8 +125,8 @@ namespace AMDiS { RobinBC::RobinBC(BoundaryType type, Operator* jOp, Operator* alphaOp, - FiniteElemSpace *rowFeSpace_, - FiniteElemSpace *colFeSpace_) + const FiniteElemSpace *rowFeSpace_, + const FiniteElemSpace *colFeSpace_) : BoundaryCondition(type, rowFeSpace_, colFeSpace_), neumannOperators(NULL), robinOperators(NULL) diff --git a/AMDiS/src/RobinBC.h b/AMDiS/src/RobinBC.h index 4bbbac48..d2a81217 100644 --- a/AMDiS/src/RobinBC.h +++ b/AMDiS/src/RobinBC.h @@ -45,21 +45,21 @@ namespace AMDiS { RobinBC(BoundaryType type, AbstractFunction > *j, AbstractFunction > *alpha, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL); + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL); /// Constructor. \f$ j \f$ and \f$ alpha \f$ are given as DOFVectors. RobinBC(BoundaryType type, DOFVectorBase *j, DOFVectorBase *alpha, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL); + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL); /// Constructor. \f$ j \f$ and \f$ alpha \f$ are given as Operator objects. RobinBC(BoundaryType type, Operator* jOp, Operator* alphaOp, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL); + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL); /// Implements BoundaryCondition::fillBoundaryCondition(); virtual void fillBoundaryCondition(DOFMatrix* matrix, @@ -96,15 +96,15 @@ namespace AMDiS { public: NeumannBC(BoundaryType type, AbstractFunction > *j, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL) + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL) : RobinBC(type, j, NULL, rowFeSpace, colFeSpace) {} NeumannBC(BoundaryType type, DOFVectorBase *j, - FiniteElemSpace *rowFeSpace, - FiniteElemSpace *colFeSpace = NULL) + const FiniteElemSpace *rowFeSpace, + const FiniteElemSpace *colFeSpace = NULL) : RobinBC(type, j, NULL, rowFeSpace, colFeSpace) {} }; diff --git a/AMDiS/src/SystemVector.cc b/AMDiS/src/SystemVector.cc index 9b644920..8399b037 100644 --- a/AMDiS/src/SystemVector.cc +++ b/AMDiS/src/SystemVector.cc @@ -17,7 +17,7 @@ namespace AMDiS { void SystemVector::serialize(std::ostream &out) { - int size = vectors.getSize(); + int size = static_cast(vectors.size()); SerUtil::serialize(out, size); for (int i = 0; i < size; i++) vectors[i]->serialize(out); @@ -26,7 +26,7 @@ namespace AMDiS { void SystemVector::deserialize(std::istream &in) { - int size, oldSize = vectors.getSize(); + int size, oldSize = static_cast(vectors.size()); SerUtil::deserialize(in, size); vectors.resize(size); for (int i = oldSize; i < size; i++) @@ -35,4 +35,14 @@ namespace AMDiS { vectors[i]->deserialize(in); } + + vector SystemVector::getComponentFeSpaces() const + { + vector result(vectors.size()); + + for (unsigned int i = 0; i < vectors.size(); i++) + result[i] = vectors[i]->getFeSpace(); + + return result; + } } diff --git a/AMDiS/src/SystemVector.h b/AMDiS/src/SystemVector.h index 47316777..79aa0a4f 100644 --- a/AMDiS/src/SystemVector.h +++ b/AMDiS/src/SystemVector.h @@ -36,57 +36,50 @@ namespace AMDiS { public: /// Constructor. SystemVector(std::string name_, - std::vector feSpace_, + std::vector feSpace_, int size) : name(name_), feSpace(feSpace_), - vectors(size) - - { - vectors.set(NULL); - } + vectors(size, NULL) + {} /// Copy Constructor. SystemVector(const SystemVector& rhs) : name(rhs.getName()), feSpace(rhs.getFeSpaces()), - vectors(rhs.getNumVectors()) - + vectors(rhs.getNumVectors(), NULL) { - for (int i = 0; i < vectors.getSize(); i++) + for (unsigned int i = 0; i < vectors.size(); i++) vectors[i] = new DOFVector(*rhs.getDOFVector(i)); } ~SystemVector() { - for (int i = 0; i < vectors.getSize(); i++) + for (unsigned int i = 0; i < vectors.size(); i++) delete vectors[i]; } - void createNewDOFVectors(std::string name) - { - for (int i = 0; i < vectors.getSize(); i++) - vectors[i] = new DOFVector(feSpace[i], "tmp"); - } - /// Sets \ref vectors[index] = vec. inline void setDOFVector(int index, DOFVector *vec) { - TEST_EXIT_DBG(index < vectors.getSize())("invalid index\n"); + TEST_EXIT_DBG(index < static_cast(vectors.size())) + ("Invalid index %d!\n", index); vectors[index] = vec; } /// Returns \ref vectors[index]. inline DOFVector *getDOFVector(int index) { - TEST_EXIT_DBG(index < vectors.getSize())("invalid index\n"); + TEST_EXIT_DBG(index < static_cast(vectors.size())) + ("Invalid index %d!\n", index); return vectors[index]; } /// Returns \ref vectors[index]. inline const DOFVector *getDOFVector(int index) const { - TEST_EXIT_DBG(index < vectors.getSize())("invalid index\n"); + TEST_EXIT_DBG(index < static_cast(vectors.size())) + ("Invalid index %d!\n", index); return vectors[index]; } @@ -99,8 +92,7 @@ namespace AMDiS { inline int getUsedSize() const { int totalSize = 0; - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + for (unsigned int i = 0; i < vectors.size(); i++) totalSize += vectors[i]->getUsedSize(); return totalSize; } @@ -108,32 +100,32 @@ namespace AMDiS { /// Returns number of contained vectors. inline int getNumVectors() const { - return vectors.getSize(); + return static_cast(vectors.size()); } inline int getSize() const { - return vectors.getSize(); + return static_cast(vectors.size()); } /// Returns the fe space for a given component. - inline FiniteElemSpace *getFeSpace(int i) const + inline const FiniteElemSpace *getFeSpace(int i) const { return feSpace[i]; } /// Returns the fe spaces for all components. - inline std::vector getFeSpaces() const + inline std::vector getFeSpaces() const { return feSpace; } - /** \brief - * Here the system vector is interpreted as one large vector. The given - * is used as a global index which indicates a local vector number and - * a local index on this vector. The operator returns this local vector - * at the local index. - */ + vector getComponentFeSpaces() const; + + /// Here the system vector is interpreted as one large vector. The given + /// is used as a global index which indicates a local vector number and + /// a local index on this vector. The operator returns this local vector + /// at the local index. inline double& operator[](int index) { int localIndex = index; @@ -162,22 +154,20 @@ namespace AMDiS { /// Sets all entries in all vectors to value. inline void set(double value) { - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + for (unsigned int i = 0; i < vectors.size(); i++) vectors[i]->set(value); } inline void setCoarsenOperation(CoarsenOperation op) { - for (int i = 0; i < static_cast(vectors.getSize()); i++) + for (int i = 0; i < static_cast(vectors.size()); i++) vectors[i]->setCoarsenOperation(op); } /// Sets all entries in all vectors to value. inline SystemVector& operator=(double value) { - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + for (int i = 0; i < static_cast(vectors.size()); i++) (*(vectors[i])) = value; return *this; } @@ -185,10 +175,11 @@ namespace AMDiS { /// Assignement operator. inline SystemVector& operator=(const SystemVector& rhs) { - TEST_EXIT_DBG(rhs.vectors.getSize() == vectors.getSize())("invalied sizes\n"); - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + TEST_EXIT_DBG(rhs.vectors.size() == vectors.size())("Invalied sizes!\n"); + + for (int i = 0; i < static_cast(vectors.size()); i++) (*(vectors[i])) = (*(rhs.getDOFVector(i))); + return *this; } @@ -198,16 +189,15 @@ namespace AMDiS { void copy(const SystemVector& rhs) { - int size = vectors.getSize(); - TEST_EXIT_DBG(size == rhs.getNumVectors())("invalid sizes\n"); + int size = static_cast(vectors.size()); + TEST_EXIT_DBG(size == rhs.getNumVectors())("Invalid sizes!\n"); for (int i = 0; i < size; i++) vectors[i]->copy(*(const_cast(rhs).getDOFVector(i))); } void interpol(std::vector >*> *f) { - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + for (unsigned int i = 0; i < vectors.size(); i++) vectors[i]->interpol((*f)[i]); } @@ -219,15 +209,14 @@ namespace AMDiS { void print() { - int size = vectors.getSize(); - for (int i = 0; i < size; i++) + for (unsigned int i = 0; i < vectors.size(); i++) vectors[i]->print(); } int calcMemoryUsage() { int result = 0; - for (int i = 0; i < static_cast(vectors.getSize()); i++) + for (unsigned int i = 0; i < vectors.size(); i++) result += vectors[i]->calcMemoryUsage(); result += sizeof(SystemVector); @@ -241,10 +230,10 @@ namespace AMDiS { std::string name; /// Finite element space. - std::vector feSpace; + std::vector feSpace; /// Local dof vectors. - Vector*> vectors; + std::vector*> vectors; }; diff --git a/AMDiS/src/parallel/ParallelCoarseSpaceMatVec.cc b/AMDiS/src/parallel/ParallelCoarseSpaceMatVec.cc index 0abb7f78..f65ada0b 100644 --- a/AMDiS/src/parallel/ParallelCoarseSpaceMatVec.cc +++ b/AMDiS/src/parallel/ParallelCoarseSpaceMatVec.cc @@ -261,7 +261,7 @@ namespace AMDiS { mpiCommGlobal.Allreduce(&sendValue, &recvAllValues, 1, MPI_INT, MPI_SUM); if (recvAllValues != 0 || alwaysCreateNnzStructure) { - vector feSpaces = getFeSpaces(seqMat); + vector feSpaces = getComponentFeSpaces(seqMat); interiorMap->setComputeMatIndex(!localMatrix); interiorMap->update(feSpaces); diff --git a/AMDiS/src/parallel/PetscSolver.cc b/AMDiS/src/parallel/PetscSolver.cc index 9b5a3848..044f7dfb 100644 --- a/AMDiS/src/parallel/PetscSolver.cc +++ b/AMDiS/src/parallel/PetscSolver.cc @@ -110,40 +110,4 @@ namespace AMDiS { VecScatterDestroy(&scatter); } - - vector PetscSolver::getFeSpaces(SystemVector *vec) - { - FUNCNAME("PetscSolver::getFeSpaces()"); - - int nComponents = vec->getSize(); - vector result(nComponents); - - for (int i = 0; i < nComponents; i++) - result[i] = vec->getFeSpace(i); - -#if (DEBUG != 0) - // === In debug mode, we test if all FE spaces of the matrix are also === - // === considered by the mesh distributor. === - checkFeSpaces(result); -#endif - - return result; - } - - - void PetscSolver::checkFeSpaces(vector& feSpaces) - { - FUNCNAME("PetscSolver::checkFeSpaces()"); - - for (unsigned int i = 0; i < feSpaces.size(); i++) { - bool found = false; - for (unsigned int j = 0; j < meshDistributor->getFeSpaces().size(); j++) - if (feSpaces[i] == meshDistributor->getFeSpace(j)) { - found = true; - break; - } - TEST_EXIT(found)("FE spaces not found in mesh distributor!\n"); - } - } - } diff --git a/AMDiS/src/parallel/PetscSolver.h b/AMDiS/src/parallel/PetscSolver.h index b745df17..27b02147 100644 --- a/AMDiS/src/parallel/PetscSolver.h +++ b/AMDiS/src/parallel/PetscSolver.h @@ -241,12 +241,6 @@ namespace AMDiS { void copyVec(Vec& originVec, Vec& destVec, vector& originIndex, vector& destIndex); - /// Checks if all given FE spaces are handled by the mesh distributor. - void checkFeSpaces(vector& feSpaces); - - /// Returns a vector with the FE spaces of all components of a system vector. - vector getFeSpaces(SystemVector *vec); - protected: MeshDistributor *meshDistributor; diff --git a/AMDiS/src/parallel/PetscSolverFeti.cc b/AMDiS/src/parallel/PetscSolverFeti.cc index e5b21db5..e788279e 100644 --- a/AMDiS/src/parallel/PetscSolverFeti.cc +++ b/AMDiS/src/parallel/PetscSolverFeti.cc @@ -1265,7 +1265,7 @@ namespace AMDiS { // === Create scatter to get solutions of all primal nodes that are === // === contained in rank's domain. === - vector feSpaces = getFeSpaces(&vec); + vector feSpaces = vec.getComponentFeSpaces(); vector globalIsIndex, localIsIndex; globalIsIndex.reserve(primalDofMap.getLocalDofs()); @@ -1352,7 +1352,7 @@ namespace AMDiS { // === Create all sets and indices. === - vector feSpaces = AMDiS::getFeSpaces(*mat); + vector feSpaces = AMDiS::getComponentFeSpaces(*mat); initialize(feSpaces); diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc index ceb79b88..ef8f7115 100644 --- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc +++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.cc @@ -68,11 +68,6 @@ namespace AMDiS { MSG(" nz unneeded: %d\n", static_cast(matInfo.nz_unneeded)); } - // === Remove Dirichlet BC DOFs. === - - // removeDirichletBcDofs(mat); - - // === Init PETSc solver. === KSPCreate(mpiCommGlobal, &kspInterior); @@ -114,7 +109,7 @@ namespace AMDiS { TEST_EXIT_DBG(coarseSpaceMap.size() == seqMat->getSize()) ("Wrong sizes %d %d\n", coarseSpaceMap.size(), seqMat->getSize()); - vector feSpaces = AMDiS::getFeSpaces(*seqMat); + vector feSpaces = AMDiS::getComponentFeSpaces(*seqMat); int nRowsRankInterior = interiorMap->getRankDofs(); int nRowsOverallInterior = interiorMap->getOverallDofs(); @@ -748,16 +743,4 @@ namespace AMDiS { } } - - void PetscSolverGlobalMatrix::removeDirichletBcDofs(Matrix *mat) - { - FUNCNAME("PetscSolverGlobalMatrix::removeDirichletBcDofs()"); - } - - - void PetscSolverGlobalMatrix::removeDirichletBcDofs(SystemVector *vec) - { - FUNCNAME("PetscSolverGlobalMatrix::removeDirichletBcDofs()"); - } - } diff --git a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h index de602bef..b7761785 100644 --- a/AMDiS/src/parallel/PetscSolverGlobalMatrix.h +++ b/AMDiS/src/parallel/PetscSolverGlobalMatrix.h @@ -79,10 +79,6 @@ namespace AMDiS { DOFVector* vec, int nRowVec, bool rankOnly = false); - void removeDirichletBcDofs(Matrix *mat); - - void removeDirichletBcDofs(SystemVector *vec); - inline void setDofVector(Vec vecInterior, DOFVector* vec, int nRowVec, bool rankOnly = false) -- GitLab