diff --git a/AMDiS/src/DOFVector.cc b/AMDiS/src/DOFVector.cc index d79dbf32d09fa6bebd6363c10e7a1db6ff8ee1c3..96ab5a49b39ee90e3755ebf7655f13bf80ffc0ef 100644 --- a/AMDiS/src/DOFVector.cc +++ b/AMDiS/src/DOFVector.cc @@ -243,7 +243,6 @@ namespace AMDiS { ("invalid basis functions"); int dow = Global::getGeo(WORLD); - int myRank = omp_get_thread_num(); int nPoints = quadFast ? quadFast->getQuadrature()->getNumPoints() : quad->getNumPoints(); static WorldVector<double> *grd = NULL; WorldVector<double> *result; @@ -261,8 +260,8 @@ namespace AMDiS { result = grd; } - const ElementVector& localVec = localVectors[myRank]; - getLocalVector(elInfo->getElement(), const_cast<ElementVector&>(localVec)); + ElementVector localVec(nBasFcts); + getLocalVector(elInfo->getElement(), localVec); mtl::dense_vector<double> grd1(dim + 1); int parts = Global::getGeo(PARTS, dim); @@ -327,7 +326,6 @@ namespace AMDiS { ("invalid basis functions"); int dow = Global::getGeo(WORLD); - int myRank = omp_get_thread_num(); int nPoints = quadFast ? quadFast->getQuadrature()->getNumPoints() : quad->getNumPoints(); static WorldVector<double> *grd = NULL; WorldVector<double> *result; @@ -345,8 +343,8 @@ namespace AMDiS { result = grd; } - const ElementVector& localVec = localVectors[myRank]; - getLocalVector(largeElInfo->getElement(), const_cast<ElementVector&>(localVec)); + ElementVector localVec(nBasFcts); + getLocalVector(largeElInfo->getElement(), localVec); const BasisFunction *basFcts = feSpace->getBasisFcts(); mtl::dense2D<double> &m = smallElInfo->getSubElemCoordsMat(basFcts->getDegree()); @@ -405,7 +403,6 @@ namespace AMDiS { Element *el = elInfo->getElement(); - int myRank = omp_get_thread_num(); int dow = Global::getGeo(WORLD); int nPoints = quadFast ? quadFast->getQuadrature()->getNumPoints() : quad->getNumPoints(); @@ -423,8 +420,8 @@ namespace AMDiS { result = vec; } - const ElementVector& localVec = localVectors[myRank]; - getLocalVector(el, const_cast<ElementVector&>(localVec)); + ElementVector localVec(nBasFcts); + getLocalVector(el, localVec); DimMat<double> D2Tmp(dim, DEFAULT_VALUE, 0.0); int parts = Global::getGeo(PARTS, dim); @@ -452,7 +449,7 @@ namespace AMDiS { } } else { const BasisFunction *basFcts = feSpace->getBasisFcts(); - DimMat<double>* D2Phi = D2Phis[omp_get_thread_num()]; + DimMat<double> D2Phi(dim, NO_INIT); for (int iq = 0; iq < nPoints; iq++) { for (int k = 0; k < parts; k++) @@ -461,11 +458,11 @@ namespace AMDiS { for (int i = 0; i < nBasFcts; i++) { WARNING("not tested after index correction\n"); - (*(basFcts->getD2Phi(i)))(quad->getLambda(iq), *D2Phi); + (*(basFcts->getD2Phi(i)))(quad->getLambda(iq), D2Phi); for (int k = 0; k < parts; k++) for (int l = 0; l < parts; l++) - D2Tmp[k][l] += localVec[i] * (*D2Phi)[k][l]; + D2Tmp[k][l] += localVec[i] * D2Phi[k][l]; } for (int i = 0; i < dow; i++) @@ -496,7 +493,7 @@ namespace AMDiS { this->set(0.0); - DegreeOfFreedom *myLocalIndices = localIndices[omp_get_thread_num()]; + std::vector<DegreeOfFreedom> myLocalIndices(nBasisFcts); ElementVector sourceLocalCoeffs(nSourceBasisFcts); if (feSpace->getMesh() == sourceFeSpace->getMesh()) { @@ -582,11 +579,11 @@ namespace AMDiS { const BasisFunction *basFcts = feSpace->getBasisFcts(); const BasisFunction *vBasFcts = vFeSpace->getBasisFcts(); - int numBasFcts = basFcts->getNumber(); + int nBasFcts = basFcts->getNumber(); int vNumBasFcts = vBasFcts->getNumber(); - if (feSpace->getMesh() == vFeSpace->getMesh()) { - DegreeOfFreedom *myLocalIndices = localIndices[omp_get_thread_num()]; + if (feSpace->getMesh() == vFeSpace->getMesh()) { + std::vector<DegreeOfFreedom> myLocalIndices(nBasFcts); mtl::dense_vector<WorldVector<double> > vLocalCoeffs(vNumBasFcts); Mesh *mesh = feSpace->getMesh(); TraverseStack stack; @@ -598,7 +595,7 @@ namespace AMDiS { basFcts->getLocalIndices(el, feSpace->getAdmin(), myLocalIndices); v->getLocalVector(el, vLocalCoeffs); - for (int i = 0; i < numBasFcts; i++) { + for (int i = 0; i < nBasFcts; i++) { if (vec[myLocalIndices[i]] == nul) { coords = basFcts->getCoords(i); vec[myLocalIndices[i]] += vBasFcts->evalUh(*coords, vLocalCoeffs, NULL) * factor; @@ -776,9 +773,8 @@ namespace AMDiS { vecAtQPs.change_dim(nPoints); - const ElementVector &localVec = localVectors[omp_get_thread_num()]; - getLocalVector(largeElInfo->getElement(), const_cast<ElementVector&>(localVec)); - + ElementVector localVec(nBasFcts); + getLocalVector(largeElInfo->getElement(), localVec); mtl::dense2D<double> &m = smallElInfo->getSubElemCoordsMat(basFcts->getDegree()); for (int i = 0; i < nPoints; i++) { diff --git a/AMDiS/src/DOFVector.h b/AMDiS/src/DOFVector.h index e93d0bc10f0e118254bf09e77379c05e2df0df07..9488d14814ddcc6b9b5d505d70da9658bb2dd333 100644 --- a/AMDiS/src/DOFVector.h +++ b/AMDiS/src/DOFVector.h @@ -64,9 +64,6 @@ namespace AMDiS { virtual ~DOFVectorBase(); - /// Creates temporary used data structures. - void createTempData(); - /** \brief * For the given element, this function returns an array of all DOFs of this * DOFVector that are defined on this element. @@ -263,15 +260,6 @@ namespace AMDiS { /// Number of basis functions of the used finite element space. int nBasFcts; - /// Are used to store temporary local dofs of an element. Thread safe. - std::vector<DegreeOfFreedom*> localIndices; - - /// Are used to store temporary local values of an element. Thread safe. - std::vector<mtl::dense_vector<T> > localVectors; - - /// Temporarly used in \ref getD2AtQPs. Thread safe. - std::vector<DimMat<double>*> D2Phis; - /// Dimension of the mesh this DOFVectorBase belongs to int dim; @@ -359,8 +347,6 @@ namespace AMDiS { this->name = rhs.name + "copy"; if (this->feSpace && this->feSpace->getAdmin()) (dynamic_cast<DOFAdmin*>(this->feSpace->getAdmin()))->addDOFIndexed(this); - - this->createTempData(); } /// Destructor diff --git a/AMDiS/src/DOFVector.hh b/AMDiS/src/DOFVector.hh index 57d7eb46039f10e14b952597098e5a4db8730604..9cc5b51c917807b30ab51554845e89cf2466d4c0 100644 --- a/AMDiS/src/DOFVector.hh +++ b/AMDiS/src/DOFVector.hh @@ -31,7 +31,6 @@ #include "AbstractFunction.h" #include "BoundaryManager.h" #include "Assembler.h" -#include "OpenMP.h" #include "Operator.h" #include "Parameters.h" #include "Traverse.h" @@ -88,36 +87,14 @@ namespace AMDiS { { nBasFcts = feSpace->getBasisFcts()->getNumber(); dim = feSpace->getMesh()->getDim(); - - createTempData(); } template<typename T> DOFVectorBase<T>::~DOFVectorBase() - { - for (int i = 0; i < static_cast<int>(localIndices.size()); i++) { - delete [] localIndices[i]; - delete D2Phis[i]; - } - } + {} - - template<typename T> - void DOFVectorBase<T>::createTempData() - { - localIndices.resize(omp_get_overall_max_threads()); - localVectors.resize(omp_get_overall_max_threads()); - D2Phis.resize(omp_get_overall_max_threads()); - - for (int i = 0; i < omp_get_overall_max_threads(); i++) { - localIndices[i] = new DegreeOfFreedom[this->nBasFcts]; - localVectors[i].change_dim(this->nBasFcts); - D2Phis[i] = new DimMat<double>(dim, NO_INIT); - } - } - - + template<typename T> DOFVector<T>::DOFVector(const FiniteElemSpace* f, std::string n) : DOFVectorBase<T>(f, n), @@ -465,9 +442,10 @@ namespace AMDiS { const_cast<BasisFunction*>(basFct)->interpol(elinfo, 0, NULL, traverseVector->interFct, NULL); - DegreeOfFreedom *myLocalIndices = this->localIndices[omp_get_thread_num()]; - basFct->getLocalIndices(const_cast<Element*>(elinfo->getElement()), admin, myLocalIndices); int nBasFcts = basFct->getNumber(); + std::vector<DegreeOfFreedom> myLocalIndices(nBasFcts); + basFct->getLocalIndices(const_cast<Element*>(elinfo->getElement()), + admin, myLocalIndices); for (int i = 0; i < nBasFcts; i++) (*traverseVector)[myLocalIndices[i]] = inter_val[i]; } @@ -991,19 +969,14 @@ namespace AMDiS { ("Element is defined on a different mesh than the DOF vector!\n"); const DOFAdmin* admin = feSpace->getAdmin(); - -#ifdef _OPENMP - std::vector<DegreeOfFreedom> localIndices(nBasFcts); - feSpace->getBasisFcts()->getLocalIndices(el, admin, &(localIndices[0])); -#else const DegreeOfFreedom *localIndices = - feSpace->getBasisFcts()->getLocalIndices(el, admin, NULL); -#endif + feSpace->getBasisFcts()->getLocalIndices(el, admin, NULL); for (int i = 0; i < nBasFcts; i++) d[i] = (*this)[localIndices[i]]; } + template<typename T> void DOFVectorBase<T>::getVecAtQPs(const ElInfo *elInfo, const Quadrature *quad, @@ -1019,9 +992,10 @@ namespace AMDiS { TEST_EXIT_DBG(!quadFast || quadFast->getBasisFunctions() == feSpace->getBasisFcts()) ("Invalid basis functions!"); - const mtl::dense_vector<T>& localVec = localVectors[omp_get_thread_num()]; - getLocalVector(elInfo->getElement(), - const_cast<mtl::dense_vector<T>& >(localVec)); + const BasisFunction *basFcts = feSpace->getBasisFcts(); + int nBasFcts = basFcts->getNumber(); + mtl::dense_vector<T> localVec(nBasFcts); + getLocalVector(elInfo->getElement(), localVec); if (quadFast) { const mtl::dense2D<double>& phi = quadFast->getPhi(); @@ -1029,9 +1003,7 @@ namespace AMDiS { vecAtQPs = phi * localVec; } else { const Quadrature *quadrature = quadFast ? quadFast->getQuadrature() : quad; - const BasisFunction *basFcts = feSpace->getBasisFcts(); int nPoints = quadrature->getNumPoints(); - int nBasFcts = basFcts->getNumber(); vecAtQPs.change_dim(nPoints); for (int i = 0; i < nPoints; i++) { diff --git a/AMDiS/src/FirstOrderAssembler.cc b/AMDiS/src/FirstOrderAssembler.cc index 40708e05a75896cdeda9ffcb5772baf2e976dbc3..f811efe13305ccc22a392781a237017031349c3d 100644 --- a/AMDiS/src/FirstOrderAssembler.cc +++ b/AMDiS/src/FirstOrderAssembler.cc @@ -207,7 +207,7 @@ namespace AMDiS { for (unsigned int i = 0; i < terms.size(); i++) (static_cast<FirstOrderTerm*>((terms[i])))->getLb(elInfo, Lb); - const mtl::dense2D<double>& phi(phiFast->getPhi()); + const mtl::dense2D<double>& phi = phiFast->getPhi(); for (int iq = 0; iq < nPoints; iq++) { Lb[iq] *= elInfo->getDet(); @@ -220,8 +220,6 @@ namespace AMDiS { mat[i][j] += factor * phi[iq][j] * dot(Lb[iq], grdPsi[i]); } } - - // TODO: Do the same performance update es in Quad01::calculateElementMatrix }