49std::string
static inline ReadTag(std::string &tagStr)
51 std::string::size_type indxBeg = tagStr.find_first_of(
'[') + 1;
52 std::string::size_type indxEnd = tagStr.find_last_of(
']') - 1;
56 (std::string(
"Error reading interface region definition:") + tagStr)
59 std::string indxStr = tagStr.substr(indxBeg, indxEnd - indxBeg + 1);
66 auto length = str.length();
67 return str.substr(1, length - 2);
73 TiXmlNode *nektar = pSession->GetElement(
"NEKTAR");
74 if (nektar ==
nullptr)
79 TiXmlNode *movement = nektar->FirstChild(
"MOVEMENT");
80 if (movement !=
nullptr)
82 bool zones = movement->FirstChild(
"ZONES") !=
nullptr;
85 TiXmlElement *zonesTag =
86 pSession->GetElement(
"NEKTAR/MOVEMENT/ZONES");
90 bool interfaces = movement->FirstChild(
"INTERFACES") !=
nullptr;
93 TiXmlElement *interfacesTag =
94 pSession->GetElement(
"NEKTAR/MOVEMENT/INTERFACES");
99 "Only one of ZONES or INTERFACES present in the MOVEMENT "
103 if (pSession->DefinesSolverInfo(
"DiffusionType"))
105 ASSERTL0(pSession->GetSolverInfo(
"DiffusionType") ==
"LDGNS",
106 "Only LDGNS is supported as the DiffusionType in "
107 "SOLVERINFO when a MOVEMENT block is defined.")
112 if (movement !=
nullptr && pSession->DefinesCmdLineArgument(
"verbose"))
114 auto comm = pSession->GetComm();
115 if (comm->TreatAsRankZero())
117 std::cout <<
"Movement Info:\n";
118 std::cout <<
"\tNum zones: " <<
m_zones.size() <<
"\n";
123 auto numEl = zone.second->GetElements().size();
128 !zone.second->GetElements().empty()
129 ? zone.second->GetElements().front()->GetShapeType()
133 if (comm->TreatAsRankZero())
135 std::cout <<
"\t- " << zone.first <<
" "
137 zone.second->GetMovementType())]
138 <<
": " << numEl <<
" "
143 if (comm->TreatAsRankZero())
145 std::cout <<
"\tNum interfaces: " <<
m_interfaces.size() <<
"\n";
151 interface.second->GetLeftInterface()->GetEdge().size();
153 interface.second->GetRightInterface()->GetEdge().size();
159 !interface.second->GetLeftInterface()->GetEdge().empty()
160 ? interface.second->GetLeftInterface()
163 ->second->GetShapeType()
165 comm->GetSpaceComm()->AllReduce(shapeTypeLeft,
168 !interface.second->GetRightInterface()->GetEdge().empty()
169 ? interface.second->GetRightInterface()
172 ->second->GetShapeType()
174 comm->GetSpaceComm()->AllReduce(shapeTypeRight,
177 if (comm->TreatAsRankZero())
179 std::cout <<
"\t- \"" << interface.first.second <<
"\": "
180 << interface.second->GetLeftInterface()->GetId()
181 <<
" (" << numLeft <<
" "
184 << interface.second->GetRightInterface()->GetId()
185 <<
" (" << numRight <<
" "
191 comm->GetSpaceComm()->Block();
192 if (comm->TreatAsRankZero())
194 std::cout << std::endl;
204 ASSERTL0(zonesTag,
"Unable to find ZONES tag in file.");
205 TiXmlElement *zonesElement = zonesTag->FirstChildElement();
208 std::string zoneType = zonesElement->Value();
213 err = zonesElement->QueryIntAttribute(
"ID", &indx);
214 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read zone ID.");
216 std::string interfaceDomainStr;
217 err = zonesElement->QueryStringAttribute(
"DOMAIN", &interfaceDomainStr);
218 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read zone domain.");
221 auto domFind = std::stoi(
ReadTag(interfaceDomainStr));
222 std::map<int, CompositeSharedPtr> domain;
223 if (domains.find(domFind) != domains.end())
225 domain = domains.at(domFind);
229 if (zoneType ==
"F" || zoneType ==
"FIXED")
232 indx, domFind, domain, coordDim));
234 else if (zoneType ==
"R" || zoneType ==
"ROTATE" ||
235 zoneType ==
"ROTATING")
237 std::string originStr;
238 err = zonesElement->QueryStringAttribute(
"ORIGIN", &originStr);
239 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read origin.");
240 std::vector<NekDouble> originVec;
246 err = zonesElement->QueryStringAttribute(
"AXIS", &axisStr);
247 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read axis.");
251 std::string angularVelStr;
252 err = zonesElement->QueryStringAttribute(
"ANGVEL", &angularVelStr);
253 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read angular velocity.");
257 pSession->GetInterpreter(), angularVelStr);
260 indx, domFind, domain, coordDim,
origin,
axis, angularVelEqn));
264 else if (zoneType ==
"T" || zoneType ==
"TRANSLATE" ||
265 zoneType ==
"TRANSLATING")
267 std::string velocityStr;
268 err = zonesElement->QueryStringAttribute(
"VELOCITY", &velocityStr);
269 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read direction.");
275 indx, domFind, domain, coordDim,
velocity));
279 else if (zoneType ==
"P" || zoneType ==
"PRESCRIBED")
281 std::string xDeformStr;
282 err = zonesElement->QueryStringAttribute(
"XDEFORM", &xDeformStr);
283 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read x deform equation.");
285 std::make_shared<LibUtilities::Equation>(
286 pSession->GetInterpreter(), xDeformStr);
288 std::string yDeformStr;
289 err = zonesElement->QueryStringAttribute(
"YDEFORM", &yDeformStr);
290 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read y deform equation.");
292 std::make_shared<LibUtilities::Equation>(
293 pSession->GetInterpreter(), yDeformStr);
295 std::string zDeformStr;
296 err = zonesElement->QueryStringAttribute(
"ZDEFORM", &zDeformStr);
297 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read z deform equation.");
299 std::make_shared<LibUtilities::Equation>(
300 pSession->GetInterpreter(), zDeformStr);
304 indx, domFind, domain, coordDim, xDeformEqn, yDeformEqn,
311 WARNINGL0(
false,
"Zone type '" + zoneType +
312 "' is unsupported. Valid types are: 'Fixed', "
313 "'Rotate', 'Translate', or 'Prescribe'.")
317 zonesElement = zonesElement->NextSiblingElement();
323 ASSERTL0(interfacesTag,
"Unable to find INTERFACES tag in file.");
324 TiXmlElement *interfaceElement = interfacesTag->FirstChildElement();
326 while (interfaceElement)
329 "INTERFACE" == (std::string)interfaceElement->Value(),
330 "Only INTERFACE tags may be present inside the INTERFACES block.")
335 err = interfaceElement->QueryStringAttribute(
"NAME", &
name);
336 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read interface name.");
337 TiXmlElement *sideElement = interfaceElement->FirstChildElement();
343 std::vector<int> cnt;
348 "In INTERFACE NAME " +
name +
349 ", only two sides may be present in each INTERFACE block.")
352 err = sideElement->QueryIntAttribute(
"ID", &indx);
353 ASSERTL0(err == TIXML_SUCCESS,
"Unable to read interface ID.");
355 std::string boundaryStr;
357 sideElement->QueryStringAttribute(
"BOUNDARY", &boundaryStr);
361 if (boundaryErr == TIXML_SUCCESS)
363 indxStr =
ReadTag(boundaryStr);
368 auto sideElVal = sideElement->ValueStr();
369 if (sideElVal ==
"LEFT" || sideElVal ==
"L")
373 else if (sideElVal ==
"RIGHT" || sideElVal ==
"R")
380 sideElement->ValueStr() +
381 " is not a valid interface side for interface "
383 name +
". Please only use LEFT or RIGHT.")
386 interfaces[cnt[cnt.size() - 1]] =
388 indx, boundaryEdge));
390 sideElement = sideElement->NextSiblingElement();
393 ASSERTL0(std::accumulate(cnt.begin(), cnt.end(), 0) == 1,
394 "You must have only one LEFT and one RIGHT side"
395 " present in interface NAME " +
400 interfaces[0], interfaces[1]));
401 interfaceElement = interfaceElement->NextSiblingElement();
412 TiXmlElement *movement =
new TiXmlElement(
"MOVEMENT");
413 root->LinkEndChild(movement);
415 TiXmlElement *zones =
new TiXmlElement(
"ZONES");
421 TiXmlElement *e =
new TiXmlElement(
label.substr(0, 1));
422 e->SetAttribute(
"ID", i.first);
424 s <<
"D[" <<
z->GetDomainID() <<
"]";
425 e->SetAttribute(
"DOMAIN", s.str());
431 auto rotate = std::static_pointer_cast<ZoneRotate>(
z);
434 e->SetAttribute(
"AXIS",
436 e->SetAttribute(
"ANGVEL",
437 rotate->GetAngularVelEqn()->GetExpression());
442 auto translate = std::static_pointer_cast<ZoneTranslate>(
z);
443 const std::vector<NekDouble> vel = translate->GetVel();
444 std::stringstream vel_s;
445 vel_s << vel[0] <<
", " << vel[1] <<
", " << vel[2];
446 e->SetAttribute(
"VELOCITY", vel_s.str());
451 auto prescribe = std::static_pointer_cast<ZonePrescribe>(
z);
454 prescribe->GetXDeformEquation()->GetExpression());
457 prescribe->GetYDeformEquation()->GetExpression());
460 prescribe->GetZDeformEquation()->GetExpression());
466 zones->LinkEndChild(e);
468 movement->LinkEndChild(zones);
470 TiXmlElement *interfaces =
new TiXmlElement(
"INTERFACES");
473 const std::string interfaceName = i.first.second;
474 TiXmlElement *e =
new TiXmlElement(
"INTERFACE");
475 e->SetAttribute(
"NAME", interfaceName);
480 TiXmlElement *left_e =
new TiXmlElement(
"L");
481 left_e->SetAttribute(
"ID", left->GetId());
482 left_e->SetAttribute(
486 e->LinkEndChild(left_e);
490 TiXmlElement *right_e =
new TiXmlElement(
"R");
491 right_e->SetAttribute(
"ID", right->GetId());
492 right_e->SetAttribute(
496 e->LinkEndChild(right_e);
498 interfaces->LinkEndChild(e);
500 movement->LinkEndChild(interfaces);
507 std::set<int> movedZoneIds;
510 if (zone.second->Move(timeStep))
512 movedZoneIds.insert(zone.first);
519 int leftId = interPair.second->GetLeftInterface()->GetId();
520 int rightId = interPair.second->GetRightInterface()->GetId();
522 if (movedZoneIds.find(leftId) != movedZoneIds.end() ||
523 movedZoneIds.find(rightId) != movedZoneIds.end())
525 m_zones[leftId]->GetMoved() =
true;
526 m_zones[rightId]->GetMoved() =
true;
#define ASSERTL0(condition, msg)
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
#define WARNINGL0(condition, msg)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
static std::string GenerateSeqString(const std::vector< T > &v)
Generate a compressed comma-separated string representation of a vector of unsigned integers.
static bool GenerateVector(const std::string &str, std::vector< T > &out)
Takes a comma-separated string and converts it to entries in a vector.
Base class for a spectral/hp element mesh.
std::map< int, std::map< int, CompositeSharedPtr > > & GetDomain()
void GetCompositeList(const std::string &compositeStr, CompositeMap &compositeVector) const
int GetSpaceDimension()
Dimension of the space (can be a 1D curve in 3D space).
void PerformMovement(NekDouble timeStep)
Movement()
Default constructor.
void ReadZones(TiXmlElement *zonesTag, MeshGraph *meshGraph, const LibUtilities::SessionReaderSharedPtr &pSession)
Read zones given TiXmlDocument.
void AddInterface(std::string name, InterfaceShPtr left, InterfaceShPtr right)
Add pair of interfaces to this data.
void ReadInterfaces(TiXmlElement *interfacesTag, MeshGraph *meshGraph)
Read interfaces given TiXmlDocument.
void AddZone(ZoneBaseShPtr zone)
Add a zone object to this Movement data.
void WriteMovement(TiXmlElement *root)
Write the MOVEMENT section of the XML file.
InterfaceCollection m_interfaces
std::map< int, ZoneBaseShPtr > m_zones
const char *const ShapeTypeMap[SIZE_ShapeType]
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< Equation > EquationSharedPtr
const std::vector< NekDouble > velocity
const NekPoint< NekDouble > origin
static std::string StripParentheses(const std::string &str)
std::shared_ptr< ZoneBase > ZoneBaseShPtr
std::shared_ptr< ZonePrescribe > ZonePrescribeShPtr
std::shared_ptr< ZoneTranslate > ZoneTranslateShPtr
std::shared_ptr< InterfacePair > InterfacePairShPtr
std::shared_ptr< ZoneRotate > ZoneRotateShPtr
std::shared_ptr< Interface > InterfaceShPtr
MovementType
Enum of zone movement type.
std::shared_ptr< ZoneFixed > ZoneFixedShPtr
static std::string ReadTag(std::string &tagStr)
const std::string MovementTypeStr[]
Map of zone movement type to movement type string.
std::map< int, CompositeSharedPtr > CompositeMap
std::vector< double > z(NPUPPER)