diff --git a/AMDiS/src/Marker.cc b/AMDiS/src/Marker.cc index 1fa00a37539644c9853c0e937d6ad660cf4af6a9..57668bf931948f0aacbb083e3daf8c771bd14b61 100644 --- a/AMDiS/src/Marker.cc +++ b/AMDiS/src/Marker.cc @@ -2,7 +2,7 @@ namespace AMDiS { - Marker *Marker::createMarker(::std::string name, int row) { + Marker *Marker::createMarker(std::string name, int row) { int strategy = 0; GET_PARAMETER(0, name + "->strategy", "%d", &strategy); @@ -28,5 +28,205 @@ namespace AMDiS { return marker; } + + void Marker::initMarking(AdaptInfo *adaptInfo, Mesh *mesh) + { + elMarkRefine = 0; + elMarkCoarsen = 0; + estSum = pow(adaptInfo->getEstSum(row == -1 ? 0 : row), p); + estMax = adaptInfo->getEstMax(row == -1 ? 0 : row); + } + + void Marker::finishMarking(AdaptInfo *adaptInfo) { + FUNCNAME("Marker::finishMarking()"); + + INFO(info, 4)("%d elements marked for refinement\n", elMarkRefine); + INFO(info, 4)("%d elements marked for coarsening\n", elMarkCoarsen); + } + + void Marker::markElement(AdaptInfo *adaptInfo, ElInfo *elInfo) { + Element *el = elInfo->getElement(); + double lError = el->getEstimation(row); + + if (adaptInfo->isRefinementAllowed(row == -1 ? 0 : row) && lError > markRLimit) { + setMark(el, adaptInfo->getRefineBisections(row == -1 ? 0 : row)); + } else { + if (adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && lError <= markCLimit) { + if (!el->getElementData()->getElementData(COARSENABLE) || + (lError + el->getCoarseningEstimation(row)) <= markCLimit) { + setMark(el, -adaptInfo->getCoarseBisections(row == -1 ? 0 : row)); + } + } + } + } + + Flag Marker::markMesh(AdaptInfo *adaptInfo, Mesh *mesh) { + initMarking(adaptInfo, mesh); + + if (!adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && + !adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) { + return 0; + } + + TraverseStack stack; + ElInfo *elInfo = NULL; + elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); + while (elInfo) { + markElement(adaptInfo, elInfo); + elInfo = stack.traverseNext(elInfo); + } + + finishMarking(adaptInfo); + + Flag markFlag; + if (elMarkRefine) + markFlag = 1; + if (elMarkCoarsen) + markFlag |= 2; + + return(markFlag); + } + + + void ESMarker::initMarking(AdaptInfo *adaptInfo, Mesh *mesh) { + FUNCNAME("ESMarker::initMarking()"); + + Marker::initMarking(adaptInfo, mesh); + + double ESThetaP = pow(ESTheta, p); + double ESThetaCP = pow(ESThetaC, p); + + double epsP = pow(adaptInfo->getSpaceTolerance(row == -1 ? 0 : row), p); + + markRLimit = ESThetaP * epsP / mesh->getNumberOfLeaves(); + markCLimit = ESThetaCP * epsP / mesh->getNumberOfLeaves(); + + INFO(info, 2) + ("start mark_limits: %.3le %.3le nt=%d\n", + markRLimit, markCLimit, mesh->getNumberOfLeaves()); + } + + Flag GERSMarker::markMesh(AdaptInfo *adaptInfo, Mesh *mesh) { + FUNCNAME("GERSMarker::markMesh()"); + + initMarking(adaptInfo, mesh); + + if (!adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && + !adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) { + return 0; + } + + GERSSum = 0.0; + + double LTheta = pow(1.0 - GERSThetaStar, p); + double improv, redfac, wanted, epsP, GERSGamma; + + epsP = pow(adaptInfo->getSpaceTolerance(row == -1 ? 0 : row), p); + + if (estSum < oldErrSum) { + improv = estSum / oldErrSum; + wanted = 0.8 * epsP / estSum; + redfac = std::min((1.0 - wanted) / (1.0 - improv), 1.0); + redfac = std::max(redfac, 0.0); + + if (redfac < 1.0) { + LTheta *= redfac; + INFO(info, 1) + ("GERS: use extrapolated theta_star = %lf\n", + pow(LTheta, 1.0 / p)); + } + } + + oldErrSum = estSum; + GERSGamma = 1.0; + + if (adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) { + if (LTheta > 0) { + do { + GERSSum = 0.0; + GERSGamma -= GERSNu; + markRLimit = GERSGamma * estMax; + + TraverseStack stack; + ElInfo *elInfo = NULL; + elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); + while (elInfo) { + markElementForRefinement(adaptInfo, elInfo); + elInfo = stack.traverseNext(elInfo); + } + + } while((GERSGamma > 0) && (GERSSum < LTheta * estSum)); + } + + INFO(info, 2) + ("GERS refinement with gamma = %.3lf\n", GERSGamma); + } + + if (adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row)) { + GERSGamma = 0.3; + LTheta = GERSThetaC * epsP; + + do { + GERSSum = 0.0; + GERSGamma -= GERSNu; + markCLimit = GERSGamma * estMax; + + TraverseStack stack; + ElInfo *elInfo = NULL; + elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); + while (elInfo) { + markElementForCoarsening(adaptInfo, elInfo); + elInfo = stack.traverseNext(elInfo); + } + + INFO(info, 6) + ("coarse loop: gamma=%.3e, sum=%.3e, limit=%.3e\n", + GERSGamma, GERSSum, LTheta); + } while(GERSSum > LTheta); + + INFO(info, 2) + ("GERS coarsening with gamma = %.3lf\n", GERSGamma); + } + + finishMarking(adaptInfo); + + Flag markFlag; + if (elMarkRefine) + markFlag = 1; + if (elMarkCoarsen) + markFlag |= 2; + + return(markFlag); + } + + void GERSMarker::markElementForRefinement(AdaptInfo *adaptInfo, ElInfo *elInfo) { + Element *el = elInfo->getElement(); + double lError = el->getEstimation(row); + + if (lError > markRLimit) { + GERSSum += lError; + setMark(el, adaptInfo->getRefineBisections(row == -1 ? 0 : row)); + } + } + + void GERSMarker::markElementForCoarsening(AdaptInfo *adaptInfo, ElInfo *elInfo) { + Element *el = elInfo->getElement(); + double lError = el->getEstimation(row); + + if (el->getMark() <= 0) { + if (el->getElementData()->getElementData(COARSENABLE)) { + lError += el->getCoarseningEstimation(row); + } + + if (lError <= markCLimit) { + GERSSum += lError; + setMark(el, -adaptInfo->getCoarseBisections(row == -1 ? 0 : row)); + } else { + setMark(el, 0); + } + } + } + + } diff --git a/AMDiS/src/Marker.h b/AMDiS/src/Marker.h index 0bb7cb883f3e1c1346852c2392731e7b9cf305e8..b644d4c2cea642aa76b3b805af02b100727df45a 100644 --- a/AMDiS/src/Marker.h +++ b/AMDiS/src/Marker.h @@ -54,7 +54,7 @@ namespace AMDiS { /** \brief * Constructor. */ - Marker(::std::string name_, int row_) + Marker(std::string name_, int row_) : name(name_), row(row_), maximumMarking(false), @@ -100,68 +100,22 @@ namespace AMDiS { /** \brief * Can be used by sub classes. Called before traversal. */ - virtual void initMarking(AdaptInfo *adaptInfo, Mesh *mesh) { - elMarkRefine = 0; - elMarkCoarsen = 0; - estSum = pow(adaptInfo->getEstSum(row == -1 ? 0 : row), p); - estMax = adaptInfo->getEstMax(row == -1 ? 0 : row); - }; + virtual void initMarking(AdaptInfo *adaptInfo, Mesh *mesh); /** \brief * Can be used by sub classes. Called after traversal. */ - virtual void finishMarking(AdaptInfo *adaptInfo) { - FUNCNAME("Marker::finishMarking()"); - - INFO(info, 4)("%d elements marked for refinement\n", elMarkRefine); - INFO(info, 4)("%d elements marked for coarsening\n", elMarkCoarsen); - }; + virtual void finishMarking(AdaptInfo *adaptInfo); /** \brief * Marks one element. */ - virtual void markElement(AdaptInfo *adaptInfo, ElInfo *elInfo) { - Element *el = elInfo->getElement(); - double lError = el->getEstimation(row); - - if (adaptInfo->isRefinementAllowed(row == -1 ? 0 : row) && lError > markRLimit) { - setMark(el, adaptInfo->getRefineBisections(row == -1 ? 0 : row)); - } else { - if (adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && lError <= markCLimit) { - if (!el->getElementData()->getElementData(COARSENABLE) || - (lError + el->getCoarseningEstimation(row)) <= markCLimit) { - setMark(el, -adaptInfo->getCoarseBisections(row == -1 ? 0 : row)); - } - } - } - }; + virtual void markElement(AdaptInfo *adaptInfo, ElInfo *elInfo); /** \brief * Marking of the mesh. */ - virtual Flag markMesh(AdaptInfo *adaptInfo, Mesh *mesh) { - initMarking(adaptInfo, mesh); - - if (!adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && - !adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) { - return 0; - } - - TraverseStack stack; - ElInfo *elInfo = NULL; - elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); - while (elInfo) { - markElement(adaptInfo, elInfo); - elInfo = stack.traverseNext(elInfo); - } - - finishMarking(adaptInfo); - - Flag markFlag; - if (elMarkRefine) markFlag = 1; - if (elMarkCoarsen) markFlag |= 2; - return(markFlag); - }; + virtual Flag markMesh(AdaptInfo *adaptInfo, Mesh *mesh); /** \brief * Sets \ref maximumMarking. @@ -180,13 +134,13 @@ namespace AMDiS { /** \brief * Creates a scalr marker depending on the strategy set in parameters. */ - static Marker *createMarker(::std::string name, int row_); + static Marker *createMarker(std::string name, int row_); protected: /** \brief * Name of the scalar marker. */ - ::std::string name; + std::string name; /** \brief * Equal to -1 for scalar problems. Component number if marker is @@ -262,7 +216,7 @@ namespace AMDiS { /** \brief * Constructor. */ - GRMarker(::std::string name_, int row_) + GRMarker(std::string name_, int row_) : Marker(name_, row_) {}; @@ -294,7 +248,7 @@ namespace AMDiS { /** \brief * Constructor. */ - MSMarker(::std::string name_, int row_) + MSMarker(std::string name_, int row_) : Marker(name_, row_), MSGamma(0.5), MSGammaC(0.1) @@ -346,7 +300,7 @@ namespace AMDiS { /** \brief * Constructor. */ - ESMarker(::std::string name_, int row_) + ESMarker(std::string name_, int row_) : Marker(name_, row_), ESTheta(0.9), ESThetaC(0.2) @@ -358,23 +312,7 @@ namespace AMDiS { /** \brief * Implementation of MarkScal::initMarking(). */ - virtual void initMarking(AdaptInfo *adaptInfo, Mesh *mesh) { - FUNCNAME("ESMarker::initMarking()"); - - Marker::initMarking(adaptInfo, mesh); - - double ESThetaP = pow(ESTheta, p); - double ESThetaCP = pow(ESThetaC, p); - - double epsP = pow(adaptInfo->getSpaceTolerance(row == -1 ? 0 : row), p); - - markRLimit = ESThetaP * epsP / mesh->getNumberOfLeaves(); - markCLimit = ESThetaCP * epsP / mesh->getNumberOfLeaves(); - - INFO(info, 2) - ("start mark_limits: %.3le %.3le nt=%d\n", - markRLimit, markCLimit, mesh->getNumberOfLeaves()); - }; + virtual void initMarking(AdaptInfo *adaptInfo, Mesh *mesh); protected: /** \brief @@ -406,7 +344,7 @@ namespace AMDiS { /** \brief * Constructor. */ - GERSMarker(::std::string name_, int row_) + GERSMarker(std::string name_, int row_) : Marker(name_, row_), oldErrSum(0.0), GERSThetaStar(0.6), @@ -421,133 +359,18 @@ namespace AMDiS { /** \brief * Implementation of Marker::markMesh(). */ - virtual Flag markMesh(AdaptInfo *adaptInfo, Mesh *mesh) { - FUNCNAME("GERSMarker::markMesh()"); - - initMarking(adaptInfo, mesh); - - if(!adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row) && - !adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) - { - return 0; - } - - GERSSum = 0.0; - - double LTheta = pow(1.0 - GERSThetaStar, p); - double improv, redfac, wanted, epsP, GERSGamma; - - epsP = pow(adaptInfo->getSpaceTolerance(row == -1 ? 0 : row), p); - - if(estSum < oldErrSum) { - improv = estSum / oldErrSum; - wanted = 0.8 * epsP / estSum; - redfac = ::std::min((1.0 - wanted) / (1.0 - improv), 1.0); - redfac = ::std::max(redfac, 0.0); - - if(redfac < 1.0) { - LTheta *= redfac; - INFO(info, 1) - ("GERS: use extrapolated theta_star = %lf\n", - pow(LTheta, 1.0 / p)); - } - } - - oldErrSum = estSum; - - GERSGamma = 1.0; - - if(adaptInfo->isRefinementAllowed(row == -1 ? 0 : row)) { - if(LTheta > 0) { - do { - GERSSum = 0.0; - - GERSGamma -= GERSNu; - markRLimit = GERSGamma * estMax; - - TraverseStack stack; - ElInfo *elInfo = NULL; - elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); - while(elInfo) { - markElementForRefinement(adaptInfo, elInfo); - elInfo = stack.traverseNext(elInfo); - } - - } while((GERSGamma > 0) && (GERSSum < LTheta * estSum)); - } - - INFO(info, 2) - ("GERS refinement with gamma = %.3lf\n", GERSGamma); - } - - if(adaptInfo->isCoarseningAllowed(row == -1 ? 0 : row)) { - GERSGamma = 0.3; - LTheta = GERSThetaC * epsP; - - do { - GERSSum = 0.0; - GERSGamma -= GERSNu; - markCLimit = GERSGamma * estMax; - - TraverseStack stack; - ElInfo *elInfo = NULL; - elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL); - while(elInfo) { - markElementForCoarsening(adaptInfo, elInfo); - elInfo = stack.traverseNext(elInfo); - } - - INFO(info, 6) - ("coarse loop: gamma=%.3e, sum=%.3e, limit=%.3e\n", - GERSGamma, GERSSum, LTheta); - } while(GERSSum > LTheta); - - INFO(info, 2) - ("GERS coarsening with gamma = %.3lf\n", GERSGamma); - } - - finishMarking(adaptInfo); - - Flag markFlag; - if(elMarkRefine) markFlag = 1; - if(elMarkCoarsen) markFlag |= 2; - return(markFlag); - }; + virtual Flag markMesh(AdaptInfo *adaptInfo, Mesh *mesh); protected: /** \brief * Refinement marking function. */ - void markElementForRefinement(AdaptInfo *adaptInfo, ElInfo *elInfo) { - Element *el = elInfo->getElement(); - double lError = el->getEstimation(row); - - if (lError > markRLimit) { - GERSSum += lError; - setMark(el, adaptInfo->getRefineBisections(row == -1 ? 0 : row)); - } - }; + void markElementForRefinement(AdaptInfo *adaptInfo, ElInfo *elInfo); /** \brief * Coarsening marking function. */ - void markElementForCoarsening(AdaptInfo *adaptInfo, ElInfo *elInfo) { - Element *el = elInfo->getElement(); - double lError = el->getEstimation(row); - - if(el->getMark() <= 0) { - if(el->getElementData()->getElementData(COARSENABLE)) { - lError += el->getCoarseningEstimation(row); - } - - if(lError <= markCLimit) { - GERSSum += lError; - setMark(el, -adaptInfo->getCoarseBisections(row == -1 ? 0 : row)); - } else { - setMark(el, 0); - } - } - }; + void markElementForCoarsening(AdaptInfo *adaptInfo, ElInfo *elInfo); protected: /** \brief