Nektar++
Classes | Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
Nektar::LibUtilities::MeshPartition Class Referenceabstract

#include <MeshPartition.h>

Inheritance diagram for Nektar::LibUtilities::MeshPartition:
Inheritance graph
[legend]
Collaboration diagram for Nektar::LibUtilities::MeshPartition:
Collaboration graph
[legend]

Classes

struct  GraphEdgeProperties
 
struct  GraphVertexProperties
 
struct  MeshComposite
 
struct  MeshCurved
 
struct  MeshEdge
 
struct  MeshElement
 
struct  MeshEntity
 
struct  MeshFace
 
struct  MeshVertex
 

Public Member Functions

 MeshPartition (const SessionReaderSharedPtr &pSession)
 
virtual ~MeshPartition ()
 
void PartitionMesh (int nParts, bool shared=false)
 
void WriteLocalPartition (SessionReaderSharedPtr &pSession)
 
void WriteAllPartitions (SessionReaderSharedPtr &pSession)
 
void PrintPartInfo (std::ostream &out)
 
void GetCompositeOrdering (CompositeOrdering &composites)
 
void GetBndRegionOrdering (BndRegionOrdering &composites)
 
void GetElementIDs (const int procid, std::vector< unsigned int > &tmp)
 

Private Types

typedef std::pair< std::string, int > MeshCurvedKey
 
typedef std::vector< unsigned int > MultiWeight
 
typedef boost::adjacency_list< boost::setS, boost::vecS, boost::undirectedS, GraphVertexProperties, boost::property< boost::edge_index_t, unsigned int, GraphEdgeProperties > > BoostGraph
 
typedef boost::subgraph< BoostGraphBoostSubGraph
 
typedef boost::graph_traits< BoostGraph >::vertex_descriptor BoostVertex
 
typedef boost::graph_traits< BoostGraph >::edge_descriptor BoostEdge
 
typedef boost::graph_traits< BoostGraph >::edge_iterator BoostEdgeIterator
 
typedef boost::graph_traits< BoostGraph >::vertex_iterator BoostVertexIterator
 
typedef boost::graph_traits< BoostGraph >::adjacency_iterator BoostAdjacencyIterator
 
typedef std::vector< unsigned int > NumModes
 
typedef std::map< std::string, NumModesNummodesPerField
 

Private Member Functions

void ReadExpansions (const SessionReaderSharedPtr &pSession)
 
void ReadGeometry (const SessionReaderSharedPtr &pSession)
 
void ReadConditions (const SessionReaderSharedPtr &pSession)
 
void WeightElements ()
 
void CreateGraph (BoostSubGraph &pGraph)
 
void PartitionGraph (BoostSubGraph &pGraph, int nParts, std::vector< BoostSubGraph > &pLocalPartition)
 
virtual void PartitionGraphImpl (int &nVerts, int &nVertConds, Nektar::Array< Nektar::OneD, int > &xadj, Nektar::Array< Nektar::OneD, int > &adjcy, Nektar::Array< Nektar::OneD, int > &vertWgt, Nektar::Array< Nektar::OneD, int > &vertSize, int &nparts, int &volume, Nektar::Array< Nektar::OneD, int > &part)=0
 
void OutputPartition (SessionReaderSharedPtr &pSession, BoostSubGraph &pGraph, TiXmlElement *pGeometry)
 
void CheckPartitions (int nParts, Array< OneD, int > &pPart)
 
int CalculateElementWeight (char elmtType, bool bndWeight, int na, int nb, int nc)
 

Private Attributes

int m_dim
 
int m_numFields
 
std::map< int, MeshVertexm_meshVertices
 
std::map< int, MeshEntitym_meshEdges
 
std::map< int, MeshEntitym_meshFaces
 
std::map< int, MeshEntitym_meshElements
 
std::map< MeshCurvedKey, MeshCurvedm_meshCurved
 
std::map< int, MeshEntitym_meshComposites
 
std::vector< unsigned int > m_domain
 
std::map< std::string, std::string > m_vertexAttributes
 
std::map< int, NummodesPerFieldm_expansions
 
std::map< std::string, int > m_fieldNameToId
 
std::map< int, MultiWeightm_vertWeights
 
BndRegionOrdering m_bndRegOrder
 
BoostSubGraph m_mesh
 
std::vector< BoostSubGraphm_localPartition
 
CommSharedPtr m_comm
 
bool m_weightingRequired
 
bool m_shared
 

Detailed Description

Definition at line 61 of file MeshPartition.h.

Member Typedef Documentation

typedef boost::graph_traits< BoostGraph >::adjacency_iterator Nektar::LibUtilities::MeshPartition::BoostAdjacencyIterator
private

Definition at line 185 of file MeshPartition.h.

typedef boost::graph_traits< BoostGraph >::edge_descriptor Nektar::LibUtilities::MeshPartition::BoostEdge
private

Definition at line 176 of file MeshPartition.h.

typedef boost::graph_traits< BoostGraph >::edge_iterator Nektar::LibUtilities::MeshPartition::BoostEdgeIterator
private

Definition at line 179 of file MeshPartition.h.

typedef boost::adjacency_list< boost::setS, boost::vecS, boost::undirectedS, GraphVertexProperties, boost::property< boost::edge_index_t, unsigned int, GraphEdgeProperties > > Nektar::LibUtilities::MeshPartition::BoostGraph
private

Definition at line 166 of file MeshPartition.h.

Definition at line 169 of file MeshPartition.h.

typedef boost::graph_traits< BoostGraph >::vertex_descriptor Nektar::LibUtilities::MeshPartition::BoostVertex
private

Definition at line 173 of file MeshPartition.h.

typedef boost::graph_traits< BoostGraph >::vertex_iterator Nektar::LibUtilities::MeshPartition::BoostVertexIterator
private

Definition at line 182 of file MeshPartition.h.

typedef std::pair<std::string, int> Nektar::LibUtilities::MeshPartition::MeshCurvedKey
private

Definition at line 128 of file MeshPartition.h.

typedef std::vector<unsigned int> Nektar::LibUtilities::MeshPartition::MultiWeight
private

Definition at line 137 of file MeshPartition.h.

typedef std::vector<unsigned int> Nektar::LibUtilities::MeshPartition::NumModes
private

Definition at line 187 of file MeshPartition.h.

typedef std::map<std::string, NumModes> Nektar::LibUtilities::MeshPartition::NummodesPerField
private

Definition at line 188 of file MeshPartition.h.

Constructor & Destructor Documentation

Nektar::LibUtilities::MeshPartition::MeshPartition ( const SessionReaderSharedPtr pSession)

Definition at line 74 of file MeshPartition.cpp.

References ReadConditions(), ReadExpansions(), and ReadGeometry().

74  :
75  m_numFields(0),
77  m_comm(pSession->GetComm()),
78  m_weightingRequired(false)
79  {
80  ReadConditions(pSession);
81  ReadGeometry(pSession);
82  ReadExpansions(pSession);
83  }
void ReadConditions(const SessionReaderSharedPtr &pSession)
std::map< std::string, int > m_fieldNameToId
void ReadExpansions(const SessionReaderSharedPtr &pSession)
void ReadGeometry(const SessionReaderSharedPtr &pSession)
Nektar::LibUtilities::MeshPartition::~MeshPartition ( )
virtual

Definition at line 85 of file MeshPartition.cpp.

86  {
87 
88  }

Member Function Documentation

int Nektar::LibUtilities::MeshPartition::CalculateElementWeight ( char  elmtType,
bool  bndWeight,
int  na,
int  nb,
int  nc 
)
private

Definition at line 1264 of file MeshPartition.cpp.

References Nektar::LibUtilities::StdSegData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdTriData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdQuadData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdHexData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdTetData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdPyrData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdPrismData::getNumberOfBndCoefficients(), Nektar::LibUtilities::StdSegData::getNumberOfCoefficients(), Nektar::LibUtilities::StdTriData::getNumberOfCoefficients(), Nektar::LibUtilities::StdQuadData::getNumberOfCoefficients(), Nektar::LibUtilities::StdHexData::getNumberOfCoefficients(), Nektar::LibUtilities::StdTetData::getNumberOfCoefficients(), Nektar::LibUtilities::StdPyrData::getNumberOfCoefficients(), and Nektar::LibUtilities::StdPrismData::getNumberOfCoefficients().

Referenced by PrintPartInfo(), and WeightElements().

1270  {
1271  int weight = 0;
1272 
1273  switch (elmtType)
1274  {
1275  case 'A':
1276  weight = bndWeight ?
1278  StdTetData ::getNumberOfCoefficients (na, nb, nc);
1279  break;
1280  case 'R':
1281  weight = bndWeight ?
1283  StdPrismData::getNumberOfCoefficients (na, nb, nc);
1284  break;
1285  case 'H':
1286  weight = bndWeight ?
1288  StdHexData ::getNumberOfCoefficients (na, nb, nc);
1289  break;
1290  case 'P':
1291  weight = bndWeight ?
1293  StdPyrData ::getNumberOfCoefficients (na, nb, nc);
1294  break;
1295  case 'Q':
1296  weight = bndWeight ?
1298  StdQuadData ::getNumberOfCoefficients (na, nb);
1299  break;
1300  case 'T':
1301  weight = bndWeight ?
1303  StdTriData ::getNumberOfCoefficients (na, nb);
1304  break;
1305  case 'S':
1306  weight = bndWeight ?
1308  StdSegData ::getNumberOfCoefficients (na);
1309  break;
1310  case 'V':
1311  weight = 1;
1312  break;
1313  default:
1314  break;
1315  }
1316 
1317  return weight;
1318  }
int getNumberOfBndCoefficients(int Na, int Nb)
Definition: ShapeType.hpp:120
int getNumberOfBndCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:159
int getNumberOfBndCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:266
int getNumberOfBndCoefficients(int Na, int Nb)
Definition: ShapeType.hpp:139
int getNumberOfBndCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:297
int getNumberOfBndCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:209
void Nektar::LibUtilities::MeshPartition::CheckPartitions ( int  nParts,
Array< OneD, int > &  pPart 
)
private

Definition at line 815 of file MeshPartition.cpp.

Referenced by PartitionGraph().

816  {
817  unsigned int i = 0;
818  unsigned int cnt = 0;
819  bool valid = true;
820 
821  // Check that every process has at least one element assigned
822  for (i = 0; i < nParts; ++i)
823  {
824  cnt = std::count(pPart.begin(), pPart.end(), i);
825  if (cnt == 0)
826  {
827  valid = false;
828  }
829  }
830 
831  // If METIS produced an invalid partition, repartition naively.
832  // Elements are assigned to processes in a round-robin fashion.
833  // It is assumed that METIS failure only occurs when the number of
834  // elements is approx. the number of processes, so this approach
835  // should not be too inefficient communication-wise.
836  if (!valid)
837  {
838  for (i = 0; i < pPart.num_elements(); ++i)
839  {
840  pPart[i] = i % nParts;
841  }
842  }
843  }
void Nektar::LibUtilities::MeshPartition::CreateGraph ( BoostSubGraph pGraph)
private

Definition at line 657 of file MeshPartition.cpp.

References Nektar::iterator, m_meshElements, m_vertWeights, and m_weightingRequired.

Referenced by PartitionMesh().

658  {
659  // Maps edge/face to first mesh element id.
660  // On locating second mesh element id, graph edge is created instead.
661  std::map<int, int> vGraphEdges;
663  int vcnt = 0;
664 
665  for (eIt = m_meshElements.begin(); eIt != m_meshElements.end();
666  ++eIt, ++vcnt)
667  {
668  BoostVertex v = boost::add_vertex(pGraph);
669  pGraph[v].id = eIt->first;
670  pGraph[v].partition = 0;
671 
673  {
674  pGraph[v].weight = m_vertWeights[eIt->first];
675  }
676 
677  // Process element entries and add graph edges
678  for (unsigned j = 0; j < eIt->second.list.size(); ++j)
679  {
680  int eId = eIt->second.list[j];
681 
682  // Look to see if we've examined this edge/face before
683  // If so, we've got both graph vertices so add edge
684  if (vGraphEdges.find(eId) != vGraphEdges.end())
685  {
686  BoostEdge e = boost::add_edge(vcnt, vGraphEdges[eId], pGraph).first;
687  pGraph[e].id = vcnt;
688  }
689  else
690  {
691  vGraphEdges[eId] = vcnt;
692  }
693  }
694  }
695  }
boost::graph_traits< BoostGraph >::vertex_descriptor BoostVertex
std::map< int, MultiWeight > m_vertWeights
std::map< int, MeshEntity > m_meshElements
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::graph_traits< BoostGraph >::edge_descriptor BoostEdge
void Nektar::LibUtilities::MeshPartition::GetBndRegionOrdering ( BndRegionOrdering composites)

Definition at line 177 of file MeshPartition.cpp.

References m_bndRegOrder.

178  {
179  bndRegs = m_bndRegOrder;
180  }
void Nektar::LibUtilities::MeshPartition::GetCompositeOrdering ( CompositeOrdering composites)

Definition at line 167 of file MeshPartition.cpp.

References Nektar::iterator, and m_meshComposites.

168  {
170  for (it = m_meshComposites.begin();
171  it != m_meshComposites.end(); ++it)
172  {
173  composites[it->first] = it->second.list;
174  }
175  }
std::map< int, MeshEntity > m_meshComposites
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void Nektar::LibUtilities::MeshPartition::GetElementIDs ( const int  procid,
std::vector< unsigned int > &  tmp 
)

Definition at line 1249 of file MeshPartition.cpp.

References ASSERTL0, m_localPartition, and m_meshElements.

1250  {
1251  BoostVertexIterator vertit, vertit_end;
1252 
1253  ASSERTL0(procid < m_localPartition.size(),"procid is less than the number of partitions");
1254 
1255  // Populate lists of elements, edges and vertices required.
1256  for ( boost::tie(vertit, vertit_end) = boost::vertices(m_localPartition[procid]);
1257  vertit != vertit_end;
1258  ++vertit)
1259  {
1260  elmtid.push_back(m_meshElements[m_localPartition[procid][*vertit].id].id);
1261  }
1262  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< int, MeshEntity > m_meshElements
std::vector< BoostSubGraph > m_localPartition
boost::graph_traits< BoostGraph >::vertex_iterator BoostVertexIterator
void Nektar::LibUtilities::MeshPartition::OutputPartition ( LibUtilities::SessionReaderSharedPtr pSession,
BoostSubGraph pGraph,
TiXmlElement *  pGeometry 
)
private

Definition at line 846 of file MeshPartition.cpp.

References Nektar::LibUtilities::MeshPartition::MeshCurved::data, Nektar::LibUtilities::MeshPartition::MeshCurved::entityid, Nektar::LibUtilities::MeshPartition::MeshCurved::entitytype, Nektar::ParseUtils::GenerateSeqVector(), Nektar::LibUtilities::MeshPartition::MeshCurved::id, Nektar::iterator, m_bndRegOrder, m_dim, m_domain, m_meshComposites, m_meshCurved, m_meshEdges, m_meshElements, m_meshFaces, m_meshVertices, m_vertexAttributes, Nektar::LibUtilities::MeshPartition::MeshCurved::npoints, and Nektar::LibUtilities::MeshPartition::MeshCurved::type.

Referenced by WriteAllPartitions(), and WriteLocalPartition().

850  {
851  // Write Geometry data
852  std::string vDim = pSession->GetElement("Nektar/Geometry")->Attribute("DIM");
853  std::string vSpace = pSession->GetElement("Nektar/Geometry")->Attribute("SPACE");
854  std::string vPart = boost::lexical_cast<std::string>(pGraph[*boost::vertices(pGraph).first].partition);
855  TiXmlElement* vElmtGeometry = new TiXmlElement("GEOMETRY");
856  vElmtGeometry->SetAttribute("DIM", vDim);
857  vElmtGeometry->SetAttribute("SPACE", vSpace);
858  vElmtGeometry->SetAttribute("PARTITION", vPart);
859 
860  TiXmlElement *vVertex = new TiXmlElement("VERTEX");
861  TiXmlElement *vEdge = new TiXmlElement("EDGE");
862  TiXmlElement *vFace = new TiXmlElement("FACE");
863  TiXmlElement *vElement = new TiXmlElement("ELEMENT");
864  TiXmlElement *vCurved = new TiXmlElement("CURVED");
865  TiXmlElement *vComposite = new TiXmlElement("COMPOSITE");
866  TiXmlElement *vDomain = new TiXmlElement("DOMAIN");
867 
868  TiXmlElement *x;
869  TiXmlText *y;
870 
871  BoostVertexIterator vertit, vertit_end;
872  int id;
873 
874  std::map<int, MeshEntity> vComposites;
875  std::map<int, MeshEntity> vElements;
876  std::map<int, MeshEntity> vEdges;
877  std::map<int, MeshEntity> vFaces;
878  std::map<int, MeshVertex> vVertices;
882 
883  // Populate lists of elements, edges and vertices required.
884  for ( boost::tie(vertit, vertit_end) = boost::vertices(pGraph);
885  vertit != vertit_end;
886  ++vertit)
887  {
888  id = pGraph[*vertit].id;
889  vElements[id] = m_meshElements[pGraph[*vertit].id];
890  }
891 
892  std::map<int, MeshEntity> * vNext = &vElements;
893  switch (m_dim)
894  {
895  case 3:
896  {
897  // Compile list of faces
898  for (vIt = vNext->begin(); vIt != vNext->end(); vIt++)
899  {
900  for (unsigned int j = 0; j < vIt->second.list.size(); ++j)
901  {
902  id = vIt->second.list[j];
903  vFaces[id] = m_meshFaces[id];
904  }
905  }
906  vNext = &vFaces;
907  }
908  case 2:
909  {
910  // Compile list of edges
911  for (vIt = vNext->begin(); vIt != vNext->end(); vIt++)
912  {
913  for (unsigned int j = 0; j < vIt->second.list.size(); ++j)
914  {
915  id = vIt->second.list[j];
916  vEdges[id] = m_meshEdges[id];
917  }
918  }
919  vNext = &vEdges;
920  }
921  case 1:
922  {
923  // Compile list of vertices
924  for (vIt = vNext->begin(); vIt != vNext->end(); vIt++)
925  {
926  for (unsigned int j = 0; j < vIt->second.list.size(); ++j)
927  {
928  id = vIt->second.list[j];
929  vVertices[id] = m_meshVertices[id];
930  }
931  }
932  }
933  }
934 
935  // Generate XML data for these mesh entities
936  for (vVertIt = vVertices.begin(); vVertIt != vVertices.end(); vVertIt++)
937  {
938  x = new TiXmlElement("V");
939  x->SetAttribute("ID", vVertIt->first);
940  std::stringstream vCoords;
941  vCoords.precision(12);
942  vCoords << std::setw(15) << vVertIt->second.x << " "
943  << std::setw(15) << vVertIt->second.y << " "
944  << std::setw(15) << vVertIt->second.z << " ";
945  y = new TiXmlText(vCoords.str());
946  x->LinkEndChild(y);
947  vVertex->LinkEndChild(x);
948  }
949 
950  // Apply transformation attributes to VERTEX section
951  for (vAttrIt = m_vertexAttributes.begin();
952  vAttrIt != m_vertexAttributes.end();
953  ++ vAttrIt)
954  {
955  vVertex->SetAttribute(vAttrIt->first, vAttrIt->second);
956  }
957 
958  if (m_dim >= 2)
959  {
960  for (vIt = vEdges.begin(); vIt != vEdges.end(); vIt++)
961  {
962  x = new TiXmlElement("E");
963  x->SetAttribute("ID", vIt->first);
964  std::stringstream vVertices;
965  vVertices << std::setw(10) << vIt->second.list[0]
966  << std::setw(10) << vIt->second.list[1] << " ";
967  y = new TiXmlText(vVertices.str());
968  x->LinkEndChild(y);
969  vEdge->LinkEndChild(x);
970  }
971  }
972 
973  if (m_dim >= 3)
974  {
975  for (vIt = vFaces.begin(); vIt != vFaces.end(); vIt++)
976  {
977  std::string vType("F");
978  vType[0] = vIt->second.type;
979  x = new TiXmlElement(vType);
980  x->SetAttribute("ID", vIt->first);
981  std::stringstream vListStr;
982  for (unsigned int i = 0; i < vIt->second.list.size(); ++i)
983  {
984  vListStr << std::setw(10) << vIt->second.list[i];
985  }
986  vListStr << " ";
987  y = new TiXmlText(vListStr.str());
988  x->LinkEndChild(y);
989  vFace->LinkEndChild(x);
990  }
991  }
992 
993  for (vIt = vElements.begin(); vIt != vElements.end(); vIt++)
994  {
995  std::string vType("T");
996  vType[0] = vIt->second.type;
997  x = new TiXmlElement(vType.c_str());
998  x->SetAttribute("ID", vIt->first);
999  std::stringstream vEdges;
1000  for (unsigned i = 0; i < vIt->second.list.size(); ++i)
1001  {
1002  vEdges << std::setw(10) << vIt->second.list[i];
1003  }
1004  vEdges << " ";
1005  y = new TiXmlText(vEdges.str());
1006  x->LinkEndChild(y);
1007  vElement->LinkEndChild(x);
1008  }
1009 
1010  if (m_dim >= 2)
1011  {
1012  std::map<MeshCurvedKey, MeshCurved>::const_iterator vItCurve;
1013  for (vItCurve = m_meshCurved.begin();
1014  vItCurve != m_meshCurved.end();
1015  ++vItCurve)
1016  {
1017  MeshCurved c = vItCurve->second;
1018 
1019  if (vEdges.find(c.entityid) != vEdges.end() ||
1020  vFaces.find(c.entityid) != vFaces.end())
1021  {
1022  x = new TiXmlElement(c.entitytype);
1023  x->SetAttribute("ID", c.id);
1024  if (c.entitytype == "E")
1025  {
1026  x->SetAttribute("EDGEID", c.entityid);
1027  }
1028  else
1029  {
1030  x->SetAttribute("FACEID", c.entityid);
1031  }
1032  x->SetAttribute("TYPE", c.type);
1033  x->SetAttribute("NUMPOINTS", c.npoints);
1034  y = new TiXmlText(c.data);
1035  x->LinkEndChild(y);
1036  vCurved->LinkEndChild(x);
1037  }
1038  }
1039  }
1040 
1041  // Generate composites section comprising only those mesh entities
1042  // which belong to this partition.
1043  for (vIt = m_meshComposites.begin(); vIt != m_meshComposites.end(); ++vIt)
1044  {
1045  bool comma = false; // Set to true after first entity output
1046  bool range = false; // True when entity IDs form a range
1047  int last_idx = -2; // Last entity ID output
1048  std::string vCompositeStr = "";
1049  for (unsigned int j = 0; j < vIt->second.list.size(); ++j)
1050  {
1051  // Based on entity type, check if in this partition
1052  switch (vIt->second.type)
1053  {
1054  case 'V':
1055  if (vVertices.find(vIt->second.list[j]) == vVertices.end())
1056  {
1057  continue;
1058  }
1059  break;
1060  case 'E':
1061  if (vEdges.find(vIt->second.list[j]) == vEdges.end())
1062  {
1063  continue;
1064  }
1065  break;
1066  case 'F':
1067  if (vFaces.find(vIt->second.list[j]) == vFaces.end())
1068  {
1069  continue;
1070  }
1071  break;
1072  default:
1073  if (vElements.find(vIt->second.list[j]) == vElements.end())
1074  {
1075  continue;
1076  }
1077  break;
1078  }
1079 
1080  // Condense consecutive entity IDs into ranges
1081  // last_idx initially -2 to avoid error for ID=0
1082  if (last_idx + 1 == vIt->second.list[j])
1083  {
1084  last_idx++;
1085  range = true;
1086  continue;
1087  }
1088  // This entity is not in range, so close previous range with
1089  // last_idx
1090  if (range)
1091  {
1092  vCompositeStr += "-" + boost::lexical_cast<std::string>(last_idx);
1093  range = false;
1094  }
1095  // Output ID, which is either standalone, or will start a
1096  // range.
1097  vCompositeStr += comma ? "," : "";
1098  vCompositeStr += boost::lexical_cast<std::string>(vIt->second.list[j]);
1099  last_idx = vIt->second.list[j];
1100  comma = true;
1101  }
1102  // If last entity is part of a range, it must be output now
1103  if (range)
1104  {
1105  vCompositeStr += "-" + boost::lexical_cast<std::string>(last_idx);
1106  }
1107 
1108  if (vCompositeStr.length() > 0)
1109  {
1110  vComposites[vIt->first] = vIt->second;
1111  x = new TiXmlElement("C");
1112  x->SetAttribute("ID", vIt->first);
1113  vCompositeStr = "X[" + vCompositeStr + "]";
1114  vCompositeStr[0] = vIt->second.type;
1115  y = new TiXmlText(vCompositeStr.c_str());
1116  x->LinkEndChild(y);
1117  vComposite->LinkEndChild(x);
1118  }
1119  }
1120 
1121  std::string vDomainListStr;
1122  bool comma = false;
1123  for (unsigned int i = 0; i < m_domain.size(); ++i)
1124  {
1125  if (vComposites.find(m_domain[i]) != vComposites.end())
1126  {
1127  vDomainListStr += comma ? "," : "";
1128  comma = true;
1129  vDomainListStr += boost::lexical_cast<std::string>(m_domain[i]);
1130  }
1131  }
1132  vDomainListStr = "C[" + vDomainListStr + "]";
1133  TiXmlText* vDomainList = new TiXmlText(vDomainListStr);
1134  vDomain->LinkEndChild(vDomainList);
1135 
1136  vElmtGeometry->LinkEndChild(vVertex);
1137  if (m_dim >= 2)
1138  {
1139  vElmtGeometry->LinkEndChild(vEdge);
1140  }
1141  if (m_dim >= 3)
1142  {
1143  vElmtGeometry->LinkEndChild(vFace);
1144  }
1145  vElmtGeometry->LinkEndChild(vElement);
1146  if (m_dim >= 2)
1147  {
1148  vElmtGeometry->LinkEndChild(vCurved);
1149  }
1150  vElmtGeometry->LinkEndChild(vComposite);
1151  vElmtGeometry->LinkEndChild(vDomain);
1152 
1153  pNektar->LinkEndChild(vElmtGeometry);
1154 
1155  if (pSession->DefinesElement("Nektar/Conditions"))
1156  {
1157  std::set<int> vBndRegionIdList;
1158  TiXmlElement* vConditions = new TiXmlElement(*pSession->GetElement("Nektar/Conditions"));
1159  TiXmlElement* vBndRegions = vConditions->FirstChildElement("BOUNDARYREGIONS");
1160  TiXmlElement* vBndConditions = vConditions->FirstChildElement("BOUNDARYCONDITIONS");
1161  TiXmlElement* vItem;
1162 
1163  if (vBndRegions)
1164  {
1165  TiXmlElement* vNewBndRegions = new TiXmlElement("BOUNDARYREGIONS");
1166  vItem = vBndRegions->FirstChildElement();
1167  while (vItem)
1168  {
1169  std::string vSeqStr = vItem->FirstChild()->ToText()->Value();
1170  std::string::size_type indxBeg = vSeqStr.find_first_of('[') + 1;
1171  std::string::size_type indxEnd = vSeqStr.find_last_of(']') - 1;
1172  vSeqStr = vSeqStr.substr(indxBeg, indxEnd - indxBeg + 1);
1173  std::vector<unsigned int> vSeq;
1174  ParseUtils::GenerateSeqVector(vSeqStr.c_str(), vSeq);
1175  std::string vListStr;
1176  bool comma = false;
1177  for (unsigned int i = 0; i < vSeq.size(); ++i)
1178  {
1179  if (vComposites.find(vSeq[i]) != vComposites.end())
1180  {
1181  vListStr += comma ? "," : "";
1182  comma = true;
1183  vListStr += boost::lexical_cast<std::string>(vSeq[i]);
1184  }
1185  }
1186  int p = atoi(vItem->Attribute("ID"));
1187 
1188  if (vListStr.length() == 0)
1189  {
1190  TiXmlElement* tmp = vItem;
1191  vItem = vItem->NextSiblingElement();
1192  vBndRegions->RemoveChild(tmp);
1193  }
1194  else
1195  {
1196  vListStr = "C[" + vListStr + "]";
1197  TiXmlText* vList = new TiXmlText(vListStr);
1198  TiXmlElement* vNewElement = new TiXmlElement("B");
1199  vNewElement->SetAttribute("ID", p);
1200  vNewElement->LinkEndChild(vList);
1201  vNewBndRegions->LinkEndChild(vNewElement);
1202  vBndRegionIdList.insert(p);
1203  vItem = vItem->NextSiblingElement();
1204  }
1205 
1206  // Store original order of boundary region.
1207  m_bndRegOrder[p] = vSeq;
1208 
1209  }
1210  vConditions->ReplaceChild(vBndRegions, *vNewBndRegions);
1211  }
1212 
1213  if (vBndConditions)
1214  {
1215  vItem = vBndConditions->FirstChildElement();
1216  while (vItem)
1217  {
1219  if ((x = vBndRegionIdList.find(atoi(vItem->Attribute("REF")))) != vBndRegionIdList.end())
1220  {
1221  vItem->SetAttribute("REF", *x);
1222  vItem = vItem->NextSiblingElement();
1223  }
1224  else
1225  {
1226  TiXmlElement* tmp = vItem;
1227  vItem = vItem->NextSiblingElement();
1228  vBndConditions->RemoveChild(tmp);
1229  }
1230  }
1231  }
1232  pNektar->LinkEndChild(vConditions);
1233  }
1234 
1235  // Distribute other sections of the XML to each process as is.
1236  TiXmlElement* vSrc = pSession->GetElement("Nektar")
1237  ->FirstChildElement();
1238  while (vSrc)
1239  {
1240  std::string vName = boost::to_upper_copy(vSrc->ValueStr());
1241  if (vName != "GEOMETRY" && vName != "CONDITIONS")
1242  {
1243  pNektar->LinkEndChild(new TiXmlElement(*vSrc));
1244  }
1245  vSrc = vSrc->NextSiblingElement();
1246  }
1247  }
std::map< std::string, std::string > m_vertexAttributes
std::map< int, MeshEntity > m_meshComposites
std::map< int, MeshVertex > m_meshVertices
std::map< MeshCurvedKey, MeshCurved > m_meshCurved
static bool GenerateSeqVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:78
std::map< int, MeshEntity > m_meshElements
std::vector< unsigned int > m_domain
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
std::map< int, MeshEntity > m_meshEdges
std::map< int, MeshEntity > m_meshFaces
boost::graph_traits< BoostGraph >::vertex_iterator BoostVertexIterator
void Nektar::LibUtilities::MeshPartition::PartitionGraph ( BoostSubGraph pGraph,
int  nParts,
std::vector< BoostSubGraph > &  pLocalPartition 
)
private

Definition at line 697 of file MeshPartition.cpp.

References CheckPartitions(), ErrorUtil::efatal, m_comm, m_shared, m_weightingRequired, NEKERROR, and PartitionGraphImpl().

Referenced by PartitionMesh().

700  {
701  int i;
702  int nGraphVerts = boost::num_vertices(pGraph);
703  int nGraphEdges = boost::num_edges(pGraph);
704 
705  // Convert boost graph into CSR format
706  BoostVertexIterator vertit, vertit_end;
707  Array<OneD, int> part(nGraphVerts,0);
708 
709  if (m_comm->GetRowComm()->TreatAsRankZero())
710  {
711  int acnt = 0;
712  int vcnt = 0;
713  int nWeight = nGraphVerts;
714  BoostAdjacencyIterator adjvertit, adjvertit_end;
715  Array<OneD, int> xadj(nGraphVerts+1,0);
716  Array<OneD, int> adjncy(2*nGraphEdges);
717  Array<OneD, int> vwgt(nWeight, 1);
718  Array<OneD, int> vsize(nGraphVerts, 1);
719 
720  for ( boost::tie(vertit, vertit_end) = boost::vertices(pGraph);
721  vertit != vertit_end;
722  ++vertit)
723  {
724  for ( boost::tie(adjvertit, adjvertit_end) = boost::adjacent_vertices(*vertit,pGraph);
725  adjvertit != adjvertit_end;
726  ++adjvertit)
727  {
728  adjncy[acnt++] = *adjvertit;
729  }
730 
731  xadj[++vcnt] = acnt;
732 
734  {
735  vwgt[vcnt-1] = pGraph[*vertit].weight[0];
736  }
737  else
738  {
739  vwgt[vcnt-1] = 1;
740  }
741  }
742 
743  // Call Metis and partition graph
744  int vol = 0;
745 
746  try
747  {
748  //////////////////////////////////////////////////////
749  // On a cartesian communicator do mesh partiotion just on the first column
750  // so there is no doubt the partitions are all the same in all the columns
751  if(m_comm->GetColumnComm()->GetRank() == 0)
752  {
753  // Attempt partitioning using METIS.
754  int ncon = 1;
755  PartitionGraphImpl(nGraphVerts, ncon, xadj, adjncy, vwgt, vsize, nParts, vol, part);
756 
757  // Check METIS produced a valid partition and fix if not.
758  CheckPartitions(nParts, part);
759  if (!m_shared)
760  {
761  // distribute among columns
762  for (i = 1; i < m_comm->GetColumnComm()->GetSize(); ++i)
763  {
764  m_comm->GetColumnComm()->Send(i, part);
765  }
766  }
767  }
768  else
769  {
770  m_comm->GetColumnComm()->Recv(0, part);
771  }
772  if (!m_shared)
773  {
774  m_comm->GetColumnComm()->Block();
775 
776  //////////////////////////////////
777  // distribute among rows
778  for (i = 1; i < m_comm->GetRowComm()->GetSize(); ++i)
779  {
780  m_comm->GetRowComm()->Send(i, part);
781  }
782  }
783  }
784  catch (...)
785  {
787  "Error in calling metis to partition graph.");
788  }
789  }
790  else
791  {
792  m_comm->GetRowComm()->Recv(0, part);
793  }
794 
795  // Create boost subgraph for this process's partitions
796  int nCols = nParts;
797  pLocalPartition.resize(nCols);
798  for (i = 0; i < nCols; ++i)
799  {
800  pLocalPartition[i] = pGraph.create_subgraph();
801  }
802 
803  // Populate subgraph
804  i = 0;
805  for ( boost::tie(vertit, vertit_end) = boost::vertices(pGraph);
806  vertit != vertit_end;
807  ++vertit, ++i)
808  {
809  pGraph[*vertit].partition = part[i];
810  boost::add_vertex(i, pLocalPartition[part[i]]);
811  }
812  }
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
Definition: ErrorUtil.hpp:158
void CheckPartitions(int nParts, Array< OneD, int > &pPart)
virtual void PartitionGraphImpl(int &nVerts, int &nVertConds, Nektar::Array< Nektar::OneD, int > &xadj, Nektar::Array< Nektar::OneD, int > &adjcy, Nektar::Array< Nektar::OneD, int > &vertWgt, Nektar::Array< Nektar::OneD, int > &vertSize, int &nparts, int &volume, Nektar::Array< Nektar::OneD, int > &part)=0
boost::graph_traits< BoostGraph >::adjacency_iterator BoostAdjacencyIterator
boost::graph_traits< BoostGraph >::vertex_iterator BoostVertexIterator
virtual void Nektar::LibUtilities::MeshPartition::PartitionGraphImpl ( int &  nVerts,
int &  nVertConds,
Nektar::Array< Nektar::OneD, int > &  xadj,
Nektar::Array< Nektar::OneD, int > &  adjcy,
Nektar::Array< Nektar::OneD, int > &  vertWgt,
Nektar::Array< Nektar::OneD, int > &  vertSize,
int &  nparts,
int &  volume,
Nektar::Array< Nektar::OneD, int > &  part 
)
privatepure virtual
void Nektar::LibUtilities::MeshPartition::PartitionMesh ( int  nParts,
bool  shared = false 
)

Definition at line 90 of file MeshPartition.cpp.

References ASSERTL0, CreateGraph(), m_localPartition, m_mesh, m_meshElements, m_shared, m_weightingRequired, PartitionGraph(), and WeightElements().

91  {
92  ASSERTL0(m_meshElements.size() >= nParts,
93  "Too few elements for this many processes.");
94  m_shared = shared;
95 
97  {
99  }
102  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< int, MeshEntity > m_meshElements
void CreateGraph(BoostSubGraph &pGraph)
std::vector< BoostSubGraph > m_localPartition
void PartitionGraph(BoostSubGraph &pGraph, int nParts, std::vector< BoostSubGraph > &pLocalPartition)
void Nektar::LibUtilities::MeshPartition::PrintPartInfo ( std::ostream &  out)

Definition at line 470 of file MeshPartition.cpp.

References ASSERTL0, CalculateElementWeight(), Nektar::iterator, m_dim, m_domain, m_expansions, m_localPartition, m_mesh, and m_meshComposites.

471  {
472  int nElmt = boost::num_vertices(m_mesh);
473  int nPart = m_localPartition.size();
474 
475  out << "# Partition information:" << std::endl;
476  out << "# No. elements : " << nElmt << std::endl;
477  out << "# No. partitions: " << nPart << std::endl;
478  out << "# ID nElmt nLocDof nBndDof" << std::endl;
479 
480  BoostVertexIterator vertit, vertit_end;
481  std::vector<int> partElmtCount(nPart, 0);
482  std::vector<int> partLocCount (nPart, 0);
483  std::vector<int> partBndCount (nPart, 0);
484 
485  std::map<int, int> elmtSizes;
486  std::map<int, int> elmtBndSizes;
487 
488  for (unsigned int i = 0; i < m_domain.size(); ++i)
489  {
490  int cId = m_domain[i];
491  NummodesPerField npf = m_expansions[cId];
492 
493  for (NummodesPerField::iterator it = npf.begin(); it != npf.end(); ++it)
494  {
495  ASSERTL0(it->second.size() == m_dim,
496  " Number of directional" \
497  " modes in expansion spec for composite id = " +
498  boost::lexical_cast<std::string>(cId) +
499  " and field " +
500  boost::lexical_cast<std::string>(it->first) +
501  " does not correspond to mesh dimension");
502 
503  int na = it->second[0];
504  int nb = it->second[1];
505  int nc = 0;
506  if (m_dim == 3)
507  {
508  nc = it->second[2];
509  }
510 
511  int weight = CalculateElementWeight(
512  m_meshComposites[cId].type, false, na, nb, nc);
513  int bndWeight = CalculateElementWeight(
514  m_meshComposites[cId].type, true, na, nb, nc);
515 
516  for (unsigned int j = 0; j < m_meshComposites[cId].list.size(); ++j)
517  {
518  int elid = m_meshComposites[cId].list[j];
519  elmtSizes[elid] = weight;
520  elmtBndSizes[elid] = bndWeight;
521  }
522  }
523  }
524 
525  for (boost::tie(vertit, vertit_end) = boost::vertices(m_mesh);
526  vertit != vertit_end; ++vertit)
527  {
528  int partId = m_mesh[*vertit].partition;
529  partElmtCount[partId]++;
530  partLocCount [partId] += elmtSizes[m_mesh[*vertit].id];
531  partBndCount [partId] += elmtBndSizes[m_mesh[*vertit].id];
532  }
533 
534  for (int i = 0; i < nPart; ++i)
535  {
536  out << i << " " << partElmtCount[i] << " " << partLocCount[i] << " " << partBndCount[i] << std::endl;
537  }
538  }
int CalculateElementWeight(char elmtType, bool bndWeight, int na, int nb, int nc)
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< int, MeshEntity > m_meshComposites
std::map< std::string, NumModes > NummodesPerField
std::map< int, NummodesPerField > m_expansions
std::vector< unsigned int > m_domain
std::vector< BoostSubGraph > m_localPartition
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
boost::graph_traits< BoostGraph >::vertex_iterator BoostVertexIterator
void Nektar::LibUtilities::MeshPartition::ReadConditions ( const SessionReaderSharedPtr pSession)
private

Definition at line 540 of file MeshPartition.cpp.

References ASSERTL0, and m_weightingRequired.

Referenced by MeshPartition().

541  {
542  if (!pSession->DefinesElement("Nektar/Conditions/SolverInfo"))
543  {
544  // No SolverInfo = no change of default action to weight
545  // mesh graph.
546  return;
547  }
548 
549  TiXmlElement* solverInfoElement =
550  pSession->GetElement("Nektar/Conditions/SolverInfo");
551 
552  TiXmlElement* solverInfo =
553  solverInfoElement->FirstChildElement("I");
554  ASSERTL0(solverInfo, "Cannot read SolverInfo tags");
555 
556  while (solverInfo)
557  {
558  // read the property name
559  ASSERTL0(solverInfo->Attribute("PROPERTY"),
560  "Missing PROPERTY attribute in solver info "
561  "section. ");
562  std::string solverProperty =
563  solverInfo->Attribute("PROPERTY");
564  ASSERTL0(!solverProperty.empty(),
565  "Solver info properties must have a non-empty "
566  "name. ");
567  // make sure that solver property is capitalised
568  std::string solverPropertyUpper =
569  boost::to_upper_copy(solverProperty);
570 
571 
572  // read the value
573  ASSERTL0(solverInfo->Attribute("VALUE"),
574  "Missing VALUE attribute in solver info section. ");
575  std::string solverValue = solverInfo->Attribute("VALUE");
576  ASSERTL0(!solverValue.empty(),
577  "Solver info properties must have a non-empty value");
578  // make sure that property value is capitalised
579  std::string propertyValueUpper =
580  boost::to_upper_copy(solverValue);
581 
582  if (solverPropertyUpper == "WEIGHTPARTITIONS")
583  {
584  if (propertyValueUpper != "UNIFORM")
585  {
586  m_weightingRequired = true;
587  }
588  return;
589  }
590  solverInfo = solverInfo->NextSiblingElement("I");
591  }
592  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void Nektar::LibUtilities::MeshPartition::ReadExpansions ( const SessionReaderSharedPtr pSession)
private

Expansiontypes will contain plenty of data, where relevant at this stage are composite ID(s) that this expansion type describes, nummodes and a list of fields that this expansion relates to. If this does not exist the variable is only set to "DefaultVar".

Definition at line 183 of file MeshPartition.cpp.

References ASSERTL0, Nektar::ParseUtils::GenerateOrderedStringVector(), Nektar::ParseUtils::GenerateOrderedVector(), Nektar::ParseUtils::GenerateSeqVector(), m_dim, m_expansions, m_fieldNameToId, and m_numFields.

Referenced by MeshPartition().

184  {
185  // Find the Expansions tag
186  TiXmlElement *expansionTypes = pSession->GetElement("Nektar/Expansions");
187 
188  // Find the Expansion type
189  TiXmlElement *expansion = expansionTypes->FirstChildElement();
190  std::string expType = expansion->Value();
191 
192 
193  if(expType != "E")
194  {
195  ASSERTL0(false,"Expansion type not defined or not supported at the moment");
196  }
197 
198  /// Expansiontypes will contain plenty of data,
199  /// where relevant at this stage are composite
200  /// ID(s) that this expansion type describes,
201  /// nummodes and a list of fields that this
202  /// expansion relates to. If this does not exist
203  /// the variable is only set to "DefaultVar".
204 
205  while (expansion)
206  {
207  std::vector<unsigned int> composite;
208  std::vector<unsigned int> nummodes;
209  std::vector<std::string> fieldName;
210 
211  const char *nModesStr = expansion->Attribute("NUMMODES");
212  ASSERTL0(nModesStr,"NUMMODES was not defined in EXPANSION section of input");
213  std::string numModesStr = nModesStr;
214  bool valid = ParseUtils::GenerateOrderedVector(numModesStr.c_str(), nummodes);
215  ASSERTL0(valid, "Unable to correctly parse the number of modes.");
216 
217  if (nummodes.size() == 1)
218  {
219  for (int i = 1; i < m_dim; i++)
220  {
221  nummodes.push_back( nummodes[0] );
222  }
223  }
224  ASSERTL0(nummodes.size() == m_dim,"Number of modes should match mesh dimension");
225 
226 
227  const char *fStr = expansion->Attribute("FIELDS");
228  if(fStr)
229  {
230  std::string fieldStr = fStr;
231  bool valid = ParseUtils::GenerateOrderedStringVector(fieldStr.c_str(),fieldName);
232  ASSERTL0(valid,"Unable to correctly parse the field string in ExpansionTypes.");
233 
234  for (int i = 0; i < fieldName.size(); ++i)
235  {
236  if (m_fieldNameToId.count(fieldName[i]) == 0)
237  {
238  int k = m_fieldNameToId.size();
239  m_fieldNameToId[ fieldName[i] ] = k;
240  m_numFields++;
241  }
242  }
243  }
244  else
245  {
246  fieldName.push_back("DefaultVar");
247  int k = m_fieldNameToId.size();
248 
249  if (m_fieldNameToId.count("DefaultVar") == 0)
250  {
251  ASSERTL0(k == 0,
252  "Omitting field variables and explicitly listing " \
253  "them in different ExpansionTypes is wrong practise");
254 
255  m_fieldNameToId[ "DefaultVar" ] = k;
256  m_numFields++;
257  }
258  }
259 
260  std::string compositeStr = expansion->Attribute("COMPOSITE");
261  ASSERTL0(compositeStr.length() > 3, "COMPOSITE must be specified in expansion definition");
262  int beg = compositeStr.find_first_of("[");
263  int end = compositeStr.find_first_of("]");
264  std::string compositeListStr = compositeStr.substr(beg+1,end-beg-1);
265  bool parseGood = ParseUtils::GenerateSeqVector(compositeListStr.c_str(), composite);
266  ASSERTL0(parseGood && !composite.empty(),
267  (std::string("Unable to read composite index range: ") + compositeListStr).c_str());
268 
269 
270  // construct mapping (field name, CompositeID) -> nummodes
271  for (int i = 0; i < composite.size(); ++i)
272  {
273  for (int j = 0; j < fieldName.size(); j++)
274  {
275  m_expansions[composite[i]][fieldName[j]] = nummodes;
276  }
277  }
278 
279  expansion = expansion->NextSiblingElement("E");
280  }
281  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static bool GenerateOrderedStringVector(const char *const str, std::vector< std::string > &vec)
Definition: ParseUtils.hpp:142
static bool GenerateOrderedVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:96
static bool GenerateSeqVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:78
std::map< std::string, int > m_fieldNameToId
std::map< int, NummodesPerField > m_expansions
void Nektar::LibUtilities::MeshPartition::ReadGeometry ( const SessionReaderSharedPtr pSession)
private

Definition at line 286 of file MeshPartition.cpp.

References ASSERTL0, Nektar::LibUtilities::MeshPartition::MeshCurved::data, Nektar::LibUtilities::MeshPartition::MeshCurved::entityid, Nektar::LibUtilities::MeshPartition::MeshCurved::entitytype, Nektar::ParseUtils::GenerateSeqVector(), Nektar::LibUtilities::MeshPartition::MeshEntity::id, Nektar::LibUtilities::MeshPartition::MeshVertex::id, Nektar::LibUtilities::MeshPartition::MeshCurved::id, Nektar::LibUtilities::MeshPartition::MeshEntity::list, m_dim, m_domain, m_meshComposites, m_meshCurved, m_meshEdges, m_meshElements, m_meshFaces, m_meshVertices, m_vertexAttributes, Nektar::LibUtilities::MeshPartition::MeshCurved::npoints, Nektar::LibUtilities::MeshPartition::MeshEntity::type, Nektar::LibUtilities::MeshPartition::MeshCurved::type, Nektar::LibUtilities::MeshPartition::MeshVertex::x, Nektar::LibUtilities::MeshPartition::MeshVertex::y, and Nektar::LibUtilities::MeshPartition::MeshVertex::z.

Referenced by MeshPartition().

287  {
288  TiXmlElement* x;
289  TiXmlElement *vGeometry, *vSubElement;
290 
291  vGeometry = pSession->GetElement("Nektar/Geometry");
292  m_dim = atoi(vGeometry->Attribute("DIM"));
293 
294  // Read mesh vertices
295  vSubElement = pSession->GetElement("Nektar/Geometry/Vertex");
296 
297  // Retrieve any VERTEX attributes specifying mesh transforms
298  std::string attr[] = {"XSCALE", "YSCALE", "ZSCALE",
299  "XMOVE", "YMOVE", "ZMOVE" };
300  for (int i = 0; i < 6; ++i)
301  {
302  const char *val = vSubElement->Attribute(attr[i].c_str());
303  if (val)
304  {
305  m_vertexAttributes[attr[i]] = std::string(val);
306  }
307  }
308 
309  x = vSubElement->FirstChildElement();
310 
311  while(x)
312  {
313  TiXmlAttribute* y = x->FirstAttribute();
314  ASSERTL0(y, "Failed to get attribute.");
315  MeshVertex v;
316  v.id = y->IntValue();
317  std::vector<std::string> vCoords;
318  std::string vCoordStr = x->FirstChild()->ToText()->Value();
319  boost::split(vCoords, vCoordStr, boost::is_any_of("\t "));
320  v.x = atof(vCoords[0].c_str());
321  v.y = atof(vCoords[1].c_str());
322  v.z = atof(vCoords[2].c_str());
323  m_meshVertices[v.id] = v;
324  x = x->NextSiblingElement();
325  }
326 
327  // Read mesh edges
328  if (m_dim >= 2)
329  {
330  vSubElement = pSession->GetElement("Nektar/Geometry/Edge");
331  ASSERTL0(vSubElement, "Cannot read edges");
332  x = vSubElement->FirstChildElement();
333  while(x)
334  {
335  TiXmlAttribute* y = x->FirstAttribute();
336  ASSERTL0(y, "Failed to get attribute.");
337  MeshEntity e;
338  e.id = y->IntValue();
339  e.type = 'E';
340  std::vector<std::string> vVertices;
341  std::string vVerticesString = x->FirstChild()->ToText()->Value();
342  boost::split(vVertices, vVerticesString, boost::is_any_of("\t "));
343  e.list.push_back(atoi(vVertices[0].c_str()));
344  e.list.push_back(atoi(vVertices[1].c_str()));
345  m_meshEdges[e.id] = e;
346  x = x->NextSiblingElement();
347  }
348  }
349 
350  // Read mesh faces
351  if (m_dim == 3)
352  {
353  vSubElement = pSession->GetElement("Nektar/Geometry/Face");
354  ASSERTL0(vSubElement, "Cannot read faces.");
355  x = vSubElement->FirstChildElement();
356  while(x)
357  {
358  TiXmlAttribute* y = x->FirstAttribute();
359  ASSERTL0(y, "Failed to get attribute.");
360  MeshEntity f;
361  f.id = y->IntValue();
362  f.type = x->Value()[0];
363  std::vector<std::string> vEdges;
364  std::string vEdgeStr = x->FirstChild()->ToText()->Value();
365  boost::split(vEdges, vEdgeStr, boost::is_any_of("\t "));
366  for (int i = 0; i < vEdges.size(); ++i)
367  {
368  f.list.push_back(atoi(vEdges[i].c_str()));
369  }
370  m_meshFaces[f.id] = f;
371  x = x->NextSiblingElement();
372  }
373  }
374 
375  // Read mesh elements
376  vSubElement = pSession->GetElement("Nektar/Geometry/Element");
377  ASSERTL0(vSubElement, "Cannot read elements.");
378  x = vSubElement->FirstChildElement();
379  while(x)
380  {
381  TiXmlAttribute* y = x->FirstAttribute();
382  ASSERTL0(y, "Failed to get attribute.");
383  MeshEntity e;
384  e.id = y->IntValue();
385  std::vector<std::string> vItems;
386  std::string vItemStr = x->FirstChild()->ToText()->Value();
387  boost::split(vItems, vItemStr, boost::is_any_of("\t "));
388  for (int i = 0; i < vItems.size(); ++i)
389  {
390  e.list.push_back(atoi(vItems[i].c_str()));
391  }
392  e.type = x->Value()[0];
393  m_meshElements[e.id] = e;
394  x = x->NextSiblingElement();
395  }
396 
397  // Read mesh curves
398  if (pSession->DefinesElement("Nektar/Geometry/Curved"))
399  {
400  vSubElement = pSession->GetElement("Nektar/Geometry/Curved");
401  x = vSubElement->FirstChildElement();
402  while(x)
403  {
404  MeshCurved c;
405  ASSERTL0(x->Attribute("ID", &c.id),
406  "Failed to get attribute ID");
407  c.type = std::string(x->Attribute("TYPE"));
408  ASSERTL0(!c.type.empty(),
409  "Failed to get attribute TYPE");
410  ASSERTL0(x->Attribute("NUMPOINTS", &c.npoints),
411  "Failed to get attribute NUMPOINTS");
412  c.data = x->FirstChild()->ToText()->Value();
413  c.entitytype = x->Value()[0];
414  if (c.entitytype == "E")
415  {
416  ASSERTL0(x->Attribute("EDGEID", &c.entityid),
417  "Failed to get attribute EDGEID");
418  }
419  else if (c.entitytype == "F")
420  {
421  ASSERTL0(x->Attribute("FACEID", &c.entityid),
422  "Failed to get attribute FACEID");
423  }
424  else
425  {
426  ASSERTL0(false, "Unknown curve type.");
427  }
428  m_meshCurved[std::make_pair(c.entitytype, c.id)] = c;
429  x = x->NextSiblingElement();
430  }
431  }
432 
433  // Read composites
434  vSubElement = pSession->GetElement("Nektar/Geometry/Composite");
435  ASSERTL0(vSubElement, "Cannot read composites.");
436  x = vSubElement->FirstChildElement();
437  while(x)
438  {
439  TiXmlAttribute* y = x->FirstAttribute();
440  ASSERTL0(y, "Failed to get attribute.");
441  MeshEntity c;
442  c.id = y->IntValue();
443  std::string vSeqStr = x->FirstChild()->ToText()->Value();
444  c.type = vSeqStr[0];
445  std::string::size_type indxBeg = vSeqStr.find_first_of('[') + 1;
446  std::string::size_type indxEnd = vSeqStr.find_last_of(']') - 1;
447  vSeqStr = vSeqStr.substr(indxBeg, indxEnd - indxBeg + 1);
448 
449  std::vector<unsigned int> vSeq;
450  ParseUtils::GenerateSeqVector(vSeqStr.c_str(), vSeq);
451 
452  for (int i = 0; i < vSeq.size(); ++i)
453  {
454  c.list.push_back(vSeq[i]);
455  }
456  m_meshComposites[c.id] = c;
457  x = x->NextSiblingElement();
458  }
459 
460  // Read Domain
461  vSubElement = pSession->GetElement("Nektar/Geometry/Domain");
462  ASSERTL0(vSubElement, "Cannot read domain");
463  std::string vSeqStr = vSubElement->FirstChild()->ToText()->Value();
464  std::string::size_type indxBeg = vSeqStr.find_first_of('[') + 1;
465  std::string::size_type indxEnd = vSeqStr.find_last_of(']') - 1;
466  vSeqStr = vSeqStr.substr(indxBeg, indxEnd - indxBeg + 1);
467  ParseUtils::GenerateSeqVector(vSeqStr.c_str(), m_domain);
468  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< std::string, std::string > m_vertexAttributes
std::map< int, MeshEntity > m_meshComposites
std::map< int, MeshVertex > m_meshVertices
std::map< MeshCurvedKey, MeshCurved > m_meshCurved
static bool GenerateSeqVector(const char *const str, std::vector< unsigned int > &vec)
Definition: ParseUtils.hpp:78
std::map< int, MeshEntity > m_meshElements
std::vector< unsigned int > m_domain
std::map< int, MeshEntity > m_meshEdges
std::map< int, MeshEntity > m_meshFaces
void Nektar::LibUtilities::MeshPartition::WeightElements ( )
private

Definition at line 609 of file MeshPartition.cpp.

References ASSERTL0, CalculateElementWeight(), Nektar::iterator, m_dim, m_domain, m_expansions, m_fieldNameToId, m_meshComposites, m_meshElements, m_numFields, and m_vertWeights.

Referenced by PartitionMesh().

610  {
611  std::vector<unsigned int> weight(m_numFields, 1);
613  for (eIt = m_meshElements.begin(); eIt != m_meshElements.end(); ++eIt)
614  {
615  m_vertWeights[eIt->first] = weight;
616  }
617 
618  for (unsigned int i = 0; i < m_domain.size(); ++i)
619  {
620  int cId = m_domain[i];
621  NummodesPerField npf = m_expansions[cId];
622 
623  for (NummodesPerField::iterator it = npf.begin(); it != npf.end(); ++it)
624  {
625  ASSERTL0(it->second.size() == m_dim,
626  " Number of directional" \
627  " modes in expansion spec for composite id = " +
628  boost::lexical_cast<std::string>(cId) +
629  " and field " +
630  boost::lexical_cast<std::string>(it->first) +
631  " does not correspond to mesh dimension");
632 
633  int na = it->second[0];
634  int nb = 0;
635  int nc = 0;
636  if (m_dim >= 2)
637  {
638  nb = it->second[1];
639  }
640  if (m_dim == 3)
641  {
642  nc = it->second[2];
643  }
644 
645  int bndWeight = CalculateElementWeight(
646  m_meshComposites[cId].type, true, na, nb, nc);
647 
648  for (unsigned int j = 0; j < m_meshComposites[cId].list.size(); ++j)
649  {
650  int elmtId = m_meshComposites[cId].list[j];
651  m_vertWeights[elmtId][m_fieldNameToId[it->first]] = bndWeight;
652  }
653  }
654  } // for i
655  }
int CalculateElementWeight(char elmtType, bool bndWeight, int na, int nb, int nc)
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
std::map< int, MeshEntity > m_meshComposites
std::map< int, MultiWeight > m_vertWeights
std::map< std::string, NumModes > NummodesPerField
std::map< std::string, int > m_fieldNameToId
std::map< int, NummodesPerField > m_expansions
std::map< int, MeshEntity > m_meshElements
std::vector< unsigned int > m_domain
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void Nektar::LibUtilities::MeshPartition::WriteAllPartitions ( LibUtilities::SessionReaderSharedPtr pSession)

Definition at line 134 of file MeshPartition.cpp.

References m_localPartition, OutputPartition(), and Nektar::LibUtilities::PortablePath().

135  {
136  for (int i = 0; i < m_localPartition.size(); ++i)
137  {
138  TiXmlDocument vNew;
139  TiXmlDeclaration * decl = new TiXmlDeclaration("1.0", "utf-8", "");
140  vNew.LinkEndChild(decl);
141 
142  TiXmlElement* vElmtNektar;
143  vElmtNektar = new TiXmlElement("NEKTAR");
144 
145  OutputPartition(pSession, m_localPartition[i], vElmtNektar);
146 
147  vNew.LinkEndChild(vElmtNektar);
148 
149  std::string dirname = pSession->GetSessionName() + "_xml";
150  fs::path pdirname(dirname);
151 
152  boost::format pad("P%1$07d.xml");
153  pad % i;
154  fs::path pFilename(pad.str());
155 
156  fs::path fullpath = pdirname / pFilename;
157 
158  if(!fs::is_directory(dirname))
159  {
160  fs::create_directory(dirname);
161  }
162 
163  vNew.SaveFile(PortablePath(fullpath));
164  }
165  }
std::vector< BoostSubGraph > m_localPartition
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
void OutputPartition(SessionReaderSharedPtr &pSession, BoostSubGraph &pGraph, TiXmlElement *pGeometry)
void Nektar::LibUtilities::MeshPartition::WriteLocalPartition ( LibUtilities::SessionReaderSharedPtr pSession)

Definition at line 104 of file MeshPartition.cpp.

References m_comm, m_localPartition, OutputPartition(), and Nektar::LibUtilities::PortablePath().

105  {
106  TiXmlDocument vNew;
107  TiXmlDeclaration * decl = new TiXmlDeclaration("1.0", "utf-8", "");
108  vNew.LinkEndChild(decl);
109 
110  TiXmlElement* vElmtNektar;
111  vElmtNektar = new TiXmlElement("NEKTAR");
112 
113  int rank = m_comm->GetRowComm()->GetRank();
114  OutputPartition(pSession, m_localPartition[rank], vElmtNektar);
115 
116  vNew.LinkEndChild(vElmtNektar);
117 
118  std::string dirname = pSession->GetSessionName() + "_xml";
119  fs::path pdirname(dirname);
120 
121  boost::format pad("P%1$07d.xml");
122  pad % rank;
123  fs::path pFilename(pad.str());
124 
125  if(!fs::is_directory(dirname))
126  {
127  fs::create_directory(dirname);
128  }
129 
130  fs::path fullpath = pdirname / pFilename;
131  vNew.SaveFile(PortablePath(fullpath));
132  }
std::vector< BoostSubGraph > m_localPartition
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
void OutputPartition(SessionReaderSharedPtr &pSession, BoostSubGraph &pGraph, TiXmlElement *pGeometry)

Member Data Documentation

BndRegionOrdering Nektar::LibUtilities::MeshPartition::m_bndRegOrder
private

Definition at line 209 of file MeshPartition.h.

Referenced by GetBndRegionOrdering(), and OutputPartition().

CommSharedPtr Nektar::LibUtilities::MeshPartition::m_comm
private

Definition at line 214 of file MeshPartition.h.

Referenced by PartitionGraph(), and WriteLocalPartition().

int Nektar::LibUtilities::MeshPartition::m_dim
private
std::vector<unsigned int> Nektar::LibUtilities::MeshPartition::m_domain
private

Definition at line 199 of file MeshPartition.h.

Referenced by OutputPartition(), PrintPartInfo(), ReadGeometry(), and WeightElements().

std::map<int, NummodesPerField> Nektar::LibUtilities::MeshPartition::m_expansions
private

Definition at line 204 of file MeshPartition.h.

Referenced by PrintPartInfo(), ReadExpansions(), and WeightElements().

std::map<std::string, int> Nektar::LibUtilities::MeshPartition::m_fieldNameToId
private

Definition at line 206 of file MeshPartition.h.

Referenced by ReadExpansions(), and WeightElements().

std::vector<BoostSubGraph> Nektar::LibUtilities::MeshPartition::m_localPartition
private
BoostSubGraph Nektar::LibUtilities::MeshPartition::m_mesh
private

Definition at line 211 of file MeshPartition.h.

Referenced by PartitionMesh(), and PrintPartInfo().

std::map<int, MeshEntity> Nektar::LibUtilities::MeshPartition::m_meshComposites
private
std::map<MeshCurvedKey, MeshCurved> Nektar::LibUtilities::MeshPartition::m_meshCurved
private

Definition at line 197 of file MeshPartition.h.

Referenced by OutputPartition(), and ReadGeometry().

std::map<int, MeshEntity> Nektar::LibUtilities::MeshPartition::m_meshEdges
private

Definition at line 194 of file MeshPartition.h.

Referenced by OutputPartition(), and ReadGeometry().

std::map<int, MeshEntity> Nektar::LibUtilities::MeshPartition::m_meshElements
private
std::map<int, MeshEntity> Nektar::LibUtilities::MeshPartition::m_meshFaces
private

Definition at line 195 of file MeshPartition.h.

Referenced by OutputPartition(), and ReadGeometry().

std::map<int, MeshVertex> Nektar::LibUtilities::MeshPartition::m_meshVertices
private

Definition at line 193 of file MeshPartition.h.

Referenced by OutputPartition(), and ReadGeometry().

int Nektar::LibUtilities::MeshPartition::m_numFields
private

Definition at line 191 of file MeshPartition.h.

Referenced by ReadExpansions(), and WeightElements().

bool Nektar::LibUtilities::MeshPartition::m_shared
private

Definition at line 217 of file MeshPartition.h.

Referenced by PartitionGraph(), and PartitionMesh().

std::map<std::string, std::string> Nektar::LibUtilities::MeshPartition::m_vertexAttributes
private

Definition at line 200 of file MeshPartition.h.

Referenced by OutputPartition(), and ReadGeometry().

std::map<int, MultiWeight> Nektar::LibUtilities::MeshPartition::m_vertWeights
private

Definition at line 207 of file MeshPartition.h.

Referenced by CreateGraph(), and WeightElements().

bool Nektar::LibUtilities::MeshPartition::m_weightingRequired
private

Definition at line 216 of file MeshPartition.h.

Referenced by CreateGraph(), PartitionGraph(), PartitionMesh(), and ReadConditions().