diff --git a/AMDiS/src/Assembler.cc b/AMDiS/src/Assembler.cc index ee9efcf6067a653ac3bbe573051a251ffb7a5b76..5d49cb97b5ce25d49a9874bec53e065a8e603855 100644 --- a/AMDiS/src/Assembler.cc +++ b/AMDiS/src/Assembler.cc @@ -38,15 +38,13 @@ namespace AMDiS { { FUNCNAME("Assembler::calculateElementMatrix()"); - if (remember && (factor != 1.0 || operat->uhOld)) { + if (remember && (factor != 1.0 || operat->uhOld)) rememberElMat = true; - } Element *el = elInfo->getElement(); - if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) { + if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) initElement(elInfo); - } if (el != lastMatEl || !operat->isOptimized()) { if (rememberElMat) @@ -84,16 +82,14 @@ namespace AMDiS { { FUNCNAME("Assembler::calculateElementMatrix()"); - if (remember && ((factor != 1.0) || (operat->uhOld))) { + if (remember && ((factor != 1.0) || (operat->uhOld))) rememberElMat = true; - } Element *el = smallElInfo->getElement(); lastVecEl = lastMatEl = NULL; - if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) { + if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) initElement(smallElInfo, largeElInfo); - } if (el != lastMatEl || !operat->isOptimized()) { if (rememberElMat) @@ -149,15 +145,13 @@ namespace AMDiS { { FUNCNAME("Assembler::calculateElementVector()"); - if (remember && factor != 1.0) { + if (remember && factor != 1.0) rememberElVec = true; - } Element *el = elInfo->getElement(); - if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) { + if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) initElement(elInfo); - } if (el != lastVecEl || !operat->isOptimized()) { if (rememberElVec) @@ -174,9 +168,9 @@ namespace AMDiS { ElementVector& vec = rememberElVec ? elementVector : userVec; if (operat->uhOld && remember) { matVecAssemble(elInfo, vec); - if (rememberElVec) { + if (rememberElVec) userVec += factor * elementVector; - } + return; } @@ -184,6 +178,7 @@ namespace AMDiS { firstOrderAssemblerGrdPsi->calculateElementVector(elInfo, vec); if (zeroOrderAssembler) zeroOrderAssembler->calculateElementVector(elInfo, vec); + if (rememberElVec) userVec += factor * elementVector; } @@ -197,15 +192,13 @@ namespace AMDiS { { FUNCNAME("Assembler::calculateElementVector()"); - if (remember && factor != 1.0) { + if (remember && factor != 1.0) rememberElVec = true; - } Element *el = mainElInfo->getElement(); - if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) { + if ((el != lastMatEl && el != lastVecEl) || !operat->isOptimized()) initElement(auxElInfo); - } if (el != lastVecEl || !operat->isOptimized()) { if (rememberElVec) @@ -228,9 +221,8 @@ namespace AMDiS { matVecAssemble(mainElInfo, auxElInfo, smallElInfo, largeElInfo, vec); } - if (rememberElVec) { + if (rememberElVec) userVec += factor * elementVector; - } return; } diff --git a/AMDiS/src/BasisFunction.cc b/AMDiS/src/BasisFunction.cc index 20d5ad0909c5c2610e417ca5947e70dda2164ace..10aef95a3b12c98329bddabe87ec32b433196bee 100644 --- a/AMDiS/src/BasisFunction.cc +++ b/AMDiS/src/BasisFunction.cc @@ -32,7 +32,7 @@ namespace AMDiS { grdTmpVec1[i] = NEW DimVec<double>(dim, DEFAULT_VALUE, 0.0); grdTmpVec2[i] = NEW DimVec<double>(dim, DEFAULT_VALUE, 0.0); } - }; + } BasisFunction::~BasisFunction() { diff --git a/AMDiS/src/CoarseningManager.cc b/AMDiS/src/CoarseningManager.cc index 3a102a2b9736a44bf51b011fcd5b5d3c4dd9dac9..69fbe11f3b323576810944799dada002b58eb573 100644 --- a/AMDiS/src/CoarseningManager.cc +++ b/AMDiS/src/CoarseningManager.cc @@ -40,30 +40,27 @@ namespace AMDiS { Element *el = el_info->getElement(); signed char mark; - if (el->getChild(0)) - { - /****************************************************************************/ - /* interior node of the tree */ - /****************************************************************************/ - mark = max(el->getChild(0)->getMark(), el->getChild(1)->getMark()); - el->setMark(std::min(mark + 1, 0)); - } - else - { - /****************************************************************************/ - /* leaf node of the tree */ - /****************************************************************************/ - if (el->getMark() < 0) el->setMark(el->getMark() - 1); - } + if (el->getChild(0)) { + /****************************************************************************/ + /* interior node of the tree */ + /****************************************************************************/ + mark = max(el->getChild(0)->getMark(), el->getChild(1)->getMark()); + el->setMark(std::min(mark + 1, 0)); + } else { + /****************************************************************************/ + /* leaf node of the tree */ + /****************************************************************************/ + if (el->getMark() < 0) + el->setMark(el->getMark() - 1); + } + return 0; } void CoarseningManager::spreadCoarsenMark() { traversePtr = this; - mesh->traverse(-1, - Mesh::CALL_EVERY_EL_POSTORDER, - spreadCoarsenMarkFunction); + mesh->traverse(-1, Mesh::CALL_EVERY_EL_POSTORDER, spreadCoarsenMarkFunction); } @@ -93,7 +90,6 @@ namespace AMDiS { Flag CoarseningManager::coarsenMesh(Mesh *aMesh) { int n_elements; - Flag flag = Mesh::CALL_EVERY_EL_POSTORDER | Mesh::FILL_NEIGH; ElInfo *el_info; mesh = aMesh; @@ -102,18 +98,23 @@ namespace AMDiS { spreadCoarsenMark(); - stack = NEW TraverseStack; + stack = new TraverseStack; do { doMore = false; - el_info = stack->traverseFirst(mesh, -1, flag); + el_info = stack->traverseFirst(mesh, -1, + Mesh::CALL_EVERY_EL_POSTORDER | Mesh::FILL_NEIGH); while (el_info) { - coarsenFunction(el_info); + int idx = el_info->getElement()->getIndex(); + + // if (idx != 2288 && idx != 2283) + coarsenFunction(el_info); + el_info = stack->traverseNext(el_info); } } while (doMore); - DELETE stack; + delete stack; cleanUpAfterCoarsen(); diff --git a/AMDiS/src/CoarseningManager2d.cc b/AMDiS/src/CoarseningManager2d.cc index e50648cc912054fc6386b17cbe67cbef01c12c6f..fb54607899b8202ed2a772b42c9229e7147f2ad9 100644 --- a/AMDiS/src/CoarseningManager2d.cc +++ b/AMDiS/src/CoarseningManager2d.cc @@ -87,15 +87,15 @@ namespace AMDiS { /****************************************************************************/ /* get new dof on el at the midpoint of the coarsening edge */ /****************************************************************************/ - if (!el->getDOF(node+2)) { - el->setDOF(node+2, mesh->getDOF(EDGE)); + if (!el->getDOF(node + 2)) { + el->setDOF(node + 2, mesh->getDOF(EDGE)); if (neigh) { - neigh->setDOF(node+2, const_cast<int*>( el->getDOF(node+2))); + neigh->setDOF(node + 2, const_cast<int*>(el->getDOF(node+2))); } } } - if (mesh->getNumberOfDOFs(EDGE) || mesh->getNumberOfDOFs(CENTER)) { + if (mesh->getNumberOfDOFs(EDGE) || mesh->getNumberOfDOFs(CENTER)) { coarsenList->addDOFParents(n_neigh); } @@ -132,7 +132,7 @@ namespace AMDiS { int CoarseningManager2d::coarsenFunction(ElInfo *el_info) { - Triangle *el = dynamic_cast<Triangle*>(const_cast<Element*>( el_info->getElement())); + Triangle *el = dynamic_cast<Triangle*>(const_cast<Element*>(el_info->getElement())); DegreeOfFreedom *edge[2]; int n_neigh, bound = 0; RCNeighbourList coarse_list(2); @@ -165,11 +165,11 @@ namespace AMDiS { /****************************************************************************/ if (el->getDOF(0,0) < el->getDOF(1,0)) { - edge[0] = const_cast<int*>( el->getDOF(0)); - edge[1] = const_cast<int*>( el->getDOF(1)); + edge[0] = const_cast<int*>(el->getDOF(0)); + edge[1] = const_cast<int*>(el->getDOF(1)); } else { - edge[1] = const_cast<int*>( el->getDOF(0)); - edge[0] = const_cast<int*>( el->getDOF(1)); + edge[1] = const_cast<int*>(el->getDOF(0)); + edge[0] = const_cast<int*>(el->getDOF(1)); } coarse_list.setElement(0, el, true); diff --git a/AMDiS/src/DOFAdmin.cc b/AMDiS/src/DOFAdmin.cc index e92bb6060deb51dba4a68b6324252a27bc1c360d..06ca719a1fd59c98a272f168b1831c158cd796c6 100755 --- a/AMDiS/src/DOFAdmin.cc +++ b/AMDiS/src/DOFAdmin.cc @@ -89,7 +89,7 @@ namespace AMDiS { ERROR_EXIT("TODO\n"); } - void DOFAdmin::freeDOFIndex(int dof) { + void DOFAdmin::freeDOFIndex(int dof) { FUNCNAME("DOFAdmin::freeDOFIndex()"); TEST_EXIT_DBG(usedCount > 0)("no dofs in use\n"); @@ -98,19 +98,17 @@ namespace AMDiS { std::list<DOFIndexedBase*>::iterator di; std::list<DOFIndexedBase*>::iterator end = dofIndexedList.end(); - for (di = dofIndexedList.begin(); di != end; ++di) { + for (di = dofIndexedList.begin(); di != end; ++di) (*di)->freeDOFContent(dof); - } std::list<DOFContainer*>::iterator dc; std::list<DOFContainer*>::iterator dcend = dofContainerList.end(); - for (dc = dofContainerList.begin(); dc != dcend; ++dc) { + for (dc = dofContainerList.begin(); dc != dcend; ++dc) (*dc)->freeDOFIndex(dof); - } dofFree[dof] = true; - + if (static_cast<int>(firstHole) > dof) firstHole = dof; @@ -118,8 +116,6 @@ namespace AMDiS { holeCount++; } - /****************************************************************************/ - int DOFAdmin::getDOFIndex() { FUNCNAME("DOFAdmin::getDOFIndex()"); @@ -160,9 +156,6 @@ namespace AMDiS { return(dof); } - - /****************************************************************************/ - void DOFAdmin::enlargeDOFLists(int minsize) { FUNCNAME("DOFAdmin::enlargeDOFLists()"); @@ -256,9 +249,6 @@ namespace AMDiS { ERROR("container not in list\n"); } - - /****************************************************************************/ - void DOFAdmin::compress(std::vector<DegreeOfFreedom> &new_dof) { FUNCNAME("DOFAdmin::compress()"); @@ -269,16 +259,14 @@ namespace AMDiS { if (holeCount < 1) return; // vector to mark used dofs - for (int i = 0; i < size; i++) { + for (int i = 0; i < size; i++) new_dof[i] = -1; - } // mark used dofs DOFIteratorBase it(this, USED_DOFS); - for (it.reset(); !it.end(); ++it) { + for (it.reset(); !it.end(); ++it) new_dof[it.getDOFIndex()] = 1; - } - + int n = 0, last = 0; for (int i = 0; i < size; i++) { /* create a MONOTONE compress */ if (new_dof[i] == 1) { @@ -290,13 +278,12 @@ namespace AMDiS { TEST_EXIT_DBG(n == usedCount)("count %d != usedCount %d\n", n, usedCount); // mark used dofs in compressed dofFree - for (int i = 0; i < n; i++) { + for (int i = 0; i < n; i++) dofFree[i] = false; - } + // mark unused dofs in compressed dofFree - for (int i = n; i < size; i++) { + for (int i = n; i < size; i++) dofFree[i] = true; - } firstHole = n; holeCount = 0; @@ -313,18 +300,14 @@ namespace AMDiS { std::list<DOFIndexedBase*>::iterator di; std::list<DOFIndexedBase*>::iterator end = dofIndexedList.end(); - for (di = dofIndexedList.begin(); di != end; ++di) { + for (di = dofIndexedList.begin(); di != end; ++di) (*di)->compressDOFIndexed(first, last, new_dof); - }; std::list<DOFContainer*>::iterator dc; std::list<DOFContainer*>::iterator endc = dofContainerList.end(); - for (dc = dofContainerList.begin(); dc != endc; dc++) { + for (dc = dofContainerList.begin(); dc != endc; dc++) (*dc)->compressDOFContainer(n, new_dof); - }; - - return; } void DOFAdmin::setNumberOfDOFs(int i,int v) { diff --git a/AMDiS/src/DOFMatrix.cc b/AMDiS/src/DOFMatrix.cc index bd206eb2cb75dc166e96f22d31aa26bd69e03c71..c1f3116bc59d5d3b3f0105727d53c42c799ba699 100644 --- a/AMDiS/src/DOFMatrix.cc +++ b/AMDiS/src/DOFMatrix.cc @@ -98,12 +98,11 @@ namespace AMDiS { typedef traits::range_generator<major, Matrix>::type cursor_type; typedef traits::range_generator<nz, cursor_type>::type icursor_type; - - for (cursor_type cursor = begin<major>(matrix), cend = end<major>(matrix); cursor != cend; ++cursor) { + + std::cout.precision(10); + for (cursor_type cursor = begin<major>(matrix), cend = end<major>(matrix); cursor != cend; ++cursor) for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) - Msg::print(" (%3d,%3d,%20.17lf)", row(*icursor), col(*icursor), value(*icursor)); - Msg::print("\n"); - } + std::cout << row(*icursor) << " " << col(*icursor) << " " << value(*icursor) << "\n"; } bool DOFMatrix::symmetric() @@ -201,7 +200,6 @@ namespace AMDiS { } } - // === Add element matrix to the global matrix using the indices mapping. === DegreeOfFreedom row, col; @@ -219,6 +217,7 @@ namespace AMDiS { for (int j = 0; j < nCol; j++) { // for all columns col = colIndices[j]; entry = elMat[i][j]; + if (add) ins[row][col]+= sign * entry; else @@ -233,33 +232,8 @@ namespace AMDiS { } void DOFMatrix::freeDOFContent(int index) - { - if (matrix.nnz() == 0) return; - - using mtl::tag::major; using mtl::tag::nz; using mtl::begin; using mtl::end; - namespace traits= mtl::traits; - typedef base_matrix_type Matrix; - - traits::row<Matrix>::type row(matrix); - traits::col<Matrix>::type col(matrix); - - typedef traits::range_generator<major, Matrix>::type cursor_type; - typedef traits::range_generator<nz, cursor_type>::type icursor_type; - - cursor_type cursor = begin<major>(matrix); - // Jump directly to corresponding row or column - cursor+= index; - - // Requires structural symmetry !!! - for (icursor_type icursor = begin<nz>(cursor), icend = end<nz>(cursor); icursor != icend; ++icursor) { - int my_row= row(*icursor), my_col= col(*icursor); - // Not very efficient (but general) - matrix.lvalue(my_row, my_col) = 0.0; // Need to call crop somewhere !!!! Peter - matrix.lvalue(my_col, my_row) = 0.0; - } - } + {} - // Should work as before void DOFMatrix::assemble(double factor, ElInfo *elInfo, const BoundaryType *bound) { FUNCNAME("DOFMatrix::assemble()"); @@ -268,13 +242,9 @@ namespace AMDiS { std::vector<Operator*>::iterator it = operators.begin(); std::vector<double*>::iterator factorIt = operatorFactor.begin(); - for (; it != operators.end(); ++it, ++factorIt) { - if ((*it)->getNeedDualTraverse() == false) { - (*it)->getElementMatrix(elInfo, - elementMatrix, - *factorIt ? **factorIt : 1.0); - } - } + for (; it != operators.end(); ++it, ++factorIt) + if ((*it)->getNeedDualTraverse() == false) + (*it)->getElementMatrix(elInfo, elementMatrix, *factorIt ? **factorIt : 1.0); addElementMatrix(factor, elementMatrix, bound, elInfo, NULL); } diff --git a/AMDiS/src/DOFMatrix.h b/AMDiS/src/DOFMatrix.h index 04becace783d19870c03bc3cd52cf3ff4820ff1c..a72f9ebd05dca232ed3e905af77eaf58110cfcbf 100644 --- a/AMDiS/src/DOFMatrix.h +++ b/AMDiS/src/DOFMatrix.h @@ -290,8 +290,6 @@ namespace AMDiS { /// Resizes \ref matrix to n rows inline void resize(int n) { TEST_EXIT_DBG(n >= 0)("Can't resize DOFMatrix to negative size\n"); - // matrix.change_dim(n, n); - //WARNING("Obsolete function without effect -- Peter\n"); } /// Returns value at logical indices a,b diff --git a/AMDiS/src/DOFVector.cc b/AMDiS/src/DOFVector.cc index f90020f5aaf82aa6bcad3b8f2da25056d80e4c88..718a3bd4cdddc61f4e597ecbd41452d6a2c9ad11 100644 --- a/AMDiS/src/DOFVector.cc +++ b/AMDiS/src/DOFVector.cc @@ -65,8 +65,8 @@ namespace AMDiS { result = grad; } else { if(result && result->getFESpace() != feSpace) { - DELETE result; - result = NEW DOFVector<WorldVector<double> >(feSpace, "gradient"); + delete result; + result = new DOFVector<WorldVector<double> >(feSpace, "gradient"); } } @@ -157,11 +157,11 @@ namespace AMDiS { if (!result) { if (vec && vec->getFESpace() != feSpace) { - DELETE vec; + delete vec; vec = NULL; } if (!vec) { - vec = NEW DOFVector<WorldVector<double> >(feSpace, "gradient"); + vec = new DOFVector<WorldVector<double> >(feSpace, "gradient"); } result = vec; } @@ -257,48 +257,42 @@ namespace AMDiS { if (quadFast) { for (int i = 0; i < nPoints; i++) { - for (int j = 0; j < dim + 1; j++) { + for (int j = 0; j < dim + 1; j++) grd1[j] = 0.0; - } - for (int j = 0; j < nBasFcts; j++) { - for (int k = 0; k < parts; k++) { + for (int j = 0; j < nBasFcts; j++) + for (int k = 0; k < parts; k++) grd1[k] += quadFast->getGradient(i, j, k) * localVec[j]; - } - } - + for (int l=0; l < dow; l++) { result[i][l] = 0.0; - for (int k = 0; k < parts; k++) { + for (int k = 0; k < parts; k++) result[i][l] += grdLambda[k][l] * grd1[k]; - } } } - } else { + + } else { const BasisFunction *basFcts = feSpace->getBasisFcts(); DimVec<double>* grdPhi = grdPhis[myRank]; for (int i = 0; i < nPoints; i++) { - for (int j = 0; j < dim + 1; j++) { + for (int j = 0; j < dim + 1; j++) grd1[j] = 0.0; - } for (int j = 0; j < nBasFcts; j++) { (*(basFcts->getGrdPhi(j)))(quad->getLambda(i), *grdPhi); - for (int k = 0; k < parts; k++) { + for (int k = 0; k < parts; k++) grd1[k] += (*grdPhi)[k] * localVec[j]; - } } for (int l = 0; l < dow; l++) { result[i][l] = 0.0; - for (int k = 0; k < parts; k++) { + for (int k = 0; k < parts; k++) result[i][l] += grdLambda[k][l] * grd1[k]; - } } } } - + return const_cast<const WorldVector<double>*>(result); } @@ -585,7 +579,7 @@ namespace AMDiS { if (feSpace->getMesh() == vFESpace->getMesh()) { DegreeOfFreedom *myLocalIndices = localIndices[omp_get_thread_num()]; - WorldVector<double> *vLocalCoeffs = NEW WorldVector<double>[vNumBasFcts]; + WorldVector<double> *vLocalCoeffs = new WorldVector<double>[vNumBasFcts]; Mesh *mesh = feSpace->getMesh(); TraverseStack stack; ElInfo *elInfo = stack.traverseFirst(mesh, -1, @@ -608,7 +602,7 @@ namespace AMDiS { elInfo = stack.traverseNext(elInfo); } - DELETE [] vLocalCoeffs; + delete [] vLocalCoeffs; } else { ERROR_EXIT("not yet for dual traverse\n"); } @@ -635,14 +629,14 @@ namespace AMDiS { result = grad; } else { if (!result) { - result = NEW WorldVector<DOFVector<double>*>; + result = new WorldVector<DOFVector<double>*>; result->set(NULL); } for (int i = 0; i < dow; i++) { if ((*result)[i] && (*result)[i]->getFESpace() != feSpace) { - DELETE (*result)[i]; - (*result)[i] = NEW DOFVector<double>(feSpace, "gradient"); + delete (*result)[i]; + (*result)[i] = new DOFVector<double>(feSpace, "gradient"); } } } @@ -720,17 +714,15 @@ namespace AMDiS { DOFVector<DegreeOfFreedom>() {} - - void DOFVectorDOF::freeDOFContent(DegreeOfFreedom dof) { + void DOFVectorDOF::freeDOFContent(DegreeOfFreedom dof) + { std::vector<DegreeOfFreedom>::iterator it; std::vector<DegreeOfFreedom>::iterator end = vec.end(); DegreeOfFreedom pos = 0; - for (it = vec.begin(); it != end; ++it, ++pos) { + for (it = vec.begin(); it != end; ++it, ++pos) if (*it == dof) *it = pos; - } } - WorldVector<DOFVector<double>*> *transform(DOFVector<WorldVector<double> > *vec, WorldVector<DOFVector<double>*> *res) { @@ -742,20 +734,17 @@ namespace AMDiS { static WorldVector<DOFVector<double>*> *result = NULL; if (!res && !result) { - result = NEW WorldVector<DOFVector<double>*>; - for (int i = 0; i < dow; i++) { - (*result)[i] = NEW DOFVector<double>(vec->getFESpace(), "transform"); - } + result = new WorldVector<DOFVector<double>*>; + for (int i = 0; i < dow; i++) + (*result)[i] = new DOFVector<double>(vec->getFESpace(), "transform"); } WorldVector<DOFVector<double>*> *r = res ? res : result; int vecSize = vec->getSize(); - for (int i = 0; i < vecSize; i++) { - for (int j = 0; j < dow; j++) { + for (int i = 0; i < vecSize; i++) + for (int j = 0; j < dow; j++) (*((*r)[j]))[i] = (*vec)[i][j]; - } - } return r; } @@ -836,12 +825,10 @@ namespace AMDiS { for (it = this->operators.begin(), factorIt = this->operatorFactor.begin(); it != this->operators.end(); - ++it, ++factorIt) { - if ((*it)->getNeedDualTraverse() == false) { + ++it, ++factorIt) + if ((*it)->getNeedDualTraverse() == false) (*it)->getElementVector(elInfo, this->elementVector, - *factorIt ? **factorIt : 1.0); - } - } + *factorIt ? **factorIt : 1.0); } addElementVector(factor, this->elementVector, bound, elInfo); diff --git a/AMDiS/src/DOFVector.hh b/AMDiS/src/DOFVector.hh index bcfb853bb02e06991d41aabf5c5f45ebbd774328..b96362f14fc4d8759e69d33c8baf1133c549476e 100644 --- a/AMDiS/src/DOFVector.hh +++ b/AMDiS/src/DOFVector.hh @@ -86,9 +86,9 @@ namespace AMDiS { for (int i = 0; i < omp_get_overall_max_threads(); i++) { localIndices[i] = GET_MEMORY(DegreeOfFreedom, this->nBasFcts); localVectors[i] = GET_MEMORY(T, this->nBasFcts); - grdPhis[i] = NEW DimVec<double>(dim, DEFAULT_VALUE, 0.0); - grdTmp[i] = NEW DimVec<double>(dim, DEFAULT_VALUE, 0.0); - D2Phis[i] = NEW DimMat<double>(dim, NO_INIT); + grdPhis[i] = new DimVec<double>(dim, DEFAULT_VALUE, 0.0); + grdTmp[i] = new DimVec<double>(dim, DEFAULT_VALUE, 0.0); + D2Phis[i] = new DimMat<double>(dim, NO_INIT); } } @@ -98,9 +98,9 @@ namespace AMDiS { for (int i = 0; i < static_cast<int>(localIndices.size()); i++) { FREE_MEMORY(localIndices[i], DegreeOfFreedom, this->nBasFcts); FREE_MEMORY(localVectors[i], T, this->nBasFcts); - DELETE grdPhis[i]; - DELETE grdTmp[i]; - DELETE D2Phis[i]; + delete grdPhis[i]; + delete grdTmp[i]; + delete D2Phis[i]; } } @@ -124,7 +124,7 @@ namespace AMDiS { if(feSpace && feSpace->getAdmin()) { (feSpace->getAdmin())->addDOFIndexed(this); } - this->boundaryManager = NEW BoundaryManager(f); + this->boundaryManager = new BoundaryManager(f); } template<typename T> @@ -134,7 +134,7 @@ namespace AMDiS { (feSpace->getAdmin())->removeDOFIndexed(this); if (this->boundaryManager) - DELETE this->boundaryManager; + delete this->boundaryManager; vec.clear(); } @@ -707,9 +707,8 @@ namespace AMDiS { if (rhs.boundaryManager) { if (this->boundaryManager) delete this->boundaryManager; - this->boundaryManager = new BoundaryManager(*rhs.boundaryManager); - // boundaryManager->setDOFVector(this); + this->boundaryManager = new BoundaryManager(*rhs.boundaryManager); } else { this->boundaryManager = NULL; } @@ -726,9 +725,8 @@ namespace AMDiS { ("pointer is NULL: %8X, %8X\n", x.getFESpace(), x.getFESpace()->getAdmin()); typename DOFVector<T>::Iterator vecIterator(dynamic_cast<DOFIndexed<T>*>(&x), USED_DOFS); - for (vecIterator.reset(); !vecIterator.end(); ++vecIterator) { + for (vecIterator.reset(); !vecIterator.end(); ++vecIterator) (*vecIterator) *= scal; - } return x; } diff --git a/AMDiS/src/FirstOrderAssembler.cc b/AMDiS/src/FirstOrderAssembler.cc index dcb19c3c27610e3be6a6bb6c72ff5190709f9552..2493b16c5dd3bdd134a31cfc6a184166d62a6f6a 100644 --- a/AMDiS/src/FirstOrderAssembler.cc +++ b/AMDiS/src/FirstOrderAssembler.cc @@ -59,7 +59,7 @@ namespace AMDiS { FirstOrderAssembler *newAssembler; // check if a new assembler is needed - for (int i = 0; i < static_cast<int>( subAssemblers->size()); i++) { + for (int i = 0; i < static_cast<int>(subAssemblers->size()); i++) { std::vector<OperatorTerm*> assTerms = *((*subAssemblers)[i]->getTerms()); sort(assTerms.begin(), assTerms.end()); @@ -70,7 +70,7 @@ namespace AMDiS { // check if all terms are pw_const bool pwConst = true; - for (int i = 0; i < static_cast<int>( opTerms.size()); i++) { + for (int i = 0; i < static_cast<int>(opTerms.size()); i++) { if (!(opTerms[i])->isPWConst()) { pwConst = false; break; diff --git a/AMDiS/src/Lagrange.cc b/AMDiS/src/Lagrange.cc index 642e420404cd7d1089fecc46add02508c1d42d89..ead066766c94c1d0d3d9c2f610efcda76eed8fe7 100644 --- a/AMDiS/src/Lagrange.cc +++ b/AMDiS/src/Lagrange.cc @@ -51,7 +51,7 @@ namespace AMDiS { { for (int i = 0; i < static_cast<int>(bary->size()); i++) { if ((*bary)[i]) { - DELETE (*bary)[i]; + delete (*bary)[i]; (*bary)[i] = NULL; } } @@ -66,7 +66,7 @@ namespace AMDiS { } } - Lagrange* newLagrange = NEW Lagrange(dim, degree); + Lagrange* newLagrange = new Lagrange(dim, degree); allBasFcts.push_back(newLagrange); return newLagrange; } @@ -76,7 +76,7 @@ namespace AMDiS { std::list<Lagrange*>::iterator it; for (it = allBasFcts.begin(); it != allBasFcts.end(); it++) { if (*it) { - DELETE *it; + delete *it; *it = NULL; } } @@ -95,14 +95,14 @@ namespace AMDiS { // for all dofs at this position for (int k = 0; k < (*nDOF)[INDEX_OF_DIM(i, dim)]; k++) { // basis function - phiDimDegree[dim][degree].push_back(NEW Phi(this, + phiDimDegree[dim][degree].push_back(new Phi(this, INDEX_OF_DIM(i ,dim),j, k)); // gradients - grdPhiDimDegree[dim][degree].push_back(NEW GrdPhi(this, + grdPhiDimDegree[dim][degree].push_back(new GrdPhi(this, INDEX_OF_DIM(i, dim), j, k)); // D2-Matrices - D2PhiDimDegree[dim][degree].push_back(NEW D2Phi(this, + D2PhiDimDegree[dim][degree].push_back(new D2Phi(this, INDEX_OF_DIM(i, dim), j, k)); } @@ -192,7 +192,7 @@ namespace AMDiS { void Lagrange::setNDOF() { if (static_cast<int>(baryDimDegree[dim][degree].size()) == 0) { - ndofDimDegree[dim][degree] = NEW DimVec<int>(dim, DEFAULT_VALUE, 0); + ndofDimDegree[dim][degree] = new DimVec<int>(dim, DEFAULT_VALUE, 0); if (degree != 0) { (*ndofDimDegree[dim][degree])[VERTEX] = 1; @@ -380,7 +380,7 @@ namespace AMDiS { } Lagrange::Phi::~Phi() { - DELETE [] vertices; + delete [] vertices; } @@ -507,7 +507,7 @@ namespace AMDiS { } Lagrange::GrdPhi::~GrdPhi() { - DELETE [] vertices; + delete [] vertices; } Lagrange::D2Phi::D2Phi(Lagrange* owner, @@ -634,7 +634,7 @@ namespace AMDiS { } Lagrange::D2Phi::~D2Phi() { - DELETE [] vertices; + delete [] vertices; } void Lagrange::createCoords(int* coordInd, @@ -644,12 +644,12 @@ namespace AMDiS { DimVec<double>* vec) { if (vec == NULL) { - vec = NEW DimVec<double>(dim, DEFAULT_VALUE, 0.0); + vec = new DimVec<double>(dim, DEFAULT_VALUE, 0.0); } if (dimIndex == numCoords - 1) { (*vec)[coordInd[dimIndex]] = double(rest) / degree; - DimVec<double>* newCoords = NEW DimVec<double>(*vec); + DimVec<double>* newCoords = new DimVec<double>(*vec); bary->push_back(newCoords); } else { for (int i = rest - 1; i >= 1; i--) { @@ -876,7 +876,7 @@ namespace AMDiS { static WorldVector<double> *inter; WorldVector<double> *rvec = - vec ? vec : (inter?inter:inter= NEW WorldVector<double>[getNumber()]); + vec ? vec : (inter ? inter : inter = new WorldVector<double>[getNumber()]); WorldVector<double> x; el_info->testFlag(Mesh::FILL_COORDS); @@ -1031,8 +1031,6 @@ namespace AMDiS { } FREE_MEMORY(wdetf_qp, double, q->getNumPoints()); - - return; } @@ -1056,7 +1054,8 @@ namespace AMDiS { (*drv)[dof_new] = (*drv)[dof0]; } - void Lagrange::refineInter1(DOFIndexed<double> *drv, RCNeighbourList* list, int n, BasisFunction* basFct) + void Lagrange::refineInter1(DOFIndexed<double> *drv, RCNeighbourList* list, int n, + BasisFunction* basFct) { FUNCNAME("Lagrange::refineInter1_1d()"); @@ -1072,7 +1071,8 @@ namespace AMDiS { (*drv)[dof_new] = 0.5 * ((*drv)[dof0] + (*drv)[dof1]); } - void Lagrange::refineInter2_1d(DOFIndexed<double> *drv, RCNeighbourList* list, int n, BasisFunction* basFct) + void Lagrange::refineInter2_1d(DOFIndexed<double> *drv, RCNeighbourList* list, int n, + BasisFunction* basFct) { FUNCNAME("Lagrange::refineInter2_1d()"); @@ -1111,7 +1111,6 @@ namespace AMDiS { cdof = el->getChild(1)->getDOF(node, n0); (*drv)[cdof] = -0.125 * (*drv)[pdof[0]] + 0.375 * (*drv)[pdof[1]] + 0.75 * (*drv)[pdof[2]]; - return; } void Lagrange::refineInter2_2d(DOFIndexed<double> *drv, @@ -1259,7 +1258,6 @@ namespace AMDiS { } } - return; } void Lagrange::refineInter3_1d(DOFIndexed<double> *drv, @@ -1355,7 +1353,6 @@ namespace AMDiS { + 0.75 *(*drv)[pdof[9]]); FREE_MEMORY(pdof, DegreeOfFreedom, basFct->getNumber()); - return; } void Lagrange::refineInter3_2d(DOFIndexed<double> *drv, @@ -1449,8 +1446,6 @@ namespace AMDiS { (0.0625 * ((*drv)[pdof[0]] - (*drv)[pdof[1]]) + 0.375 * (*drv)[pdof[3]] - 0.125 * (*drv)[pdof[6]] + 0.1875 * (-(*drv)[pdof[7]] + (*drv)[pdof[8]]) + 0.75 * (*drv)[pdof[9]]); - - return; } void Lagrange::refineInter3_3d(DOFIndexed<double> *drv, @@ -1461,10 +1456,12 @@ namespace AMDiS { Element *el; const DegreeOfFreedom *cd; DegreeOfFreedom pd[20], cdi; - int i, typ, lr_set, node0, n0; + int i, typ, lr_set, node0, n0; const DOFAdmin *admin; - if (n < 1) return; + if (n < 1) + return; + el = list->getElement(0); typ = list->getType(0); @@ -1673,8 +1670,6 @@ namespace AMDiS { } } } - - return; } void Lagrange::refineInter4_1d(DOFIndexed<double> *drv, @@ -1683,9 +1678,10 @@ namespace AMDiS { { FUNCNAME("Lagrange::refineInter4_1d"); - if (n < 1) return; + if (n < 1) + return; - Element *el; + Element *el; DegreeOfFreedom *pdof = GET_MEMORY(DegreeOfFreedom, basFct->getNumber()); const DegreeOfFreedom *cdof; const DOFAdmin *admin; @@ -1830,7 +1826,6 @@ namespace AMDiS { (*drv)[cdof[14]] = (*drv)[pdof[13]]; FREE_MEMORY(pdof, DegreeOfFreedom, basFct->getNumber()); - return; } void Lagrange::refineInter4_2d(DOFIndexed<double> *drv, @@ -1843,11 +1838,11 @@ namespace AMDiS { const DegreeOfFreedom *cdof; const DOFAdmin *admin; - if (n < 1) return; - el = list->getElement(0); + if (n < 1) + return; + el = list->getElement(0); admin = drv->getFESpace()->getAdmin(); - basFct->getLocalIndices(el, admin, pdof); /****************************************************************************/ @@ -1922,7 +1917,8 @@ namespace AMDiS { - 0.03125*(*drv)[pdof[11]] + 0.75*(*drv)[pdof[14]]); (*drv)[cdof[14]] = (*drv)[pdof[13]]; - if (n <= 1) return; + if (n <= 1) + return; /****************************************************************************/ /* adjust neighbour values */ @@ -1982,8 +1978,6 @@ namespace AMDiS { + 0.09375*(*drv)[pdof[9]] - 0.046875*(*drv)[pdof[10]] - 0.03125*(*drv)[pdof[11]] + 0.75*(*drv)[pdof[14]]); (*drv)[cdof[14]] = (*drv)[pdof[13]]; - - return; } void Lagrange::refineInter4_3d(DOFIndexed<double> *drv, @@ -1993,15 +1987,14 @@ namespace AMDiS { FUNCNAME("Lagrange::refineInter4_3d"); DegreeOfFreedom pd[35]; const DegreeOfFreedom *cd; - Element *el; - int i, typ, lr_set; - const DOFAdmin *admin; + int i, lr_set; - if (n < 1) return; - el = list->getElement(0); - typ = list->getType(0); + if (n < 1) + return; - admin = drv->getFESpace()->getAdmin(); + Element *el = list->getElement(0); + int typ = list->getType(0); + const DOFAdmin *admin = drv->getFESpace()->getAdmin(); basFct->getLocalIndices(el, admin, pd); @@ -2545,28 +2538,23 @@ namespace AMDiS { } } } - - return; } - // ============ coarseRestrict fcts ================== void Lagrange::coarseRestr0(DOFIndexed<double> *drv, RCNeighbourList *list, int n, BasisFunction* basFct) { - FUNCNAME("Lagrange::coarseRestr0"); - Element *el; - DegreeOfFreedom dof_new, dof0, dof1; - int n0; + FUNCNAME("Lagrange::coarseRestr0()"); - if (n < 1) return; + if (n < 1) + return; - n0 = drv->getFESpace()->getAdmin()->getNumberOfPreDOFs(CENTER); - el = list->getElement(0); - dof0 = el->getDOF(0,n0); /* 1st endpoint of refinement edge */ - dof1 = el->getDOF(1,n0); /* 2nd endpoint of refinement edge */ - dof_new = el->getChild(0)->getDOF(basFct->getDim(), n0); + int n0 = drv->getFESpace()->getAdmin()->getNumberOfPreDOFs(CENTER); + Element* el = list->getElement(0); + int dof0 = el->getDOF(0,n0); /* 1st endpoint of refinement edge */ + int dof1 = el->getDOF(1,n0); /* 2nd endpoint of refinement edge */ + int dof_new = el->getChild(0)->getDOF(basFct->getDim(), n0); /* newest vertex is DIM */ (*drv)[dof0] += 0.5*(*drv)[dof_new]; (*drv)[dof1] += 0.5*(*drv)[dof_new]; @@ -2576,64 +2564,56 @@ namespace AMDiS { RCNeighbourList *list, int n, BasisFunction* basFct) { - FUNCNAME("Lagrange::coarseRestr1"); - Element *el; - DegreeOfFreedom dof_new, dof0, dof1; - int n0; + FUNCNAME("Lagrange::coarseRestr1()"); - if (n < 1) return; + if (n < 1) + return; - n0 = drv->getFESpace()->getAdmin()->getNumberOfPreDOFs(VERTEX); - el = list->getElement(0); - dof0 = el->getDOF(0,n0); /* 1st endpoint of refinement edge */ - dof1 = el->getDOF(1,n0); /* 2nd endpoint of refinement edge */ - dof_new = el->getChild(0)->getDOF(basFct->getDim(), n0); - /* newest vertex is DIM */ - (*drv)[dof0] += 0.5*(*drv)[dof_new]; - (*drv)[dof1] += 0.5*(*drv)[dof_new]; + int n0 = drv->getFESpace()->getAdmin()->getNumberOfPreDOFs(VERTEX); + Element *el = list->getElement(0); + + // 1st endpoint of refinement edge + DegreeOfFreedom dof0 = el->getDOF(0, n0); + // 2nd endpoint of refinement edge + DegreeOfFreedom dof1 = el->getDOF(1, n0); + DegreeOfFreedom dof_new = el->getChild(0)->getDOF(basFct->getDim(), n0); + // newest vertex is DIM + (*drv)[dof0] += 0.5 * (*drv)[dof_new]; + (*drv)[dof1] += 0.5 * (*drv)[dof_new]; } void Lagrange::coarseRestr2_2d(DOFIndexed<double> *drv, RCNeighbourList *list, int n, BasisFunction* basFct) { - FUNCNAME("Lagrange::coarseRestr2_2d"); - Element *el; - int node, n0; - DegreeOfFreedom cdof2, cdof3, cdof4; - const DegreeOfFreedom *pdof; - const DOFAdmin *admin; + FUNCNAME("Lagrange::coarseRestr2_2d()"); - if (n < 1) return; - el = list->getElement(0); - - if (!drv->getFESpace()) - { - ERROR("no fe_space in dof_real_vec\n"); - return; - } - else if (!drv->getFESpace()->getBasisFcts()) - { - ERROR("no basis functions in fe_space %s\n", - drv->getFESpace()->getName().c_str()); - return; - } + if (n < 1) + return; - admin = drv->getFESpace()->getAdmin(); + if (!drv->getFESpace()) { + ERROR("no fe_space in dof_real_vec\n"); + return; + } else if (!drv->getFESpace()->getBasisFcts()) { + ERROR("no basis functions in fe_space %s\n", drv->getFESpace()->getName().c_str()); + return; + } - pdof = basFct->getLocalIndices(el, admin, NULL); + Element *el = list->getElement(0); + const DOFAdmin *admin = drv->getFESpace()->getAdmin(); + const DegreeOfFreedom *pdof = basFct->getLocalIndices(el, admin, NULL); /****************************************************************************/ /* contributions of dofs located on child[0] */ /****************************************************************************/ - node = drv->getFESpace()->getMesh()->getNode(VERTEX); - n0 = admin->getNumberOfPreDOFs(VERTEX); - cdof2 = el->getChild(0)->getDOF(node+2, n0); + int node = drv->getFESpace()->getMesh()->getNode(VERTEX); + int n0 = admin->getNumberOfPreDOFs(VERTEX); + DegreeOfFreedom cdof2 = el->getChild(0)->getDOF(node+2, n0); node = drv->getFESpace()->getMesh()->getNode(EDGE); n0 = admin->getNumberOfPreDOFs(EDGE); - cdof3 = el->getChild(0)->getDOF(node, n0); - cdof4 = el->getChild(0)->getDOF(node+1, n0); + DegreeOfFreedom cdof3 = el->getChild(0)->getDOF(node, n0); + DegreeOfFreedom cdof4 = el->getChild(0)->getDOF(node+1, n0); (*drv)[pdof[0]] += 0.375*(*drv)[cdof3] - 0.125*(*drv)[cdof4]; (*drv)[pdof[1]] += -0.125*((*drv)[cdof3] + (*drv)[cdof4]); @@ -2651,28 +2631,26 @@ namespace AMDiS { (*drv)[pdof[1]] += 0.375*(*drv)[cdof4]; (*drv)[pdof[5]] += 0.75*(*drv)[cdof4]; - if (n > 1) - { - el = list->getElement(1); - pdof = basFct->getLocalIndices(el, admin, NULL); - - /****************************************************************************/ - /* first set those values not effected by previous element */ - /****************************************************************************/ - - cdof4 = el->getChild(0)->getDOF(node+1, n0); - (*drv)[pdof[3]] += 0.5*(*drv)[cdof4]; - (*drv)[pdof[4]] += 0.5*(*drv)[cdof4]; - - /****************************************************************************/ - /* and now, update the values in the refinement edge */ - /****************************************************************************/ - - (*drv)[pdof[0]] += -0.125*(*drv)[cdof4]; - (*drv)[pdof[1]] += -0.125*(*drv)[cdof4]; - (*drv)[pdof[5]] += 0.25*(*drv)[cdof4]; - } - return; + if (n > 1) { + el = list->getElement(1); + pdof = basFct->getLocalIndices(el, admin, NULL); + + /****************************************************************************/ + /* first set those values not effected by previous element */ + /****************************************************************************/ + + cdof4 = el->getChild(0)->getDOF(node+1, n0); + (*drv)[pdof[3]] += 0.5*(*drv)[cdof4]; + (*drv)[pdof[4]] += 0.5*(*drv)[cdof4]; + + /****************************************************************************/ + /* and now, update the values in the refinement edge */ + /****************************************************************************/ + + (*drv)[pdof[0]] += -0.125*(*drv)[cdof4]; + (*drv)[pdof[1]] += -0.125*(*drv)[cdof4]; + (*drv)[pdof[5]] += 0.25*(*drv)[cdof4]; + } } @@ -2690,17 +2668,14 @@ namespace AMDiS { if (n < 1) return; el = list->getElement(0); - if (!drv->getFESpace()) - { - ERROR("no fe_space in dof_real_vec\n"); - return; - } - else if (!drv->getFESpace()->getBasisFcts()) - { - ERROR("no basis functions in fe_space %s\n", - drv->getFESpace()->getName().c_str()); - return; - } + if (!drv->getFESpace()) { + ERROR("no fe_space in dof_real_vec\n"); + return; + } else if (!drv->getFESpace()->getBasisFcts()) { + ERROR("no basis functions in fe_space %s\n", + drv->getFESpace()->getName().c_str()); + return; + } admin = drv->getFESpace()->getAdmin(); @@ -2786,7 +2761,6 @@ namespace AMDiS { break; } } - return; } @@ -2803,17 +2777,14 @@ namespace AMDiS { if (n < 1) return; el = list->getElement(0); - if (!drv->getFESpace()) - { - ERROR("no fe_space in dof_real_vec\n"); - return; - } - else if (!drv->getFESpace()->getBasisFcts()) - { - ERROR("no basis functions in fe_space %s\n", - drv->getFESpace()->getName().c_str()); - return; - } + if (!drv->getFESpace()) { + ERROR("no fe_space in dof_real_vec\n"); + return; + } else if (!drv->getFESpace()->getBasisFcts()) { + ERROR("no basis functions in fe_space %s\n", + drv->getFESpace()->getName().c_str()); + return; + } admin = drv->getFESpace()->getAdmin(); @@ -2902,8 +2873,6 @@ namespace AMDiS { (*drv)[pdof[7]] -= 0.1875*(*drv)[dof9]; (*drv)[pdof[8]] += 0.1875*(*drv)[dof9]; (*drv)[pdof[9]] += 0.75*(*drv)[dof9]; - - return; } void Lagrange::coarseRestr3_2d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -3013,8 +2982,6 @@ namespace AMDiS { (*drv)[pdof[7]] -= 0.1875*(*drv)[dof9]; (*drv)[pdof[8]] += 0.1875*(*drv)[dof9]; (*drv)[pdof[9]] += 0.75*(*drv)[dof9]; - - return; } void Lagrange::coarseRestr3_3d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -3290,7 +3257,6 @@ namespace AMDiS { } } } - return; } void Lagrange::coarseRestr4_1d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -3464,8 +3430,6 @@ namespace AMDiS { (*drv)[pdof[13]] += (*drv)[cdof[14]] + 0.9375*(*drv)[cdof[12]] + 0.375*(*drv)[cdof[13]]; (*drv)[pdof[14]] += 0.75*(*drv)[cdof[13]]; - - return; } void Lagrange::coarseRestr4_2d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -3639,8 +3603,6 @@ namespace AMDiS { (*drv)[pdof[13]] += (*drv)[cdof[14]] + 0.9375*(*drv)[cdof[12]] + 0.375*(*drv)[cdof[13]]; (*drv)[pdof[14]] += 0.75*(*drv)[cdof[13]]; - - return; } void Lagrange::coarseRestr4_3d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -4334,8 +4296,6 @@ namespace AMDiS { } } } - - return; } // ====== coarseInter functions ====== @@ -4344,25 +4304,22 @@ namespace AMDiS { int n, BasisFunction* basFct) { FUNCNAME("Lagrange::coarseInter0"); - Element *el; DegreeOfFreedom cdof, pdof; const DOFAdmin *admin; const Mesh *mesh = NULL; if (n < 1) return; - el = list->getElement(0); - if (!drv->getFESpace()) - { + Element* el = list->getElement(0); + + if (!drv->getFESpace()) { ERROR("no fe_space in dof_real_vec\n"); return; - } - else if (!drv->getFESpace()->getBasisFcts()) - { + } else if (!drv->getFESpace()->getBasisFcts()) { ERROR("no basis functions in fe_space %s\n", drv->getFESpace()->getName().c_str()); return; - } + } admin = drv->getFESpace()->getAdmin(); mesh = drv->getFESpace()->getMesh(); @@ -4375,8 +4332,6 @@ namespace AMDiS { admin->getNumberOfPreDOFs(CENTER)); pdof = el->getDOF(mesh->getNode(CENTER)+2, admin->getNumberOfPreDOFs(CENTER)); (*drv)[pdof] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter2_1d(DOFIndexed<double> *drv, RCNeighbourList *list, @@ -4415,8 +4370,6 @@ namespace AMDiS { admin->getNumberOfPreDOFs(VERTEX)); pdof = el->getDOF(mesh->getNode(CENTER), admin->getNumberOfPreDOFs(CENTER)); (*drv)[pdof] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter2_2d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4454,8 +4407,6 @@ namespace AMDiS { admin->getNumberOfPreDOFs(VERTEX)); pdof = el->getDOF(mesh->getNode(EDGE)+2, admin->getNumberOfPreDOFs(EDGE)); (*drv)[pdof] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter2_3d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4488,8 +4439,6 @@ namespace AMDiS { admin->getNumberOfPreDOFs(VERTEX)); pdof = el->getDOF(mesh->getNode(EDGE), admin->getNumberOfPreDOFs(EDGE)); (*drv)[pdof] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter3_1d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4581,8 +4530,6 @@ namespace AMDiS { (*drv)[el->getDOF(mesh->getNode(CENTER), admin->getNumberOfPreDOFs(CENTER))] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter3_2d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4674,8 +4621,6 @@ namespace AMDiS { (*drv)[el->getDOF(mesh->getNode(CENTER), admin->getNumberOfPreDOFs(CENTER))] = (*drv)[cdof]; - - return; } void Lagrange::coarseInter3_3d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4776,7 +4721,6 @@ namespace AMDiS { break; } } - return; } void Lagrange::coarseInter4_1d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4853,8 +4797,6 @@ namespace AMDiS { cdof = basFct->getLocalIndices(el->getChild(1), admin, NULL); (*drv)[pdof[13]] = (*drv)[cdof[14]]; - - return; } void Lagrange::coarseInter4_2d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -4931,8 +4873,6 @@ namespace AMDiS { cdof = basFct->getLocalIndices(el->getChild(1), admin, NULL); (*drv)[pdof[13]] = (*drv)[cdof[14]]; - - return; } void Lagrange::coarseInter4_3d(DOFIndexed<double> *drv, RCNeighbourList* list, @@ -5069,8 +5009,6 @@ namespace AMDiS { } } } - - return; } } diff --git a/AMDiS/src/Lagrange.h b/AMDiS/src/Lagrange.h index 1887cedfee4e111957b5b168d5fe07964026ab2c..97621586cc266ce55d8a08960187de1531a8fe4c 100644 --- a/AMDiS/src/Lagrange.h +++ b/AMDiS/src/Lagrange.h @@ -108,61 +108,43 @@ namespace AMDiS { GeoIndex position, int positionIndex, int nodeIndex, int** vertices); - /** \brief - * Implements BasisFunction::refineInter - */ - inline void refineInter(DOFIndexed<double> *drv, - RCNeighbourList* list, - int n) + /// Implements BasisFunction::refineInter + inline void refineInter(DOFIndexed<double> *drv, RCNeighbourList* list, int n) { if (refineInter_fct) (*refineInter_fct)(drv, list, n, this); } - /** \brief - * Implements BasisFunction::coarseRestrict - */ - inline void coarseRestr(DOFIndexed<double> *drv, - RCNeighbourList* list, - int n) + /// Implements BasisFunction::coarseRestrict + inline void coarseRestr(DOFIndexed<double> *drv, RCNeighbourList* list, int n) { if (coarseRestr_fct) (*coarseRestr_fct)(drv, list, n, this); } - /** \brief - * Implements BasisFunction::coarseInter - */ - inline void coarseInter(DOFIndexed<double> *drv, - RCNeighbourList* list, - int n) + /// Implements BasisFunction::coarseInter + inline void coarseInter(DOFIndexed<double> *drv, RCNeighbourList* list, int n) { if (coarseInter_fct) (*coarseInter_fct)(drv, list, n, this); } - /** \brief - * Implements BasisFunction::getLocalIndices(). - */ + /// Implements BasisFunction::getLocalIndices(). const DegreeOfFreedom *getLocalIndices(const Element*, const DOFAdmin*, DegreeOfFreedom*) const; - + /// void getLocalIndicesVec(const Element*, const DOFAdmin*, Vector<DegreeOfFreedom>*) const; - /** \brief - * Implements BasisFunction::l2ScpFctBas - */ + /// Implements BasisFunction::l2ScpFctBas void l2ScpFctBas(Quadrature* q, AbstractFunction<double, WorldVector<double> >* f, DOFVector<double>* fh); - /** \brief - * Implements BasisFunction::l2ScpFctBas - */ + /// Implements BasisFunction::l2ScpFctBas void l2ScpFctBas(Quadrature* q, AbstractFunction<WorldVector<double>, WorldVector<double> >* f, DOFVector<WorldVector<double> >* fh); diff --git a/AMDiS/src/Operator.cc b/AMDiS/src/Operator.cc index 910ea062615151ba89b0db10a436180e3dc38e84..1d1e7752669c4c1b4c133beba3cc64c9a5c983e6 100644 --- a/AMDiS/src/Operator.cc +++ b/AMDiS/src/Operator.cc @@ -1210,8 +1210,9 @@ namespace AMDiS { void MatrixGradient_SOT::getLALt(const ElInfo *elInfo, int nPoints, DimMat<double> **LALt) const { const DimVec<WorldVector<double> > &Lambda = elInfo->getGrdLambda(); - for (int iq = 0; iq < nPoints; iq++) + for (int iq = 0; iq < nPoints; iq++) { lalt(Lambda, (*f)(gradAtQPs[iq]), (*LALt[iq]), symmetric, 1.0); + } } void FctGradient_SOT::getLALt(const ElInfo *elInfo, int nPoints, DimMat<double> **LALt) const { @@ -1453,12 +1454,11 @@ namespace AMDiS { double *result, double fac) const { - for (int iq = 0; iq < nPoints; iq++) { + for (int iq = 0; iq < nPoints; iq++) result[iq] += fac * (*f)(vecAtQPs[iq], gradAtQPs[iq], coordsAtQPs[iq]) * uhAtQP[iq]; - } } void Vec2AndGradAtQP_ZOT::eval(int nPoints, @@ -1468,27 +1468,25 @@ namespace AMDiS { double *result, double fac) const { - for (int iq = 0; iq < nPoints; iq++) { + for (int iq = 0; iq < nPoints; iq++) result[iq] += fac * (*f)(vecAtQPs1[iq], gradAtQPs[iq], vecAtQPs2[iq]) * uhAtQP[iq]; - } } void Vec2AndGradAtQP_ZOT::getC(const ElInfo *, int nPoints, std::vector<double> &C) const { - for (int iq = 0; iq < nPoints; iq++) { + for (int iq = 0; iq < nPoints; iq++) C[iq] += (*f)(vecAtQPs1[iq], gradAtQPs[iq], vecAtQPs2[iq]); - } } void VecAndGradAtQP_ZOT::getC(const ElInfo *, int nPoints, - std::vector<double> &C) const { - for (int iq = 0; iq < nPoints; iq++) { + std::vector<double> &C) const + { + for (int iq = 0; iq < nPoints; iq++) C[iq] += (*f)(vecAtQPs[iq], gradAtQPs[iq]); - } } void VecAndGradAtQP_ZOT::eval(int nPoints, @@ -1499,8 +1497,7 @@ namespace AMDiS { double fac) const { for (int iq = 0; iq < nPoints; iq++) - result[iq] += - fac * (*f)(vecAtQPs[iq], gradAtQPs[iq]) * uhAtQP[iq]; + result[iq] += fac * (*f)(vecAtQPs[iq], gradAtQPs[iq]) * uhAtQP[iq]; } void VecAndGradVecAtQP_ZOT::getC(const ElInfo *, int nPoints, diff --git a/AMDiS/src/Parameters.cc b/AMDiS/src/Parameters.cc index 359b2527c5e618f98176eafaa8881c5e91d4f830..3c654269df9ede43b55f0e4e18503c19d0b6e16b 100644 --- a/AMDiS/src/Parameters.cc +++ b/AMDiS/src/Parameters.cc @@ -159,7 +159,7 @@ namespace AMDiS { va_end(arg); - return(count); + return count; } int Parameters::getGlobalParameter(int flag, const std::string& key, std::string* param) @@ -211,9 +211,6 @@ namespace AMDiS { } singlett->inputFile.close(); - - - return; } const std::string& Parameters::getActFile(const std::string& aFilename) @@ -221,8 +218,6 @@ namespace AMDiS { FUNCNAME("Parameters::getActFile()"); std::list< std::string >::iterator i; - std::string actfile; - for (i = filenames.begin(); i != filenames.end(); i++) if (aFilename==*i) break; @@ -240,25 +235,21 @@ namespace AMDiS { param tmp(allParam[i]); allParam[i] = allParam[j]; allParam[j] = tmp; - return; } void Parameters::qsort(int left, int right) { - int i, last; - - if (left >= right) return; + if (left >= right) return; - swap(left, (left+right)/2); - last = left; - for (i = left+1; i <= right; i++) + swap(left, (left + right) / 2); + int last = left; + for (int i = left + 1; i <= right; i++) if (allParam[i].key.compare(allParam[left].key) < 0) swap(++last, i); swap(left, last); - qsort(left, last-1); - qsort(last+1, right); - return; + qsort(left, last - 1); + qsort(last + 1, right); } @@ -407,7 +398,7 @@ namespace AMDiS { return(""); } - return(parameter); + return parameter; } void Parameters::addParam(const std::string& key, @@ -449,17 +440,14 @@ namespace AMDiS { newPar.funcName.assign(fname); newPar.lineNo = nLine; allParam.insert(it, newPar); - - return; } void Parameters::initIntern() { - if (singlett == NULL) { - singlett = NEW Parameters; - } + if (singlett == NULL) + singlett = new Parameters; } void Parameters::printParameters() @@ -468,8 +456,6 @@ namespace AMDiS { initIntern(); - std::cout << "SIZE: " << singlett->allParam.size() << std::endl; - std::vector<param>::iterator it; for (it = singlett->allParam.begin(); it != singlett->allParam.end(); it++) { MSG("%s: %s\n", (*it).key.data(), (*it).parameters.data()); @@ -630,18 +616,15 @@ namespace AMDiS { printParameters(); Global::init(); - return; } void Parameters::readArgv(int argc, char **argv) { FUNCNAME("Parameters::readArgv()"); - for (int i = 0; i < argc; i++) { - if (strcmp("-rs", argv[i]) == 0) { + for (int i = 0; i < argc; i++) + if (strcmp("-rs", argv[i]) == 0) ADD_PARAMETER(0, "argv->rs", argv[i + 1]); - } - } } int Parameters::initFuncName(const char * call_fct, const char * call_file, @@ -651,7 +634,7 @@ namespace AMDiS { param_call_file = call_file; param_call_line = call_line; - return(1); + return 1; } int Parameters::binSearch(const std::string& key, int n_keys) @@ -765,9 +748,8 @@ namespace AMDiS { void Parameters::clear() { - if (singlett != NULL) { - DELETE singlett; - } + if (singlett != NULL) + delete singlett; } int Parameters::param::operator==(const param& aParam)const{ diff --git a/AMDiS/src/ProblemScal.cc b/AMDiS/src/ProblemScal.cc index 5306b30f58d1bbbcafee931b750de6db474a2cca..13c0a9ff577b7d9e3c070ca9d56b181041e655cc 100644 --- a/AMDiS/src/ProblemScal.cc +++ b/AMDiS/src/ProblemScal.cc @@ -30,33 +30,32 @@ namespace AMDiS { { FUNCNAME("ProblemScal::~ProblemScal()"); - for (int i = 0; i < static_cast<int>(fileWriters.size()); i++) { - DELETE fileWriters[i]; - } + for (int i = 0; i < static_cast<int>(fileWriters.size()); i++) + delete fileWriters[i]; if (systemMatrix) - DELETE systemMatrix; + delete systemMatrix; if (rhs) - DELETE rhs; + delete rhs; if (solution) - DELETE solution; + delete solution; if (estimator) - DELETE estimator; + delete estimator; if (marker) - DELETE marker; + delete marker; if (solver) - DELETE solver; + delete solver; if (mesh) - DELETE mesh; + delete mesh; FiniteElemSpace::clear(); Lagrange::clear(); } - void ProblemScal::writeFiles(AdaptInfo *adaptInfo, bool force) { - for (int i = 0; i < static_cast<int>(fileWriters.size()); i++) { + void ProblemScal::writeFiles(AdaptInfo *adaptInfo, bool force) + { + for (int i = 0; i < static_cast<int>(fileWriters.size()); i++) fileWriters[i]->writeFiles(adaptInfo, force); - } } void ProblemScal::interpolInitialSolution(AbstractFunction<double, WorldVector<double> > *fct) @@ -170,20 +169,20 @@ namespace AMDiS { TEST_EXIT(dim)("no problem dimension specified!\n"); // create the mesh - mesh = NEW Mesh(meshName, dim); + mesh = new Mesh(meshName, dim); switch (dim) { case 1: - coarseningManager = NEW CoarseningManager1d(); - refinementManager = NEW RefinementManager1d(); + coarseningManager = new CoarseningManager1d(); + refinementManager = new RefinementManager1d(); break; case 2: - coarseningManager = NEW CoarseningManager2d(); - refinementManager = NEW RefinementManager2d(); + coarseningManager = new CoarseningManager2d(); + refinementManager = new RefinementManager2d(); break; case 3: - coarseningManager = NEW CoarseningManager3d(); - refinementManager = NEW RefinementManager3d(); + coarseningManager = new CoarseningManager3d(); + refinementManager = new RefinementManager3d(); break; default: ERROR_EXIT("invalid dim!\n"); @@ -445,9 +444,9 @@ namespace AMDiS { void ProblemScal::createMatricesAndVectors() { // === create vectors and system matrix === - systemMatrix = NEW DOFMatrix(feSpace, feSpace, "system matrix"); - rhs = NEW DOFVector<double>(feSpace, "rhs"); - solution = NEW DOFVector<double>(feSpace, "solution"); + systemMatrix = new DOFMatrix(feSpace, feSpace, "system matrix"); + rhs = new DOFVector<double>(feSpace, "rhs"); + solution = new DOFVector<double>(feSpace, "solution"); solution->setCoarsenOperation(COARSE_INTERPOL); solution->set(0.0); /* initialize u_h ! */ @@ -469,7 +468,7 @@ namespace AMDiS { void ProblemScal::createEstimator() { // create and set leaf data prototype - mesh->setElementDataPrototype(NEW LeafDataEstimatable(NEW LeafDataCoarsenable)); + mesh->setElementDataPrototype(new LeafDataEstimatable(new LeafDataCoarsenable)); // === create estimator === std::string estimatorType("0"); @@ -496,13 +495,12 @@ namespace AMDiS { void ProblemScal::createFileWriter() { - fileWriters.push_back(NEW FileWriter(name + "->output", mesh, solution)); + fileWriters.push_back(new FileWriter(name + "->output", mesh, solution)); int writeSerialization = 0; GET_PARAMETER(0, name + "->output->write serialization", "%d", &writeSerialization); - if (writeSerialization) { - fileWriters.push_back(NEW Serializer<ProblemScal>(this)); - } + if (writeSerialization) + fileWriters.push_back(new Serializer<ProblemScal>(this)); } void ProblemScal::estimate(AdaptInfo *adaptInfo) @@ -516,16 +514,9 @@ namespace AMDiS { TIME_USED(first,clock())); adaptInfo->setEstSum(estimator->getErrorSum(), 0); - - adaptInfo-> - setTimeEstSum(estimator->getTimeEst(), 0); - - adaptInfo-> - setEstMax(estimator->getErrorMax(), 0); - - adaptInfo-> - setTimeEstMax(estimator->getTimeEstMax(), 0); - + adaptInfo->setTimeEstSum(estimator->getTimeEst(), 0); + adaptInfo->setEstMax(estimator->getErrorMax(), 0); + adaptInfo->setTimeEstMax(estimator->getTimeEstMax(), 0); } else { WARNING("no estimator\n"); } @@ -558,6 +549,7 @@ namespace AMDiS { // Correct dimensionality of matrix base_matrix.change_dim(feSpace->getAdmin()->getUsedSize(), feSpace->getAdmin()->getUsedSize()); + set_to_zero(base_matrix); // Reuse old sparsity information (if available) or default systemMatrix->startInsertion(nnz_per_row); diff --git a/AMDiS/src/ProblemVec.cc b/AMDiS/src/ProblemVec.cc index 1226a8e731c893dc65da71884be55a1528f7102e..3d379374b77b576b9003fa18b50b4dbc2a522e5d 100644 --- a/AMDiS/src/ProblemVec.cc +++ b/AMDiS/src/ProblemVec.cc @@ -57,9 +57,8 @@ namespace AMDiS { TEST_EXIT(meshes.size() == 1)("Daran muss ich noch arbeiten!\n"); componentMeshes.resize(nComponents); - for (int i = adoptProblem->getNumComponents(); i < nComponents; i++) { + for (int i = adoptProblem->getNumComponents(); i < nComponents; i++) componentMeshes[i] = componentMeshes[0]; - } } } @@ -579,11 +578,14 @@ namespace AMDiS { { FUNCNAME("ProblemVec::refineMesh()"); + // if (adaptInfo->getTimestepNumber() == 9) return 0; + int nMeshes = static_cast<int>(meshes.size()); Flag refineFlag = 0; - for (int i = 0; i < nMeshes; i++) { - refineFlag |= refinementManager->refineMesh(meshes[i]); - } + for (int i = 0; i < nMeshes; i++) + if (adaptInfo->isRefinementAllowed(i)) + refineFlag |= refinementManager->refineMesh(meshes[i]); + return refineFlag; } @@ -593,13 +595,10 @@ namespace AMDiS { int nMeshes = static_cast<int>(meshes.size()); Flag coarsenFlag = 0; - for (int i = 0; i < nMeshes; i++) { - if (adaptInfo->isCoarseningAllowed(i)) { + for (int i = 0; i < nMeshes; i++) + if (adaptInfo->isCoarseningAllowed(i)) coarsenFlag |= coarseningManager->coarsenMesh(meshes[i]); - WARNING("coarsening for component %d no allowed\n", i); - } - } return coarsenFlag; } @@ -636,9 +635,8 @@ namespace AMDiS { double wtime = omp_get_wtime(); #endif - for (int i = 0; i < static_cast<int>(meshes.size()); i++) { + for (int i = 0; i < static_cast<int>(meshes.size()); i++) meshes[i]->dofCompress(); - } Flag assembleFlag = flag | diff --git a/AMDiS/src/RCNeighbourList.h b/AMDiS/src/RCNeighbourList.h index 7e49f23004185fb80ae9beaae5737393f4b133ed..797c96065d70c7ec4cca18bc55a5b6aa64bf4ff5 100644 --- a/AMDiS/src/RCNeighbourList.h +++ b/AMDiS/src/RCNeighbourList.h @@ -48,32 +48,28 @@ namespace AMDiS { class RCNeighbourList { public: - /** \brief - * Constructs a RCNeighbourList of size maxEdgeNeigh - */ + /// Constructs a RCNeighbourList of size maxEdgeNeigh RCNeighbourList(int maxEdgeNeigh); RCNeighbourList() {} - /** \brief - * Destructor - */ + /// Destructor virtual ~RCNeighbourList(); - /** \brief - * Sets flag of \ref rclist[i] = true - */ - inline void setCoarsePatch(int i) {rclist[i]->flag=true;}; + /// Sets flag of \ref rclist[i] = true + inline void setCoarsePatch(int i) { + rclist[i]->flag = true; + } - /** \brief - * Sets flag of \ref rclist[i] = f - */ - inline void setCoarsePatch(int i,bool f) {rclist[i]->flag=f;}; + /// Sets flag of \ref rclist[i] = f + inline void setCoarsePatch(int i, bool f) { + rclist[i]->flag = f; + } - /** \brief - * Returns \ref rclist[i].flag - */ - inline const bool isPatchCoarse(int i) const {return rclist[i]->flag;}; + /// Returns \ref rclist[i].flag + inline const bool isPatchCoarse(int i) const { + return rclist[i]->flag; + } /** \brief * If \ref rclist[i].neigh[j] is not a NULL pointer @@ -93,27 +89,21 @@ namespace AMDiS { return rclist[i]->neigh[j]? rclist[i]->neigh[j]->el : NULL; } - /** \brief - * Returns \ref rclist[i].el - */ + /// Returns \ref rclist[i].el inline Element* getElement(int i) const { if (static_cast<int>(rclist.size()) <= i) return NULL; return rclist[i]->el; } - /** \brief - * Sets \ref rclist[i].el to el and \ref rclist[i].flag to cp. - */ + /// Sets \ref rclist[i].el to el and \ref rclist[i].flag to cp. inline const Element* setElement(int i, const Element* el, bool cp = false) { rclist[i]->el = const_cast<Element*>(el); rclist[i]->flag = cp; return el; } - /** \brief - * Returns \ref rclist[i].elType - */ + /// Returns \ref rclist[i].elType inline int getType(int i) const { return rclist[i]->elType; } @@ -122,33 +112,28 @@ namespace AMDiS { rclist[i]->elType = type; } - /** \brief - * If patch can be coarsend return true, else false and reset the element - * marks - */ + /// If patch can be coarsend return true, else false and reset the element marks. virtual bool doCoarsePatch(int n_neigh); - /** \brief - * Sets \ref rclist[i].oppVertex[j] = k - */ - inline void setOppVertex(int i,int j,int k) {rclist[i]->oppVertex[j]=k;}; + /// Sets \ref rclist[i].oppVertex[j] = k + inline void setOppVertex(int i,int j,int k) { + rclist[i]->oppVertex[j] = k; + } - /** \brief - * Returns \ref rclist[i].oppVertex[j] - */ - inline int getOppVertex(int i, int j) { return rclist[i]->oppVertex[j]; }; + /// Returns \ref rclist[i].oppVertex[j] + inline int getOppVertex(int i, int j) { + return rclist[i]->oppVertex[j]; + } - /** \brief - * Sets \ref rclist[i].elType = t - */ - inline void setElType(int i,unsigned char t) {rclist[i]->elType=t;}; + /// Sets \ref rclist[i].elType = t + inline void setElType(int i,unsigned char t) { + rclist[i]->elType = t; + } - /** \brief - * Sets \ref coarseningManager = cm - */ + /// Sets \ref coarseningManager = cm inline void setCoarseningManager(CoarseningManager *cm) { coarseningManager = cm; - }; + } /** \brief * Fills \ref rclist[i].neigh and \ref rclist[i].oppVertex infos ( 0 <= i @@ -168,14 +153,10 @@ namespace AMDiS { */ void addDOFParents(int n_neigh); - /** \brief - * Removes DOFs during refinement (3d) - */ + /// Removes DOFs during refinement (3d) void removeDOFParent(int index); - /** \brief - * Removes DOFs during refinement (2d) - */ + /// Removes DOFs during refinement (2d) void removeDOFParents(int n_neigh); RCNeighbourList *periodicSplit(DegreeOfFreedom *edge[2], @@ -184,20 +165,14 @@ namespace AMDiS { int *n_neigh_periodic); protected: - /** \brief - * Information about one Element of the patch - */ + /// Information about one Element of the patch class RCListElement { public: - /** \brief - * Pointer to the Element - */ + /// Pointer to the Element Element* el; - /** \brief - * This is the no-th entry in the RCNeighbourList - */ + /// This is the no-th entry in the RCNeighbourList int no; /** \brief @@ -213,9 +188,7 @@ namespace AMDiS { */ RCListElement* neigh[2]; - /** \brief - * opp vertex[0/1] the opposite vertex of neigh[0/1] (only 3d) - */ + /// opp vertex[0/1] the opposite vertex of neigh[0/1] (only 3d) int oppVertex[2]; /** \brief @@ -226,21 +199,12 @@ namespace AMDiS { * RCListElement vector (only 3d) */ unsigned char elType; - - // /** \brief - // * Stores whether coordinate interpolation must be done. - // */ - // bool newCoords; }; - /** \brief - * Refinement/coarsening patch - */ + /// Refinement/coarsening patch std::vector<RCListElement*> rclist; - /** \brief - * Pointer to the CoarseningManager - */ + /// Pointer to the CoarseningManager CoarseningManager *coarseningManager; }; diff --git a/AMDiS/src/RefinementManager.cc b/AMDiS/src/RefinementManager.cc index 2ad6ea82787ab3dbd6803599437b997a2a5f11ca..7971f5281f0cf1c26df28e61d5cb0c8876e62cb4 100644 --- a/AMDiS/src/RefinementManager.cc +++ b/AMDiS/src/RefinementManager.cc @@ -34,8 +34,7 @@ namespace AMDiS { return static_cast<Flag>(0); globalMark = mark; - aMesh->traverse(-1, - Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_BOUND, + aMesh->traverse(-1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS | Mesh::FILL_BOUND, globalRefineFct); return refineMesh(aMesh); } diff --git a/AMDiS/src/RefinementManager2d.cc b/AMDiS/src/RefinementManager2d.cc index ef9681c6dd2361dae66b11e864c72ab249e44a4a..945e819f9bf3647ad4335a9583655d0b31de44f5 100644 --- a/AMDiS/src/RefinementManager2d.cc +++ b/AMDiS/src/RefinementManager2d.cc @@ -19,20 +19,17 @@ namespace AMDiS { { FUNCNAME("RefinementManager::refineFunction()"); - int n_neigh; bool bound = false; - DegreeOfFreedom *edge[2]; - - RCNeighbourList* refineList = NEW RCNeighbourList(2); + RCNeighbourList* refineList = new RCNeighbourList(2); if (el_info->getElement()->getMark() <= 0) { - DELETE refineList; + delete refineList; return el_info; /* element may not be refined */ } refineList->setElement(0, el_info->getElement()); - n_neigh = 1; + int n_neigh = 1; if (el_info->getProjection(0) && el_info->getProjection(0)->getType() == VOLUME_PROJECTION) { @@ -44,11 +41,11 @@ namespace AMDiS { /****************************************************************************/ if (el_info->getElement()->getDOF(0,0) < el_info->getElement()->getDOF(1, 0)) { - edge[0] = const_cast<int*>( el_info->getElement()->getDOF(0)); - edge[1] = const_cast<int*>( el_info->getElement()->getDOF(1)); + edge[0] = const_cast<int*>(el_info->getElement()->getDOF(0)); + edge[1] = const_cast<int*>(el_info->getElement()->getDOF(1)); } else { - edge[1] = const_cast<int*>( el_info->getElement()->getDOF(0)); - edge[0] = const_cast<int*>( el_info->getElement()->getDOF(1)); + edge[1] = const_cast<int*>(el_info->getElement()->getDOF(0)); + edge[0] = const_cast<int*>(el_info->getElement()->getDOF(1)); } /****************************************************************************/ @@ -89,11 +86,10 @@ namespace AMDiS { newDOF = refinePatch(edge, periodicList, n_neigh_periodic, bound); - DELETE periodicList; + delete periodicList; - if (firstNewDOF == -1) { + if (firstNewDOF == -1) firstNewDOF = newDOF; - } if (lastNewDOF != -1) { for (it = mesh->getPeriodicAssociations().begin(); it != end; ++it) { @@ -138,7 +134,7 @@ namespace AMDiS { /* and now refine the patch */ /****************************************************************************/ - DELETE refineList; + delete refineList; return el_info; } @@ -155,7 +151,7 @@ namespace AMDiS { projector = el_info->getProjection(2); if (el->getFirstChild() && projector && (!el->isNewCoordSet())) { - WorldVector<double> *new_coord = NEW WorldVector<double>; + WorldVector<double> *new_coord = new WorldVector<double>; for (int j = 0; j < dow; j++) (*new_coord)[j] = (el_info->getCoord(0)[j] + el_info->getCoord(1)[j])*0.5; @@ -264,8 +260,7 @@ namespace AMDiS { } - void RefinementManager2d::bisectTriangle(Triangle *el, - DegreeOfFreedom* newDOFs[3]) + void RefinementManager2d::bisectTriangle(Triangle *el, DegreeOfFreedom* newDOFs[3]) { FUNCNAME("RefinementManager2d::bisectTriangle()"); TEST_EXIT_DBG(mesh)("no mesh!\n"); @@ -274,8 +269,8 @@ namespace AMDiS { child[0] = dynamic_cast<Triangle*>(mesh->createNewElement(el)); child[1] = dynamic_cast<Triangle*>(mesh->createNewElement(el)); - - signed char newMark = max(0, el->getMark()-1); + + signed char newMark = max(0, el->getMark() - 1); child[0]->setMark(newMark); child[1]->setMark(newMark); el->setMark(0); @@ -303,8 +298,8 @@ namespace AMDiS { /****************************************************************************/ for (int i_child = 0; i_child < 2; i_child++) { - child[i_child]->setDOF(i_child, const_cast<int*>( el->getDOF(2))); - child[i_child]->setDOF(1-i_child, const_cast<int*>( el->getDOF(i_child))); + child[i_child]->setDOF(i_child, const_cast<int*>(el->getDOF(2))); + child[i_child]->setDOF(1 - i_child, const_cast<int*>(el->getDOF(i_child))); } /****************************************************************************/ diff --git a/AMDiS/src/SecondOrderAssembler.cc b/AMDiS/src/SecondOrderAssembler.cc index 1a9e6ed80d14ec98a84d3213f228521b9b108a05..acd1e4d71ff6c889d6b275136ee62eaa7682c117 100644 --- a/AMDiS/src/SecondOrderAssembler.cc +++ b/AMDiS/src/SecondOrderAssembler.cc @@ -45,7 +45,7 @@ namespace AMDiS { sort(opTerms.begin(), opTerms.end()); // check if a new assembler is needed - for (int i = 0; i < static_cast<int>( subAssemblers->size()); i++) { + for (int i = 0; i < static_cast<int>(subAssemblers->size()); i++) { std::vector<OperatorTerm*> assTerms = *((*subAssemblers)[i]->getTerms()); sort(assTerms.begin(), assTerms.end()); @@ -56,7 +56,7 @@ namespace AMDiS { // check if all terms are pw_const bool pwConst = true; - for (int i = 0; i < static_cast<int>( op->secondOrder[myRank].size()); i++) { + for (int i = 0; i < static_cast<int>(op->secondOrder[myRank].size()); i++) { if (!op->secondOrder[myRank][i]->isPWConst()) { pwConst = false; break; @@ -204,6 +204,10 @@ namespace AMDiS { for (int i = 0; i < static_cast<int>(terms[myRank].size()); i++) { (static_cast<SecondOrderTerm*>(terms[myRank][i]))->getLALt(elInfo, nPoints, LALt); } + +// for (int i = 0; i < nPoints; i++) { +// LALt[i]->print(); +// } VectorOfFixVecs< DimVec<double> > *grdPsi, *grdPhi; diff --git a/AMDiS/src/SolverMatrix.h b/AMDiS/src/SolverMatrix.h index 8e6fef89e7e19176f84be86e73bae4ed49d9896d..4d40eda505dad4d047b95f5551b2a40fb0c72df1 100644 --- a/AMDiS/src/SolverMatrix.h +++ b/AMDiS/src/SolverMatrix.h @@ -68,10 +68,11 @@ namespace AMDiS { block_starts[i+1]= block_starts[i] + A[i][i]->getFESpace()->getAdmin()->getUsedSize(); matrix.change_dim(block_starts[ns], block_starts[ns]); + set_to_zero(matrix); DOFMatrix::inserter_type ins(matrix); for (int rb= 0; rb < ns; ++rb) for (int cb= 0; cb < ns; ++cb) - if (A[rb][cb]) + if (A[rb][cb]) ins[block_starts[rb]][block_starts[cb]] << A[rb][cb]->getBaseMatrix(); } diff --git a/AMDiS/src/SubAssembler.cc b/AMDiS/src/SubAssembler.cc index f1c33085295b9b7f5f78717f3c6549f310690a82..3c238fd857ac2865af56e3f0dfa643252c43604f 100644 --- a/AMDiS/src/SubAssembler.cc +++ b/AMDiS/src/SubAssembler.cc @@ -94,15 +94,13 @@ namespace AMDiS { // set values at QPs invalid std::map<const DOFVectorBase<double>*, ValuesAtQPs*>::iterator it1; - for (it1 = valuesAtQPs.begin(); it1 != valuesAtQPs.end(); ++it1) { + for (it1 = valuesAtQPs.begin(); it1 != valuesAtQPs.end(); ++it1) ((*it1).second)->valid = false; - } // set gradients at QPs invalid std::map<const DOFVectorBase<double>*, GradientsAtQPs*>::iterator it2; - for (it2 = gradientsAtQPs.begin(); it2 != gradientsAtQPs.end(); ++it2) { + for (it2 = gradientsAtQPs.begin(); it2 != gradientsAtQPs.end(); ++it2) ((*it2).second)->valid = false; - } int myRank = omp_get_thread_num(); // calls initElement of each term @@ -255,9 +253,9 @@ namespace AMDiS { Quadrature *localQuad = quad ? quad : quadrature; - if (!gradientsAtQPs[vec]) { + if (!gradientsAtQPs[vec]) gradientsAtQPs[vec] = new GradientsAtQPs; - } + gradientsAtQPs[vec]->values.resize(localQuad->getNumPoints()); WorldVector<double> *values = gradientsAtQPs[vec]->values.getValArray(); diff --git a/AMDiS/src/ZeroOrderAssembler.cc b/AMDiS/src/ZeroOrderAssembler.cc index e3d26969f2f70e743ce1892668e40a68ffaad4d2..1e5210d75dd91f50d235b56269485c0d39223a65 100644 --- a/AMDiS/src/ZeroOrderAssembler.cc +++ b/AMDiS/src/ZeroOrderAssembler.cc @@ -29,8 +29,9 @@ namespace AMDiS { int myRank = omp_get_thread_num(); // check if a assembler is needed at all - if (op->zeroOrder[myRank].size() == 0) + if (op->zeroOrder[myRank].size() == 0) { return NULL; + } ZeroOrderAssembler *newAssembler; @@ -57,7 +58,7 @@ namespace AMDiS { // check if all terms are pw_const bool pwConst = true; - for (int i = 0; i < static_cast<int>( op->zeroOrder[myRank].size()); i++) { + for (int i = 0; i < static_cast<int>(op->zeroOrder[myRank].size()); i++) { if (!op->zeroOrder[myRank][i]->isPWConst()) { pwConst = false; break; @@ -284,7 +285,7 @@ namespace AMDiS { } c[0] *= elInfo->getDet(); - + if (symmetric) { for (int i = 0; i < nRow; i++) { mat[i][i] += c[0] * q00->getValue(i,i);