From 715117bf73e220252d67f6cbdcc8975a75d3b195 Mon Sep 17 00:00:00 2001 From: Thomas Witkowski <thomas.witkowski@gmx.de> Date: Fri, 13 Nov 2009 14:20:32 +0000 Subject: [PATCH] Removed several static functions from DOFVector. --- AMDiS/src/DOFVector.h | 23 ---- AMDiS/src/DOFVector.hh | 237 +++++++++++++++++++---------------------- 2 files changed, 111 insertions(+), 149 deletions(-) diff --git a/AMDiS/src/DOFVector.h b/AMDiS/src/DOFVector.h index 32a75d1f..89c7df7e 100644 --- a/AMDiS/src/DOFVector.h +++ b/AMDiS/src/DOFVector.h @@ -562,22 +562,6 @@ namespace AMDiS { DOFVector<WorldVector<T> >* getRecoveryGradient(DOFVector<WorldVector<T> >*) const; - protected: - /// Used by Int while mesh traversal - static int Int_fct(ElInfo* elinfo); - - /// Used by L1Norm while mesh traversal - static int L1Norm_fct(ElInfo* elinfo); - - /// Used by L2Norm while mesh traversal - static int L2Norm_fct(ElInfo* elinfo); - - /// Used by H1Norm while mesh traversal - static int H1Norm_fct(ElInfo* elinfo); - - /// Used by DoubleWell while mesh traversal - static int DoubleWell_fct(ElInfo* elinfo); - protected: /// Data container std::vector<T> vec; @@ -590,13 +574,6 @@ namespace AMDiS { /// Used for mesh traversal static DOFVector<T> *traverseVector; - - protected: - /// Used while calculating vector norms - static FastQuadrature *quad_fast; - - /// Stores the last calculated vector norm - static double norm; }; diff --git a/AMDiS/src/DOFVector.hh b/AMDiS/src/DOFVector.hh index 3feec5cd..e78ea500 100644 --- a/AMDiS/src/DOFVector.hh +++ b/AMDiS/src/DOFVector.hh @@ -138,12 +138,6 @@ namespace AMDiS { template<typename T> DOFVector<T> * DOFVector<T>::traverseVector = NULL; - template<typename T> - FastQuadrature *DOFVector<T>::quad_fast = NULL; - - template<typename T> - double DOFVector<T>::norm = 0.0; - template<typename T> void DOFVectorBase<T>::addElementVector(T factor, const ElementVector &elVec, @@ -455,29 +449,28 @@ namespace AMDiS { q = Quadrature::provideQuadrature(this->dim, deg); } - quad_fast = + FastQuadrature *quadFast = FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), *q, INIT_PHI); - norm = 0.0; - traverseVector = const_cast<DOFVector<T>*>(this); - mesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET, Int_fct); - - return norm; - } - - template<typename T> - int DOFVector<T>::Int_fct(ElInfo *elinfo) - { - double det = elinfo->getDet(); - const T* uh_vec = traverseVector->getVecAtQPs(elinfo, NULL, quad_fast, NULL); - int numPoints = quad_fast->getNumPoints(); - double normT = 0.0; - for (int iq = 0; iq < numPoints; iq++) - normT += quad_fast->getWeight(iq) * (uh_vec[iq]); + double result = 0.0; + int nPoints = quadFast->getNumPoints(); + std::vector<T> uh_vec(this->feSpace->getBasisFcts()->getNumber()); + TraverseStack stack; + ElInfo *elInfo = + stack.traverseFirst(mesh, -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET); + while (elInfo) { + double det = elInfo->getDet(); + double normT = 0.0; + this->getVecAtQPs(elInfo, NULL, quadFast, &(uh_vec[0])); + for (int iq = 0; iq < nPoints; iq++) + normT += quadFast->getWeight(iq) * (uh_vec[iq]); + result += det * normT; - norm += det * normT; + elInfo = stack.traverseNext(elInfo); + } - return 0; + return result; } template<typename T> @@ -485,122 +478,114 @@ namespace AMDiS { { FUNCNAME("DOFVector::L1Norm()"); + Mesh* mesh = this->feSpace->getMesh(); + if (!q) { int deg = 2 * this->feSpace->getBasisFcts()->getDegree(); q = Quadrature::provideQuadrature(this->dim, deg); } - quad_fast = FastQuadrature::provideFastQuadrature(this-> feSpace->getBasisFcts(), - *q, INIT_PHI); - norm = 0.0; - traverseVector = const_cast<DOFVector<T>*>(this); - Mesh* mesh = this->feSpace->getMesh(); - - mesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET, - L1Norm_fct); - - return norm; - } - - template<typename T> - int DOFVector<T>::L1Norm_fct(ElInfo *elinfo) - { - double det = elinfo->getDet(); - const T* uh_vec = traverseVector->getVecAtQPs(elinfo, NULL, quad_fast, NULL); - int numPoints = quad_fast->getNumPoints(); - double normT = 0.0; + FastQuadrature *quadFast = + FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), *q, INIT_PHI); - for (int iq = 0; iq < numPoints; iq++) - normT += quad_fast->getWeight(iq) * abs(uh_vec[iq]); + double result = 0.0; + int nPoints = quadFast->getNumPoints(); + std::vector<T> uh_vec(this->feSpace->getBasisFcts()->getNumber()); + TraverseStack stack; + ElInfo *elInfo = + stack.traverseFirst(mesh, -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET); + while (elInfo) { + double det = elInfo->getDet(); + double normT = 0.0; + this->getVecAtQPs(elInfo, NULL, quadFast, &(uh_vec[0])); + for (int iq = 0; iq < nPoints; iq++) + normT += quadFast->getWeight(iq) * abs(uh_vec[iq]); + result += det * normT; - norm += det * normT; + elInfo = stack.traverseNext(elInfo); + } - return 0; + return result; } - template<typename T> double DOFVector<T>::L2NormSquare(Quadrature* q) const { FUNCNAME("DOFVector::L2NormSquare()"); + Mesh* mesh = this->feSpace->getMesh(); + if (!q) { int deg = 2 * this->feSpace->getBasisFcts()->getDegree(); q = Quadrature::provideQuadrature(this->dim, deg); } - quad_fast = FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), - *q, INIT_PHI); - - norm = 0.0; - traverseVector = const_cast<DOFVector<T>*>(this); - Mesh* mesh = this->feSpace->getMesh(); - mesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET, - L2Norm_fct); - - return norm; - } - - template<typename T> - int DOFVector<T>::L2Norm_fct(ElInfo *elinfo) - { - double det = elinfo->getDet(); - const T *uh_vec = traverseVector->getVecAtQPs(elinfo, NULL, quad_fast, NULL); - int numPoints = quad_fast->getNumPoints(); - double normT = 0.0; - - for (int iq = 0; iq < numPoints; iq++) - normT += quad_fast->getWeight(iq) * sqr(uh_vec[iq]); - - norm += det * normT; + FastQuadrature *quadFast = + FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), *q, INIT_PHI); - return 0; - } + double result = 0.0; + int nPoints = quadFast->getNumPoints(); + std::vector<T> uh_vec(this->feSpace->getBasisFcts()->getNumber()); + TraverseStack stack; + ElInfo *elInfo = + stack.traverseFirst(mesh, -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET); + while (elInfo) { + double det = elInfo->getDet(); + double normT = 0.0; + this->getVecAtQPs(elInfo, NULL, quadFast, &(uh_vec[0])); + for (int iq = 0; iq < nPoints; iq++) + normT += quadFast->getWeight(iq) * sqr(uh_vec[iq]); + result += det * normT; - template<typename T> - int DOFVector<T>::H1Norm_fct(ElInfo *elinfo) - { - double det = elinfo->getDet(); - const WorldVector<T> *grduh_vec = - traverseVector->getGrdAtQPs(elinfo, NULL, quad_fast, NULL); - int dimOfWorld = Global::getGeo(WORLD); - int numPoints = quad_fast->getNumPoints(); - double normT = 0.0; - - for (int iq = 0; iq < numPoints; iq++) { - double norm2 = 0.0; - for (int j = 0; j < dimOfWorld; j++) - norm2 += sqr(grduh_vec[iq][j]); - - normT += quad_fast->getWeight(iq) * norm2; + elInfo = stack.traverseNext(elInfo); } - norm += det * normT; - return 0; + return result; } - - template<typename T> + template<typename T> double DOFVector<T>::H1NormSquare(Quadrature *q) const { FUNCNAME("DOFVector::H1NormSquare()"); + Mesh* mesh = this->feSpace->getMesh(); + if (!q) { int deg = 2 * this->feSpace->getBasisFcts()->getDegree() - 2; q = Quadrature::provideQuadrature(this->dim, deg); } - quad_fast = FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), - *q, INIT_GRD_PHI); + FastQuadrature *quadFast = + FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), *q, INIT_GRD_PHI); + + double result = 0.0; + int nPoints = quadFast->getNumPoints(); + int dimOfWorld = Global::getGeo(WORLD); + std::vector<WorldVector<T> > grduh_vec(this->feSpace->getBasisFcts()->getNumber()); + TraverseStack stack; + ElInfo *elInfo = + stack.traverseFirst(mesh, -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET | Mesh::FILL_GRD_LAMBDA); + while (elInfo) { + double det = elInfo->getDet(); + double normT = 0.0; + this->getGrdAtQPs(elInfo, NULL, quadFast, &(grduh_vec[0])); + + for (int iq = 0; iq < nPoints; iq++) { + double norm2 = 0.0; + for (int j = 0; j < dimOfWorld; j++) + norm2 += sqr(grduh_vec[iq][j]); + normT += quadFast->getWeight(iq) * norm2; + } - norm = 0.0; - traverseVector = const_cast<DOFVector<T>*>(this); - Mesh *mesh = this->feSpace->getMesh(); + result += det * normT; - mesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET | - Mesh::FILL_GRD_LAMBDA, H1Norm_fct); + elInfo = stack.traverseNext(elInfo); + } - return norm; + return result; } template<typename T> @@ -1081,35 +1066,35 @@ namespace AMDiS { double DOFVector<T>::DoubleWell(Quadrature* q) const { FUNCNAME("DOFVector::DoubleWell()"); - + + Mesh* mesh = this->feSpace->getMesh(); + if (!q) { int deg = 2 * this->feSpace->getBasisFcts()->getDegree(); q = Quadrature::provideQuadrature(this->dim, deg); } - quad_fast = FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), - *q, INIT_PHI); - norm = 0.0; - traverseVector = const_cast<DOFVector<T>*>(this); - Mesh* mesh = this->feSpace->getMesh(); - mesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET, - DoubleWell_fct); + FastQuadrature *quadFast = + FastQuadrature::provideFastQuadrature(this->feSpace->getBasisFcts(), *q, INIT_PHI); + + double result = 0.0; + int nPoints = quadFast->getNumPoints(); + std::vector<T> uh_vec(this->feSpace->getBasisFcts()->getNumber()); + TraverseStack stack; + ElInfo *elInfo = + stack.traverseFirst(mesh, -1, + Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_DET); + while (elInfo) { + double det = elInfo->getDet(); + double normT = 0.0; + this->getVecAtQPs(elInfo, NULL, quadFast, &(uh_vec[0])); + for (int iq = 0; iq < nPoints; iq++) + normT += quadFast->getWeight(iq) * sqr(uh_vec[iq]) * sqr(1.0 - uh_vec[iq]); + result += det * normT; - return norm; - } + elInfo = stack.traverseNext(elInfo); + } - template<typename T> - int DOFVector<T>::DoubleWell_fct(ElInfo *elinfo) - { - const T *uh_vec = traverseVector->getVecAtQPs(elinfo, NULL, quad_fast, NULL); - int nPoints = quad_fast->getNumPoints(); - double normT = 0.0; - for (int iq = 0; iq < nPoints; iq++) - normT += quad_fast->getWeight(iq) * sqr(uh_vec[iq]) * sqr(1.0 - uh_vec[iq]); - - norm += elinfo->getDet() * normT; - - return 0; + return result; } - } -- GitLab