diff --git a/AMDiS/src/Boundary.h b/AMDiS/src/Boundary.h index bf3567de13032d9c82bdfa018755fdc0c2cf2a60..022d3153601da6543e6179c2ac114a81abe752f9 100644 --- a/AMDiS/src/Boundary.h +++ b/AMDiS/src/Boundary.h @@ -45,116 +45,6 @@ namespace AMDiS { */ typedef signed int BoundaryType; - // /** \ingroup Triangulation - // * \brief - // * Holds information about the type of boundary associated to an edge/face, - // * and how new vertices are projected to the boundary in the case of curved - // * boundaries. - // - // class Boundary - // { - // public: - // MEMORY_MANAGED(Boundary); - - // /** \brief - // * constructor - // - // Boundary(BoundaryType type=0) { - // bound = type; - // }; - - // /** \brief - // * copy constructor - // - // Boundary(const Boundary& old) { bound = old.getBound(); }; - - // /** \brief - // * destructor - // - // virtual ~Boundary() {}; - - // /** \brief - // * assignment operator - // - // Boundary& operator=(const Boundary& old) { - // if (this!=&old) bound = old.getBound(); - // return *this; - // }; - - // /** \brief - // * Returns - // * -true: if a new vertex should be projected to a curved boundary - // * -false: otherwise - // - // virtual bool interpolateBoundary() { return false; }; - - // /** \brief - // * Projection to the curved boundary - // - // virtual void interpolateBoundary(WorldVector<double>& ) {}; - - // /** \brief - // * Returns \ref bound - // - // inline const BoundaryType getBound() const { return bound; }; - - // /** \brief - // * Returns - // * -true: is \ref bound is INTERIOR - // * -false: otherwise - // - // inline const bool isInterior() const {return (bound == INTERIOR);}; - - // /** \brief - // * Returns - // * -true: is \ref bound is DIRICHLET - // * -false: otherwise - // - // inline const bool isDirichlet() const {return (bound == DIRICHLET);}; - - // /** \brief - // * Returns - // * -true: is \ref bound is NEUMANN - // * -false: otherwise - // - // inline const bool isNeumann() const {return (bound == NEUMANN);}; - - // /** \brief - // * Returns - // * -true: is \ref bound is ROBIN - // * -false: otherwise - // - // inline const bool isRobin() const {return (bound == ROBIN);}; - - // /** \brief - // * Returns the new value of \ref bound with respect to its old value and - // * the value of bound_. - // - // BoundaryType newVal(const BoundaryType bound_); - - // /** \brief - // * Returns the Boundary with given type from \ref boundaryMap. - // - // static Boundary* getBoundary(BoundaryType type); - - // /** \brief - // * Adds Boundary b to \ref boundaryMap. - // - // //static void addBoundary(Boundary *b); - - // protected: - // /** \brief - // * type of this boundary - // - // BoundaryType bound; - - // protected: - // /** \brief - // * stl map of all existing boundaries. - // - // static ::std::map<BoundaryType, Boundary*> boundaryMap; - // }; - BoundaryType newBound(BoundaryType oldBound, BoundaryType newBound); } diff --git a/AMDiS/src/BoundaryCondition.h b/AMDiS/src/BoundaryCondition.h index d07e2d6385a3ff761630d0115bd856b42853ea5f..3f9d85dfe7f994de51abd0c80bbc9e437bf52ba2 100644 --- a/AMDiS/src/BoundaryCondition.h +++ b/AMDiS/src/BoundaryCondition.h @@ -43,7 +43,7 @@ namespace AMDiS { * Sub class of BoundaryCondition. Local boundary conditions are filled * while mesh traversal. */ - class BoundaryCondition //: public BoundaryCondition + class BoundaryCondition { public: /** \brief @@ -121,7 +121,9 @@ namespace AMDiS { */ virtual double boundResidual(ElInfo *elInfo, DOFMatrix *matrix, - const DOFVectorBase<double> *dv) { return 0.0; }; + const DOFVectorBase<double> *dv) { + return 0.0; + }; /** \brief * Returns whether the condition must be treated as dirichlet condition diff --git a/AMDiS/src/BoundaryManager.cc b/AMDiS/src/BoundaryManager.cc index 0a515bbd9397038754c5f75047fe486165989906..7baa7276cdcec95c6751317be6696389044930fc 100644 --- a/AMDiS/src/BoundaryManager.cc +++ b/AMDiS/src/BoundaryManager.cc @@ -15,8 +15,8 @@ namespace AMDiS { { double result = 0; ::std::map<BoundaryType, BoundaryCondition*>::iterator it; - for(it = localBCs.begin(); it != localBCs.end(); ++it) { - if((*it).second) + for (it = localBCs.begin(); it != localBCs.end(); ++it) { + if ((*it).second) result += (*it).second->boundResidual(elInfo, matrix, dv); } return result; @@ -35,7 +35,7 @@ namespace AMDiS { ::std::map<BoundaryType, BoundaryCondition*>::iterator it; - if(localBCs.size() > 0) { + if (localBCs.size() > 0) { // get boundaries of all DOFs localBound = basisFcts->getBound(elInfo, NULL); @@ -45,18 +45,18 @@ namespace AMDiS { admin, NULL); // apply non dirichlet boundary conditions - for(it = localBCs.begin(); it != localBCs.end(); ++it) { - if((*it).second) { - if(!(*it).second->isDirichlet()) { + for (it = localBCs.begin(); it != localBCs.end(); ++it) { + if ((*it).second) { + if (!(*it).second->isDirichlet()) { (*it).second->fillBoundaryCondition(vec, elInfo, dofIndices, localBound, nBasFcts); } } } // apply dirichlet boundary conditions - for(it = localBCs.begin(); it != localBCs.end(); ++it) { - if((*it).second) { - if((*it).second->isDirichlet()) { + for (it = localBCs.begin(); it != localBCs.end(); ++it) { + if ((*it).second) { + if ((*it).second->isDirichlet()) { (*it).second->fillBoundaryCondition(vec, elInfo, dofIndices, localBound, nBasFcts); } } @@ -97,7 +97,7 @@ namespace AMDiS { for (it = localBCs.begin(); it != localBCs.end(); ++it) { if ((*it).second) { if ((*it).second->isDirichlet()) { - (*it).second->fillBoundaryCondition(mat, elInfo, dofIndices, localBound, nBasFcts); + (*it).second->fillBoundaryCondition(mat, elInfo, dofIndices, localBound, nBasFcts); } } } @@ -165,16 +165,16 @@ namespace AMDiS { void BoundaryManager::exitVector(DOFVectorBase<double> *vector) { ::std::map<BoundaryType, BoundaryCondition*>::iterator it; - for(it = localBCs.begin(); it != localBCs.end(); ++it) { - if((*it).second) { - if(!(*it).second->isDirichlet()) { + for (it = localBCs.begin(); it != localBCs.end(); ++it) { + if ((*it).second) { + if (!(*it).second->isDirichlet()) { (*it).second->exitVector(vector); } } } - for(it = localBCs.begin(); it != localBCs.end(); ++it) { - if((*it).second) { - if((*it).second->isDirichlet()) { + for (it = localBCs.begin(); it != localBCs.end(); ++it) { + if ((*it).second) { + if ((*it).second->isDirichlet()) { (*it).second->exitVector(vector); } } diff --git a/AMDiS/src/DOFMatrix.cc b/AMDiS/src/DOFMatrix.cc index af868e410538632da3aa12337a9acc6930a38266..9a8beef21397913425f2eb7230a48e1fd21a30fc 100644 --- a/AMDiS/src/DOFMatrix.cc +++ b/AMDiS/src/DOFMatrix.cc @@ -250,7 +250,7 @@ namespace AMDiS { if (condition && condition->isDirichlet()) { MatrixRow *matrixRow = &(matrix[row]); - if(coupleMatrix) { + if (coupleMatrix) { matrixRow->resize(0); } else { matrixRow->resize(1); @@ -447,9 +447,9 @@ namespace AMDiS { ElementMatrix *DOFMatrix::assemble(double factor, ElInfo *elInfo, const BoundaryType *bound, Operator *op) { - FUNCNAME("DOFMatrix::assemble"); + FUNCNAME("DOFMatrix::assemble()"); - if(!op && operators.size() == 0) { + if (!op && operators.size() == 0) { //WARNING("no operator\n"); return NULL; } @@ -459,12 +459,12 @@ namespace AMDiS { elementMatrix = operat->getAssembler()->initElementMatrix(elementMatrix, elInfo); - if(op) { + if (op) { op->getElementMatrix(elInfo, elementMatrix); } else { ::std::vector<Operator*>::iterator it; ::std::vector<double*>::iterator factorIt; - for(it = operators.begin(), factorIt = operatorFactor.begin(); + for (it = operators.begin(), factorIt = operatorFactor.begin(); it != operators.end(); ++it, ++factorIt) { diff --git a/AMDiS/src/DirichletBC.cc b/AMDiS/src/DirichletBC.cc index 24658e6335269610bf721d0f22d594bd2d906122..3b4a924f953db31839fe81122f93d309e0cced58 100644 --- a/AMDiS/src/DirichletBC.cc +++ b/AMDiS/src/DirichletBC.cc @@ -8,7 +8,9 @@ namespace AMDiS { DirichletBC::DirichletBC(BoundaryType type, DOFVectorBase<double> *vec) - : BoundaryCondition(type, vec->getFESpace()), f(NULL), dofVec(vec) + : BoundaryCondition(type, vec->getFESpace()), + f(NULL), + dofVec(vec) {} void DirichletBC::fillBoundaryCondition(DOFMatrix* matrix, @@ -18,10 +20,9 @@ namespace AMDiS { int nBasFcts) { FUNCNAME("DirichletBC::fillBoundaryCondition()"); + TEST_EXIT(matrix->getRowFESpace() == rowFESpace) ("invalid row fe space\n"); - TEST_EXIT(matrix->getColFESpace() == colFESpace) - ("invalid col fe space\n"); } void DirichletBC::fillBoundaryCondition(DOFVectorBase<double>* vector, @@ -31,20 +32,20 @@ namespace AMDiS { int nBasFcts) { FUNCNAME("DirichletBC::fillBoundaryCondition()"); + TEST_EXIT(vector->getFESpace() == rowFESpace) ("invalid row fe space\n"); const BasisFunction *basFcts = rowFESpace->getBasisFcts(); - int i; - for(i=0; i < nBasFcts; i++) { - if(localBound[i] == boundaryType) { - if(f) { + for (int i = 0; i < nBasFcts; i++) { + if (localBound[i] == boundaryType) { + if (f) { DimVec<double> *coords = basFcts->getCoords(i); const WorldVector<double> *worldCoords = elInfo->coordToWorld(*coords, NULL); double fAtCoords = f ? (*f)(*worldCoords) : 0; (*vector)[dofIndices[i]] = fAtCoords; } - if(dofVec) { + if (dofVec) { (*vector)[dofIndices[i]] = (*dofVec)[dofIndices[i]]; } } diff --git a/AMDiS/src/DirichletBC.h b/AMDiS/src/DirichletBC.h index 040930b9a8c3297d88138413db6cbe6035be9ee8..4079124e51eef8b059a1e30afb0080d113e47a49 100644 --- a/AMDiS/src/DirichletBC.h +++ b/AMDiS/src/DirichletBC.h @@ -51,14 +51,17 @@ namespace AMDiS { */ DirichletBC(BoundaryType type, AbstractFunction<double, WorldVector<double> > *fct, - FiniteElemSpace *rowFESpace_) - : BoundaryCondition(type, rowFESpace_, NULL), f(fct), dofVec(NULL) + FiniteElemSpace *rowFESpace, + FiniteElemSpace *colFESpace = NULL) + : BoundaryCondition(type, rowFESpace, colFESpace), + f(fct), + dofVec(NULL) {}; /** \brief * Constructor. */ - DirichletBC(BoundaryType type, + DirichletBC(BoundaryType type, DOFVectorBase<double> *vec); /** \brief @@ -84,9 +87,13 @@ namespace AMDiS { */ double boundResidual(ElInfo*, DOFMatrix *, - const DOFVectorBase<double>*) { return 0.0; }; + const DOFVectorBase<double>*) { + return 0.0; + }; - bool isDirichlet() { return true; }; + bool isDirichlet() { + return true; + }; inline AbstractFunction<double, WorldVector<double> > *getF() { return f; diff --git a/AMDiS/src/Element.h b/AMDiS/src/Element.h index 6f43acad83f28c86ff3daec30ea1386b07faf54f..fa7985e6fa78aee84f0e9158d99ada41e93a1a9b 100644 --- a/AMDiS/src/Element.h +++ b/AMDiS/src/Element.h @@ -134,17 +134,23 @@ namespace AMDiS { * Returns true if Element is a leaf element (\ref child[0] == NULL), returns * false otherwise. */ - inline const bool isLeaf() const { return (child[0]==NULL); }; + inline const bool isLeaf() const { + return (child[0]==NULL); + }; /** \brief * Returns \ref dof[i][j] which is the j-th DOF of the i-th node of Element. */ - const DegreeOfFreedom getDOF(int i,int j) const { return dof[i][j];}; + const DegreeOfFreedom getDOF(int i, int j) const { + return dof[i][j]; + }; /** \brief * Returns \ref dof[i] which is a pointer to the DOFs of the i-th node. */ - const DegreeOfFreedom* getDOF(int i) const {return dof[i];}; + const DegreeOfFreedom* getDOF(int i) const { + return dof[i]; + }; /** \brief * Returns a pointer to the DOFs of this Element @@ -156,7 +162,9 @@ namespace AMDiS { /** \brief * Returns \ref mesh of Element */ - inline Mesh* getMesh() const { return mesh; }; + inline Mesh* getMesh() const { + return mesh; + }; /** \brief * Returns \ref elementData's error estimation, if Element is a leaf element diff --git a/AMDiS/src/GMResSolver.hh b/AMDiS/src/GMResSolver.hh index 876179bd93fe1fa7f67c149738c026095ce85313..fec30262d0a605ef9edc8a0fb8b5dbecee428221 100644 --- a/AMDiS/src/GMResSolver.hh +++ b/AMDiS/src/GMResSolver.hh @@ -311,6 +311,14 @@ namespace AMDiS { double old_res = -1.0; +// ::std::cout << "TEST!" << ::std::endl; + +// matVec->matVec(NoTranspose, *x, *r); +// *r -= *b; +// r->print(); + + + if (norm(b) < TOL) { INFO(this->info, 2)("b == 0, x = 0 is the solution of the linear system\n"); setValue(*x, 0.0); diff --git a/AMDiS/src/Operator.cc b/AMDiS/src/Operator.cc index 499b42c57e752b7713d19ad15f4039bdddb8eee3..2d1eaed2b3d021c0248e27b0b04c17b84b92b783 100644 --- a/AMDiS/src/Operator.cc +++ b/AMDiS/src/Operator.cc @@ -216,7 +216,7 @@ namespace AMDiS { ElementMatrix *userMat, double factor) { - if(!assembler) { + if (!assembler) { initAssembler(NULL, NULL, NULL, NULL); } @@ -227,7 +227,7 @@ namespace AMDiS { ElementVector *userVec, double factor) { - if(!assembler) { + if (!assembler) { initAssembler(NULL, NULL, NULL, NULL); } diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc index 89e17255465b8e4f336d6a7d8ebbd75e8a248bf0..d154f6b2b63324c08e823e1d14271642120dad64 100644 --- a/AMDiS/src/ProblemVec.cc +++ b/AMDiS/src/ProblemVec.cc @@ -809,7 +809,7 @@ namespace AMDiS { { FUNCNAME("ProblemVec::addMatrixOperator()"); - if(!(*systemMatrix_)[i][j]) { + if (!(*systemMatrix_)[i][j]) { TEST_EXIT(i != j)("should have been created already\n"); (*systemMatrix_)[i][j] = NEW DOFMatrix(componentSpaces_[i], componentSpaces_[j], @@ -845,8 +845,10 @@ namespace AMDiS { (*systemMatrix_)[system][i]->getBoundaryManager()->addBoundaryCondition(dirichlet); } } + if (rhs_) rhs_->getDOFVector(system)->getBoundaryManager()->addBoundaryCondition(dirichlet); + if (solution_) solution_->getDOFVector(system)->getBoundaryManager()->addBoundaryCondition(dirichlet); } @@ -860,7 +862,7 @@ namespace AMDiS { new NeumannBC(type, n, componentSpaces_[row], componentSpaces_[col]); - if(rhs_) + if (rhs_) rhs_->getDOFVector(row)->getBoundaryManager()->addBoundaryCondition(neumann); } @@ -874,9 +876,10 @@ namespace AMDiS { new RobinBC(type, n, r, componentSpaces_[row], componentSpaces_[col]); - if(rhs_) + if (rhs_) rhs_->getDOFVector(row)->getBoundaryManager()->addBoundaryCondition(robin); - if(systemMatrix_ && (*systemMatrix_)[row][col]) { + + if (systemMatrix_ && (*systemMatrix_)[row][col]) { (*systemMatrix_)[row][col]->getBoundaryManager()->addBoundaryCondition(robin); } } diff --git a/AMDiS/src/ProblemVec.h b/AMDiS/src/ProblemVec.h index 11222f5b529434820e0f036c788d3d84cd34cdec..c27a8223c5f946ae0b7362ad15df87ae71d20370 100644 --- a/AMDiS/src/ProblemVec.h +++ b/AMDiS/src/ProblemVec.h @@ -446,7 +446,9 @@ namespace AMDiS { /** \brief * Sets \ref solver_. */ - inline void setSolver(OEMSolver<SystemVector>* sol) { solver_ = sol; }; + inline void setSolver(OEMSolver<SystemVector>* sol) { + solver_ = sol; + }; /** \brief * Sets \ref leftPrecon_. diff --git a/AMDiS/src/SystemVector.h b/AMDiS/src/SystemVector.h index 3048f31b779987804d27c09bfc7eadc60064cc3d..acca3e7b044174630008e4dc448b6b00c81beb62 100644 --- a/AMDiS/src/SystemVector.h +++ b/AMDiS/src/SystemVector.h @@ -139,12 +139,7 @@ namespace AMDiS { } }; - virtual ~SystemVector() { - // int i, size = vectors.getSize(); - // for(i = 0; i < size; i++) { - // if(vectors[i]) DELETE vectors[i]; - // } - }; + virtual ~SystemVector() {}; /** \brief * Sets \ref vectors[index] = vec. @@ -298,8 +293,8 @@ namespace AMDiS { }; void print() { - int i, size = vectors.getSize(); - for(i = 0; i < size; i++) { + int size = vectors.getSize(); + for (int i = 0; i < size; i++) { vectors[i]->print(); } };