Added: protected method SIMinput::parseTopologySet used by the

<refine> and <raisorder> tags, which also allows specifying a
topology set defining the patch(es) to refine

src/SIM/SIM1D.C

@@ -52,19 +52,9 @@ bool SIM1D::parseGeometryTag (const TiXmlElement* elem)
 
   if (!strcasecmp(elem->Value(),"refine"))
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM1D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM1D* pch = nullptr;
     RealArray xi;
@@ -72,7 +62,7 @@ bool SIM1D::parseGeometryTag (const TiXmlElement* elem)
     {
       int addu = 0;
       utl::getAttribute(elem,"u",addu);
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM1D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j <<" "<< addu << std::endl;
@@ -80,7 +70,7 @@ bool SIM1D::parseGeometryTag (const TiXmlElement* elem)
         }
     }
     else
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM1D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j
@@ -94,24 +84,14 @@ bool SIM1D::parseGeometryTag (const TiXmlElement* elem)
 
   else if (!strcasecmp(elem->Value(),"raiseorder"))
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM1D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM1D* pch = nullptr;
     int addu = 0;
     utl::getAttribute(elem,"u",addu);
-    for (int j = lowpatch; j <= uppatch; j++)
+    for (int j : patches)
       if ((pch = dynamic_cast<ASM1D*>(this->getPatch(j,true))))
       {
         IFEM::cout <<"\tRaising order of P"<< j <<" "<< addu << std::endl;
@@ -462,7 +442,7 @@ bool SIM1D::addConstraint (int patch, int lndx, int ldim, int dirs, int code,
     IFEM::cout << " V"<< lndx;
   IFEM::cout <<" in direction(s) "<< dirs;
   if (code != 0) IFEM::cout <<" code = "<< code;
-  if (basis > 1) IFEM::cout <<" basis = "<< int(basis);
+  if (basis > 1) IFEM::cout <<" basis = "<< (int)basis;
 #if SP_DEBUG > 1
   std::cout << std::endl;
 #endif

src/SIM/SIM2D.C

@@ -71,7 +71,8 @@ bool SIM2D::addConnection (int master, int slave, int mIdx,
                            int sIdx, int orient, int basis,
                            bool coordCheck, int dim)
 {
-  if (orient < 0 || orient > 1) {
+  if (orient < 0 || orient > 1)
+  {
     std::cerr <<" *** SIM2D::addConnection: Invalid orientation "<< orient <<"."
               << std::endl;
     return false;
@@ -79,10 +80,10 @@ bool SIM2D::addConnection (int master, int slave, int mIdx,
 
   int lmaster = this->getLocalPatchIndex(master);
   int lslave = this->getLocalPatchIndex(slave);
-
   if (lmaster > 0 && lslave > 0)
   {
-    if (dim > 0) {
+    if (dim < 1) return false;
+
     IFEM::cout <<"\tConnecting P"<< slave <<" E"<< sIdx
                <<" to P"<< master <<" E"<< mIdx
                <<" reversed? "<< orient << std::endl;
@@ -101,7 +102,6 @@ bool SIM2D::addConnection (int master, int slave, int mIdx,
       if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck))
         return false;
   }
-  }
   else
     adm.dd.ghostConnections.insert(DomainDecomposition::Interface{master, slave,
                                                                   mIdx, sIdx,
@@ -117,19 +117,9 @@ bool SIM2D::parseGeometryTag (const TiXmlElement* elem)
 
   if (!strcasecmp(elem->Value(),"refine") && !isRefined)
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM2D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM2D* pch = nullptr;
     RealArray xi;
@@ -138,7 +128,7 @@ bool SIM2D::parseGeometryTag (const TiXmlElement* elem)
       int addu = 0, addv = 0;
       utl::getAttribute(elem,"u",addu);
       utl::getAttribute(elem,"v",addv);
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM2D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j
@@ -149,9 +139,10 @@ bool SIM2D::parseGeometryTag (const TiXmlElement* elem)
     }
     else
     {
+      // Non-uniform (graded) refinement
       int dir = 1;
       utl::getAttribute(elem,"dir",dir);
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM2D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j <<" dir="<< dir
@@ -166,25 +157,15 @@ bool SIM2D::parseGeometryTag (const TiXmlElement* elem)
 
   else if (!strcasecmp(elem->Value(),"raiseorder") && !isRefined)
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM2D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM2D* pch = nullptr;
     int addu = 0, addv = 0;
     utl::getAttribute(elem,"u",addu);
     utl::getAttribute(elem,"v",addv);
-    for (int j = lowpatch; j <= uppatch; j++)
+    for (int j : patches)
       if ((pch = dynamic_cast<ASM2D*>(this->getPatch(j,true))))
       {
         IFEM::cout <<"\tRaising order of P"<< j
@@ -653,7 +634,7 @@ bool SIM2D::addConstraint (int patch, int lndx, int ldim, int dirs, int code,
   IFEM::cout <<" in direction(s) "<< dirs;
   if (lndx < 0) IFEM::cout << (project ? " (local projected)" : " (local)");
   if (code != 0) IFEM::cout <<" code = "<< abs(code);
-  if (basis > 1) IFEM::cout <<" basis = "<< int(basis);
+  if (basis > 1) IFEM::cout <<" basis = "<< (int)basis;
 #if SP_DEBUG > 1
   std::cout << std::endl;
 #endif

src/SIM/SIM3D.C

@@ -49,7 +49,8 @@ bool SIM3D::addConnection (int master, int slave, int mIdx,
                            int sIdx, int orient, int basis,
                            bool coordCheck, int dim)
 {
-  if (orient < 0 || orient > 7) {
+  if (orient < 0 || orient > 7)
+  {
     std::cerr <<" *** SIM3D::addConnection: Invalid orientation "<< orient <<"."
               << std::endl;
     return false;
@@ -59,7 +60,8 @@ bool SIM3D::addConnection (int master, int slave, int mIdx,
   int lslave = this->getLocalPatchIndex(slave);
   if (lmaster > 0 && lslave > 0)
   {
-    if (dim > 1) {
+    if (dim < 2) return false;
+
     IFEM::cout <<"\tConnecting P"<< slave <<" F"<< sIdx
                <<" to P"<< master <<" F"<< mIdx
                <<" orient "<< orient << std::endl;
@@ -78,7 +80,7 @@ bool SIM3D::addConnection (int master, int slave, int mIdx,
       if (!spch->connectPatch(sIdx,*mpch,mIdx,orient,b,coordCheck))
         return false;
   }
-  } else
+  else
     adm.dd.ghostConnections.insert(DomainDecomposition::Interface{master, slave,
                                                                   mIdx, sIdx,
                                                                   orient, dim, basis});
@@ -93,19 +95,9 @@ bool SIM3D::parseGeometryTag (const TiXmlElement* elem)
 
   if (!strcasecmp(elem->Value(),"refine") && !isRefined)
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM3D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM3D* pch = nullptr;
     RealArray xi;
@@ -115,7 +107,7 @@ bool SIM3D::parseGeometryTag (const TiXmlElement* elem)
       utl::getAttribute(elem,"u",addu);
       utl::getAttribute(elem,"v",addv);
       utl::getAttribute(elem,"w",addw);
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM3D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j
@@ -130,7 +122,7 @@ bool SIM3D::parseGeometryTag (const TiXmlElement* elem)
       // Non-uniform (graded) refinement
       int dir = 1;
       utl::getAttribute(elem,"dir",dir);
-      for (int j = lowpatch; j <= uppatch; j++)
+      for (int j : patches)
         if ((pch = dynamic_cast<ASM3D*>(this->getPatch(j,true))))
         {
           IFEM::cout <<"\tRefining P"<< j <<" dir="<< dir
@@ -145,26 +137,16 @@ bool SIM3D::parseGeometryTag (const TiXmlElement* elem)
 
   else if (!strcasecmp(elem->Value(),"raiseorder") && !isRefined)
   {
-    int lowpatch = 1, uppatch = 1;
+    IntVec patches;
-    if (utl::getAttribute(elem,"patch",lowpatch))
+    if (!this->parseTopologySet(elem,patches))
-      uppatch = lowpatch;
-    if (utl::getAttribute(elem,"lowerpatch",lowpatch))
-      uppatch = myModel.size();
-    utl::getAttribute(elem,"upperpatch",uppatch);
-
-    if (lowpatch < 1 || uppatch > nGlPatches)
-    {
-      std::cerr <<" *** SIM3D::parse: Invalid patch indices, lower="
-                << lowpatch <<" upper="<< uppatch << std::endl;
       return false;
-    }
 
     ASM3D* pch = nullptr;
     int addu = 0, addv = 0, addw = 0;
     utl::getAttribute(elem,"u",addu);
     utl::getAttribute(elem,"v",addv);
     utl::getAttribute(elem,"w",addw);
-    for (int j = lowpatch; j <= uppatch; j++)
+    for (int j : patches)
       if ((pch = dynamic_cast<ASM3D*>(this->getPatch(j,true))))
       {
         IFEM::cout <<"\tRaising order of P"<< j
@@ -745,6 +727,7 @@ ASMbase* SIM3D::readPatch (std::istream& isp, int pchInd,
     else
       pch->idx = myModel.size();
   }
+
   if (checkRHSys && pch)
     if (dynamic_cast<ASM3D*>(pch)->checkRightHandSystem())
       IFEM::cout <<"\tSwapped."<< std::endl;

src/SIM/SIMinput.C

@@ -25,6 +25,7 @@
 #include "tinyxml.h"
 #include <fstream>
 #include <sstream>
+#include <numeric>
 
 
 bool SIMinput::parseGeometryTag (const TiXmlElement* elem)
@@ -419,6 +420,16 @@ bool SIMinput::parseICTag (const TiXmlElement* elem)
 }
 
 
+bool SIMinput::parseLinSolTag (const TiXmlElement* elem)
+{
+  if (!strcasecmp(elem->Value(),"class"))
+    if (elem->FirstChild())
+      opt.setLinearSolver(elem->FirstChild()->Value());
+
+  return true;
+}
+
+
 static bool noDumpDataYet = true; //!< To read only once in adaptive loops
 
 bool SIMinput::parseOutputTag (const TiXmlElement* elem)
@@ -515,9 +526,9 @@ bool SIMinput::parse (const TiXmlElement* elem)
 
 int SIMinput::parseMaterialSet (const TiXmlElement* elem, int mindex)
 {
-  std::string set;
+  std::string setName;
-  utl::getAttribute(elem,"set",set);
+  utl::getAttribute(elem,"set",setName);
-  int code = this->getUniquePropertyCode(set,0);
+  int code = this->getUniquePropertyCode(setName,0);
 
   if (code == 0)
     utl::getAttribute(elem,"code",code);
@@ -529,6 +540,51 @@ int SIMinput::parseMaterialSet (const TiXmlElement* elem, int mindex)
 }
 
 
+bool SIMinput::parseTopologySet (const TiXmlElement* elem, IntVec& patches)
+{
+  std::string setName;
+  if (utl::getAttribute(elem,"set",setName))
+  {
+    auto tit = myEntitys.find(setName);
+    if (tit == myEntitys.end())
+    {
+      std::cerr <<" *** SIMinput::parseTopologySet: Undefined topology set \""
+                << setName <<"\"."<< std::endl;
+      return false;
+    }
+
+    patches.clear();
+    for (const auto& top : tit->second)
+      if (top.idim == (short int)this->getNoParamDim())
+        patches.push_back(top.patch);
+    if (!patches.empty())
+      return true;
+
+    std::cerr <<" *** SIMinput::parseTopologySet: Invalid topology set \""
+              << setName <<"\" (no patches in this set)."<< std::endl;
+    return false;
+  }
+
+  int lowpatch = 1, uppatch = 1;
+  if (utl::getAttribute(elem,"patch",lowpatch))
+    uppatch = lowpatch;
+  if (utl::getAttribute(elem,"lowerpatch",lowpatch))
+    uppatch = myModel.size();
+  utl::getAttribute(elem,"upperpatch",uppatch);
+
+  if (lowpatch < 1 || uppatch > nGlPatches)
+  {
+    std::cerr <<" *** SIMinput::parseTopologySet: Invalid patch indices, lower="
+              << lowpatch <<" upper="<< uppatch << std::endl;
+    return false;
+  }
+
+  patches.resize(uppatch-lowpatch+1);
+  std::iota(patches.begin(),patches.end(),lowpatch);
+  return true;
+}
+
+
 const char** SIMinput::getPrioritizedTags () const
 {
   // Tags to be parsed first, and in the order specified
@@ -764,16 +820,6 @@ bool SIMinput::parse (char* keyWord, std::istream& is)
 }
 
 
-bool SIMinput::parseLinSolTag (const TiXmlElement* elem)
-{
-  if (!strcasecmp(elem->Value(),"class"))
-    if (elem->FirstChild())
-      opt.setLinearSolver(elem->FirstChild()->Value());
-
-  return true;
-}
-
-
 bool SIMinput::createFEMmodel (char resetNumb)
 {
   if (resetNumb)
@@ -852,7 +898,7 @@ bool SIMinput::createPropertySet (const std::string& setName, int pc)
   }
 
   // Create the actual property objects that are used during simulation
-  for (auto top : tit->second)
+  for (const auto& top : tit->second)
     myProps.push_back(Property(Property::UNDEFINED,pc,
                                top.patch,top.idim,top.item));
 
@@ -1033,7 +1079,7 @@ bool SIMinput::setInitialCondition (SIMdependency* fieldHolder,
   std::map<std::string,PatchVec> basisMap;
 
   // Loop over the initial conditions
-  for (auto it : info)
+  for (const auto& it : info)
   {
     // Do we have this field?
     Vector* field = fieldHolder->getField(it.sim_field);
@@ -1083,7 +1129,7 @@ bool SIMinput::setInitialCondition (SIMdependency* fieldHolder,
   }
 
   // Clean up basis patches
-  for (auto itb : basisMap)
+  for (const auto& itb : basisMap)
     for (size_t i = 0; i < itb.second.size(); i++)
       delete itb.second[i];
 
@@ -1098,10 +1144,10 @@ bool SIMinput::setInitialConditions (SIMdependency* fieldHolder)
     fieldHolder = this;
 
   bool result = true;
-  for (auto it : myICs)
+  for (const auto& it : myICs)
     if (it.first != "nofile")
       result &= this->setInitialCondition(fieldHolder,it.first,it.second);
-    else for (auto ic : it.second)
+    else for (const auto& ic : it.second)
     {
       // Do we have this field?
       Vector* field = fieldHolder->getField(ic.sim_field);
@@ -1123,8 +1169,8 @@ bool SIMinput::setInitialConditions (SIMdependency* fieldHolder)
 
 bool SIMinput::hasIC (const std::string& name) const
 {
-  for (auto it : myICs)
+  for (const auto& it : myICs)
-    for (auto ic : it.second)
+    for (const auto& ic : it.second)
       if (ic.sim_field.find(name) == 0)
         return true;
 

src/SIM/SIMinput.h

@@ -115,6 +115,10 @@ protected:
   //! \param[in] mindex Index into problem-dependent material property container
   //! \return The property code to be associated with the material
   int parseMaterialSet(const TiXmlElement* elem, int mindex);
+  //! \brief Parses the "set" attribute of a refine/raiseorder XML-tag.
+  //! \param[in] elem The XML element extract the set name from
+  //! \param[in] patches List of patch indices of the specified set
+  bool parseTopologySet(const TiXmlElement* elem, std::vector<int>& patches);
 
   //! \brief Creates a set of Property objects.
   //! \param[in] setName Name of the topology set the property is defined on