Nektar++
MeshGraphIOXmlCompressed.cpp
Go to the documentation of this file.
1////////////////////////////////////////////////////////////////////////////////
2//
3// File: MeshGraphIOXmlCompressed.cpp
4//
5// For more information, please see: http://www.nektar.info/
6//
7// The MIT License
8//
9// Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10// Department of Aeronautics, Imperial College London (UK), and Scientific
11// Computing and Imaging Institute, University of Utah (USA).
12//
13// Permission is hereby granted, free of charge, to any person obtaining a
14// copy of this software and associated documentation files (the "Software"),
15// to deal in the Software without restriction, including without limitation
16// the rights to use, copy, modify, merge, publish, distribute, sublicense,
17// and/or sell copies of the Software, and to permit persons to whom the
18// Software is furnished to do so, subject to the following conditions:
19//
20// The above copyright notice and this permission notice shall be included
21// in all copies or substantial portions of the Software.
22//
23// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24// OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26// THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29// DEALINGS IN THE SOFTWARE.
30//
31// Description:
32//
33//
34////////////////////////////////////////////////////////////////////////////////
35
37
41
44
45// These are required for the Write(...) and Import(...) functions.
46#include <boost/algorithm/string/predicate.hpp>
47#include <boost/archive/iterators/base64_from_binary.hpp>
48#include <boost/archive/iterators/binary_from_base64.hpp>
49#include <boost/archive/iterators/transform_width.hpp>
50#include <boost/format.hpp>
51#include <boost/iostreams/copy.hpp>
52#include <boost/iostreams/filter/zlib.hpp>
53#include <boost/iostreams/filtering_stream.hpp>
54
55#include <tinyxml.h>
56
58{
59
62 "XmlCompressed", MeshGraphIOXmlCompressed::create,
63 "IO with compressed Xml geometry");
64
66{
67 PointGeomMap &vertSet = m_meshGraph->GetAllPointGeoms();
68 int spaceDimension = m_meshGraph->GetSpaceDimension();
69
70 // Now read the vertices
71 TiXmlElement *element = m_xmlGeom->FirstChildElement("VERTEX");
72 ASSERTL0(element, "Unable to find mesh VERTEX tag in file.");
73
74 NekDouble xscale, yscale, zscale;
75
76 // check to see if any scaling parameters are in
77 // attributes and determine these values
78 LibUtilities::Interpreter expEvaluator;
79 const char *xscal = element->Attribute("XSCALE");
80 if (!xscal)
81 {
82 xscale = 1.0;
83 }
84 else
85 {
86 std::string xscalstr = xscal;
87 int expr_id = expEvaluator.DefineFunction("", xscalstr);
88 xscale = expEvaluator.Evaluate(expr_id);
89 }
90
91 const char *yscal = element->Attribute("YSCALE");
92 if (!yscal)
93 {
94 yscale = 1.0;
95 }
96 else
97 {
98 std::string yscalstr = yscal;
99 int expr_id = expEvaluator.DefineFunction("", yscalstr);
100 yscale = expEvaluator.Evaluate(expr_id);
101 }
102
103 const char *zscal = element->Attribute("ZSCALE");
104 if (!zscal)
105 {
106 zscale = 1.0;
107 }
108 else
109 {
110 std::string zscalstr = zscal;
111 int expr_id = expEvaluator.DefineFunction("", zscalstr);
112 zscale = expEvaluator.Evaluate(expr_id);
113 }
114
115 NekDouble xmove, ymove, zmove;
116
117 // check to see if any moving parameters are in
118 // attributes and determine these values
119
120 const char *xmov = element->Attribute("XMOVE");
121 if (!xmov)
122 {
123 xmove = 0.0;
124 }
125 else
126 {
127 std::string xmovstr = xmov;
128 int expr_id = expEvaluator.DefineFunction("", xmovstr);
129 xmove = expEvaluator.Evaluate(expr_id);
130 }
131
132 const char *ymov = element->Attribute("YMOVE");
133 if (!ymov)
134 {
135 ymove = 0.0;
136 }
137 else
138 {
139 std::string ymovstr = ymov;
140 int expr_id = expEvaluator.DefineFunction("", ymovstr);
141 ymove = expEvaluator.Evaluate(expr_id);
142 }
143
144 const char *zmov = element->Attribute("ZMOVE");
145 if (!zmov)
146 {
147 zmove = 0.0;
148 }
149 else
150 {
151 std::string zmovstr = zmov;
152 int expr_id = expEvaluator.DefineFunction("", zmovstr);
153 zmove = expEvaluator.Evaluate(expr_id);
154 }
155
156 std::string IsCompressed;
157 element->QueryStringAttribute("COMPRESSED", &IsCompressed);
158
159 if (boost::iequals(IsCompressed,
161 {
162 // Extract the vertex body
163 TiXmlNode *vertexChild = element->FirstChild();
164 ASSERTL0(vertexChild, "Unable to extract the data from the compressed "
165 "vertex tag.");
166
167 std::string vertexStr;
168 if (vertexChild->Type() == TiXmlNode::TINYXML_TEXT)
169 {
170 vertexStr += vertexChild->ToText()->ValueStr();
171 }
172
173 std::vector<SpatialDomains::MeshVertex> vertData;
175 vertData);
176
177 int indx;
178 NekDouble xval, yval, zval;
179 for (int i = 0; i < vertData.size(); ++i)
180 {
181 indx = vertData[i].id;
182 xval = vertData[i].x;
183 yval = vertData[i].y;
184 zval = vertData[i].z;
185
186 xval = xval * xscale + xmove;
187 yval = yval * yscale + ymove;
188 zval = zval * zscale + zmove;
189
191 spaceDimension, indx, xval, yval, zval));
192
193 vert->SetGlobalID(indx);
194 vertSet[indx] = vert;
195 }
196 }
197 else
198 {
199 ASSERTL0(false, "Compressed formats do not match. Expected :" +
201 " but got " + IsCompressed);
202 }
203}
204
206{
207 auto &curvedEdges = m_meshGraph->GetCurvedEdges();
208 auto &curvedFaces = m_meshGraph->GetCurvedFaces();
209 int spaceDimension = m_meshGraph->GetSpaceDimension();
210
211 // check to see if any scaling parameters are in
212 // attributes and determine these values
213 TiXmlElement *element = m_xmlGeom->FirstChildElement("VERTEX");
214 ASSERTL0(element, "Unable to find mesh VERTEX tag in file.");
215
216 NekDouble xscale, yscale, zscale;
217
218 LibUtilities::Interpreter expEvaluator;
219 const char *xscal = element->Attribute("XSCALE");
220 if (!xscal)
221 {
222 xscale = 1.0;
223 }
224 else
225 {
226 std::string xscalstr = xscal;
227 int expr_id = expEvaluator.DefineFunction("", xscalstr);
228 xscale = expEvaluator.Evaluate(expr_id);
229 }
230
231 const char *yscal = element->Attribute("YSCALE");
232 if (!yscal)
233 {
234 yscale = 1.0;
235 }
236 else
237 {
238 std::string yscalstr = yscal;
239 int expr_id = expEvaluator.DefineFunction("", yscalstr);
240 yscale = expEvaluator.Evaluate(expr_id);
241 }
242
243 const char *zscal = element->Attribute("ZSCALE");
244 if (!zscal)
245 {
246 zscale = 1.0;
247 }
248 else
249 {
250 std::string zscalstr = zscal;
251 int expr_id = expEvaluator.DefineFunction("", zscalstr);
252 zscale = expEvaluator.Evaluate(expr_id);
253 }
254
255 NekDouble xmove, ymove, zmove;
256
257 // check to see if any moving parameters are in
258 // attributes and determine these values
259
260 const char *xmov = element->Attribute("XMOVE");
261 if (!xmov)
262 {
263 xmove = 0.0;
264 }
265 else
266 {
267 std::string xmovstr = xmov;
268 int expr_id = expEvaluator.DefineFunction("", xmovstr);
269 xmove = expEvaluator.Evaluate(expr_id);
270 }
271
272 const char *ymov = element->Attribute("YMOVE");
273 if (!ymov)
274 {
275 ymove = 0.0;
276 }
277 else
278 {
279 std::string ymovstr = ymov;
280 int expr_id = expEvaluator.DefineFunction("", ymovstr);
281 ymove = expEvaluator.Evaluate(expr_id);
282 }
283
284 const char *zmov = element->Attribute("ZMOVE");
285 if (!zmov)
286 {
287 zmove = 0.0;
288 }
289 else
290 {
291 std::string zmovstr = zmov;
292 int expr_id = expEvaluator.DefineFunction("", zmovstr);
293 zmove = expEvaluator.Evaluate(expr_id);
294 }
295
296 /// Look for elements in CURVE block.
297 TiXmlElement *field = m_xmlGeom->FirstChildElement("CURVED");
298
299 if (!field) // return if no curved entities
300 {
301 return;
302 }
303
304 std::string IsCompressed;
305 field->QueryStringAttribute("COMPRESSED", &IsCompressed);
306
307 if (IsCompressed.size() == 0)
308 {
309 // this could be that the curved tag is empty
310 // in this case we dont want to read it
311 return;
312 }
313
314 ASSERTL0(boost::iequals(IsCompressed,
316 "Compressed formats do not match. Expected :" +
318 IsCompressed);
319
320 std::vector<SpatialDomains::MeshCurvedInfo> edginfo;
321 std::vector<SpatialDomains::MeshCurvedInfo> facinfo;
323
324 // read edge, face info and curved poitns.
325 TiXmlElement *x = field->FirstChildElement();
326 while (x)
327 {
328 const char *entitytype = x->Value();
329 // read in edge or face info
330 if (boost::iequals(entitytype, "E"))
331 {
332 // read in data
333 std::string elmtStr;
334 TiXmlNode *child = x->FirstChild();
335
336 if (child->Type() == TiXmlNode::TINYXML_TEXT)
337 {
338 elmtStr += child->ToText()->ValueStr();
339 }
340
342 edginfo);
343 }
344 else if (boost::iequals(entitytype, "F"))
345 {
346 // read in data
347 std::string elmtStr;
348 TiXmlNode *child = x->FirstChild();
349
350 if (child->Type() == TiXmlNode::TINYXML_TEXT)
351 {
352 elmtStr += child->ToText()->ValueStr();
353 }
354
356 facinfo);
357 }
358 else if (boost::iequals(entitytype, "DATAPOINTS"))
359 {
360 int32_t id;
361 ASSERTL0(x->Attribute("ID", &id),
362 "Failed to get ID from PTS section");
363 cpts.id = id;
364
365 // read in data
366 std::string elmtStr;
367
368 TiXmlElement *DataIdx = x->FirstChildElement("INDEX");
369 ASSERTL0(DataIdx, "Cannot read data index tag in compressed "
370 "curved section");
371
372 TiXmlNode *child = DataIdx->FirstChild();
373 if (child->Type() == TiXmlNode::TINYXML_TEXT)
374 {
375 elmtStr = child->ToText()->ValueStr();
376 }
377
379 cpts.index);
380
381 TiXmlElement *DataPts = x->FirstChildElement("POINTS");
382 ASSERTL0(DataPts, "Cannot read data pts tag in compressed "
383 "curved section");
384
385 child = DataPts->FirstChild();
386 if (child->Type() == TiXmlNode::TINYXML_TEXT)
387 {
388 elmtStr = child->ToText()->ValueStr();
389 }
390
392 cpts.pts);
393 }
394 else
395 {
396 ASSERTL0(false, "Unknown tag in curved section");
397 }
398 x = x->NextSiblingElement();
399 }
400
401 // rescale (x,y,z) points;
402 for (int i = 0; i < cpts.pts.size(); ++i)
403 {
404 cpts.pts[i].x = xscale * cpts.pts[i].x + xmove;
405 cpts.pts[i].y = yscale * cpts.pts[i].y + ymove;
406 cpts.pts[i].z = zscale * cpts.pts[i].z + zmove;
407 }
408
409 for (int i = 0; i < edginfo.size(); ++i)
410 {
411 int edgeid = edginfo[i].entityid;
413
415 edgeid, ptype = (LibUtilities::PointsType)edginfo[i].ptype));
416
417 // load points
418 int offset = edginfo[i].ptoffset;
419 for (int j = 0; j < edginfo[i].npoints; ++j)
420 {
421 int idx = cpts.index[offset + j];
422
424 spaceDimension, edginfo[i].id, cpts.pts[idx].x, cpts.pts[idx].y,
425 cpts.pts[idx].z));
426 curve->m_points.push_back(vert);
427 }
428
429 curvedEdges[edgeid] = curve;
430 }
431
432 for (int i = 0; i < facinfo.size(); ++i)
433 {
434 int faceid = facinfo[i].entityid;
436
438 faceid, ptype = (LibUtilities::PointsType)facinfo[i].ptype));
439
440 int offset = facinfo[i].ptoffset;
441 for (int j = 0; j < facinfo[i].npoints; ++j)
442 {
443 int idx = cpts.index[offset + j];
444
446 spaceDimension, facinfo[i].id, cpts.pts[idx].x, cpts.pts[idx].y,
447 cpts.pts[idx].z));
448 curve->m_points.push_back(vert);
449 }
450
451 curvedFaces[faceid] = curve;
452 }
453}
454
456{
457 auto &segGeoms = m_meshGraph->GetAllSegGeoms();
458 auto &curvedEdges = m_meshGraph->GetCurvedEdges();
459 int spaceDimension = m_meshGraph->GetSpaceDimension();
460
461 CurveMap::iterator it;
462
463 /// Look for elements in ELEMENT block.
464 TiXmlElement *field = m_xmlGeom->FirstChildElement("EDGE");
465
466 ASSERTL0(field, "Unable to find EDGE tag in file.");
467
468 std::string IsCompressed;
469 field->QueryStringAttribute("COMPRESSED", &IsCompressed);
470
471 ASSERTL0(boost::iequals(IsCompressed,
473 "Compressed formats do not match. Expected :" +
475 IsCompressed);
476 // Extract the edge body
477 TiXmlNode *edgeChild = field->FirstChild();
478 ASSERTL0(edgeChild, "Unable to extract the data from "
479 "the compressed edge tag.");
480
481 std::string edgeStr;
482 if (edgeChild->Type() == TiXmlNode::TINYXML_TEXT)
483 {
484 edgeStr += edgeChild->ToText()->ValueStr();
485 }
486
487 std::vector<SpatialDomains::MeshEdge> edgeData;
489
490 int indx;
491 for (int i = 0; i < edgeData.size(); ++i)
492 {
493 indx = edgeData[i].id;
494 PointGeomSharedPtr vertices[2] = {
495 m_meshGraph->GetVertex(edgeData[i].v0),
496 m_meshGraph->GetVertex(edgeData[i].v1)};
497 SegGeomSharedPtr edge;
498
499 it = curvedEdges.find(indx);
500 if (it == curvedEdges.end())
501 {
503 indx, spaceDimension, vertices);
504 }
505 else
506 {
508 indx, spaceDimension, vertices, it->second);
509 }
510 segGeoms[indx] = edge;
511 }
512}
513
515{
516 auto &curvedFaces = m_meshGraph->GetCurvedFaces();
517 auto &triGeoms = m_meshGraph->GetAllTriGeoms();
518 auto &quadGeoms = m_meshGraph->GetAllQuadGeoms();
519
520 /// Look for elements in FACE block.
521 TiXmlElement *field = m_xmlGeom->FirstChildElement("FACE");
522
523 ASSERTL0(field, "Unable to find FACE tag in file.");
524
525 /// All faces are of the form: "<? ID="#"> ... </?>", with
526 /// ? being an element type (either Q or T).
527 /// They might be in compressed format and so then need upacking.
528
529 TiXmlElement *element = field->FirstChildElement();
530 CurveMap::iterator it;
531
532 while (element)
533 {
534 std::string elementType(element->ValueStr());
535
536 ASSERTL0(elementType == "Q" || elementType == "T",
537 (std::string("Unknown 3D face type: ") + elementType).c_str());
538
539 std::string IsCompressed;
540 element->QueryStringAttribute("COMPRESSED", &IsCompressed);
541
542 ASSERTL0(
543 boost::iequals(IsCompressed,
545 "Compressed formats do not match. Expected :" +
547 IsCompressed);
548
549 // Extract the face body
550 TiXmlNode *faceChild = element->FirstChild();
551 ASSERTL0(faceChild, "Unable to extract the data from "
552 "the compressed face tag.");
553
554 std::string faceStr;
555 if (faceChild->Type() == TiXmlNode::TINYXML_TEXT)
556 {
557 faceStr += faceChild->ToText()->ValueStr();
558 }
559
560 int indx;
561 if (elementType == "T")
562 {
563 std::vector<SpatialDomains::MeshTri> faceData;
565 faceData);
566
567 for (int i = 0; i < faceData.size(); ++i)
568 {
569 indx = faceData[i].id;
570
571 /// See if this face has curves.
572 it = curvedFaces.find(indx);
573
574 /// Create a TriGeom to hold the new definition.
576 m_meshGraph->GetSegGeom(faceData[i].e[0]),
577 m_meshGraph->GetSegGeom(faceData[i].e[1]),
578 m_meshGraph->GetSegGeom(faceData[i].e[2])};
579
580 TriGeomSharedPtr trigeom;
581 if (it == curvedFaces.end())
582 {
583 trigeom =
585 }
586 else
587 {
589 indx, edges, it->second);
590 }
591 trigeom->SetGlobalID(indx);
592 triGeoms[indx] = trigeom;
593 }
594 }
595 else if (elementType == "Q")
596 {
597 std::vector<SpatialDomains::MeshQuad> faceData;
599 faceData);
600
601 for (int i = 0; i < faceData.size(); ++i)
602 {
603 indx = faceData[i].id;
604
605 /// See if this face has curves.
606 it = curvedFaces.find(indx);
607
608 /// Create a QuadGeom to hold the new definition.
610 m_meshGraph->GetSegGeom(faceData[i].e[0]),
611 m_meshGraph->GetSegGeom(faceData[i].e[1]),
612 m_meshGraph->GetSegGeom(faceData[i].e[2]),
613 m_meshGraph->GetSegGeom(faceData[i].e[3])};
614
615 QuadGeomSharedPtr quadgeom;
616 if (it == curvedFaces.end())
617 {
618 quadgeom =
620 }
621 else
622 {
624 indx, edges, it->second);
625 }
626 quadgeom->SetGlobalID(indx);
627 quadGeoms[indx] = quadgeom;
628 }
629 }
630 /// Keep looking
631 element = element->NextSiblingElement();
632 }
633}
634
636{
637 auto &curvedEdges = m_meshGraph->GetCurvedEdges();
638 auto &segGeoms = m_meshGraph->GetAllSegGeoms();
639 int spaceDimension = m_meshGraph->GetSpaceDimension();
640
641 TiXmlElement *field = nullptr;
642
643 /// Look for elements in ELEMENT block.
644 field = m_xmlGeom->FirstChildElement("ELEMENT");
645
646 ASSERTL0(field, "Unable to find ELEMENT tag in file.");
647
648 /// All elements are of the form: "<S ID = n> ... </S>", with
649 /// ? being the element type.
650
651 TiXmlElement *segment = field->FirstChildElement("S");
652 CurveMap::iterator it;
653
654 while (segment)
655 {
656 std::string IsCompressed;
657 segment->QueryStringAttribute("COMPRESSED", &IsCompressed);
658 ASSERTL0(
659 boost::iequals(IsCompressed,
661 "Compressed formats do not match. Expected :" +
663 IsCompressed);
664
665 // Extract the face body
666 TiXmlNode *child = segment->FirstChild();
667 ASSERTL0(child, "Unable to extract the data from "
668 "the compressed face tag.");
669
670 std::string str;
671 if (child->Type() == TiXmlNode::TINYXML_TEXT)
672 {
673 str += child->ToText()->ValueStr();
674 }
675
676 int indx;
677
678 std::vector<SpatialDomains::MeshEdge> data;
680
681 for (int i = 0; i < data.size(); ++i)
682 {
683 indx = data[i].id;
684
685 /// See if this face has curves.
686 it = curvedEdges.find(indx);
687
688 PointGeomSharedPtr vertices[2] = {
689 m_meshGraph->GetVertex(data[i].v0),
690 m_meshGraph->GetVertex(data[i].v1)};
692
693 if (it == curvedEdges.end())
694 {
696 indx, spaceDimension, vertices);
697 seg->SetGlobalID(indx); // Set global mesh id
698 }
699 else
700 {
702 indx, spaceDimension, vertices, it->second);
703 seg->SetGlobalID(indx); // Set global mesh id
704 }
705 seg->SetGlobalID(indx);
706 segGeoms[indx] = seg;
707 }
708 /// Keep looking for additional segments
709 segment = segment->NextSiblingElement("S");
710 }
711}
712
714{
715 auto &curvedFaces = m_meshGraph->GetCurvedFaces();
716 auto &triGeoms = m_meshGraph->GetAllTriGeoms();
717 auto &quadGeoms = m_meshGraph->GetAllQuadGeoms();
718
719 /// Look for elements in ELEMENT block.
720 TiXmlElement *field = m_xmlGeom->FirstChildElement("ELEMENT");
721
722 ASSERTL0(field, "Unable to find ELEMENT tag in file.");
723
724 // Set up curve map for curved elements on an embedded manifold.
725 CurveMap::iterator it;
726
727 /// All elements are of the form: "<? ID="#"> ... </?>", with
728 /// ? being the element type.
729
730 TiXmlElement *element = field->FirstChildElement();
731
732 while (element)
733 {
734 std::string elementType(element->ValueStr());
735
736 ASSERTL0(
737 elementType == "Q" || elementType == "T",
738 (std::string("Unknown 2D element type: ") + elementType).c_str());
739
740 std::string IsCompressed;
741 element->QueryStringAttribute("COMPRESSED", &IsCompressed);
742
743 ASSERTL0(
744 boost::iequals(IsCompressed,
746 "Compressed formats do not match. Expected :" +
748 IsCompressed);
749
750 // Extract the face body
751 TiXmlNode *faceChild = element->FirstChild();
752 ASSERTL0(faceChild, "Unable to extract the data from "
753 "the compressed face tag.");
754
755 std::string faceStr;
756 if (faceChild->Type() == TiXmlNode::TINYXML_TEXT)
757 {
758 faceStr += faceChild->ToText()->ValueStr();
759 }
760
761 int indx;
762 if (elementType == "T")
763 {
764 std::vector<SpatialDomains::MeshTri> faceData;
766 faceData);
767
768 for (int i = 0; i < faceData.size(); ++i)
769 {
770 indx = faceData[i].id;
771
772 /// See if this face has curves.
773 it = curvedFaces.find(indx);
774
775 /// Create a TriGeom to hold the new definition.
777 m_meshGraph->GetSegGeom(faceData[i].e[0]),
778 m_meshGraph->GetSegGeom(faceData[i].e[1]),
779 m_meshGraph->GetSegGeom(faceData[i].e[2])};
780
781 TriGeomSharedPtr trigeom;
782 if (it == curvedFaces.end())
783 {
784 trigeom =
786 }
787 else
788 {
790 indx, edges, it->second);
791 }
792 trigeom->SetGlobalID(indx);
793 triGeoms[indx] = trigeom;
794 }
795 }
796 else if (elementType == "Q")
797 {
798 std::vector<SpatialDomains::MeshQuad> faceData;
800 faceData);
801
802 for (int i = 0; i < faceData.size(); ++i)
803 {
804 indx = faceData[i].id;
805
806 /// See if this face has curves.
807 it = curvedFaces.find(indx);
808
809 /// Create a QuadGeom to hold the new definition.
811 m_meshGraph->GetSegGeom(faceData[i].e[0]),
812 m_meshGraph->GetSegGeom(faceData[i].e[1]),
813 m_meshGraph->GetSegGeom(faceData[i].e[2]),
814 m_meshGraph->GetSegGeom(faceData[i].e[3])};
815
816 QuadGeomSharedPtr quadgeom;
817 if (it == curvedFaces.end())
818 {
819 quadgeom =
821 }
822 else
823 {
825 indx, edges, it->second);
826 }
827 quadgeom->SetGlobalID(indx);
828 quadGeoms[indx] = quadgeom;
829 }
830 }
831 /// Keep looking
832 element = element->NextSiblingElement();
833 }
834}
835
837{
838 auto &tetGeoms = m_meshGraph->GetAllTetGeoms();
839 auto &pyrGeoms = m_meshGraph->GetAllPyrGeoms();
840 auto &prismGeoms = m_meshGraph->GetAllPrismGeoms();
841 auto &hexGeoms = m_meshGraph->GetAllHexGeoms();
842
843 /// Look for elements in ELEMENT block.
844 TiXmlElement *field = m_xmlGeom->FirstChildElement("ELEMENT");
845
846 ASSERTL0(field, "Unable to find ELEMENT tag in file.");
847
848 /// All elements are of the form: "<? ID="#"> ... </?>", with
849 /// ? being the element type.
850
851 TiXmlElement *element = field->FirstChildElement();
852
853 while (element)
854 {
855 std::string elementType(element->ValueStr());
856
857 // A - tet, P - pyramid, R - prism, H - hex
858 ASSERTL0(
859 elementType == "A" || elementType == "P" || elementType == "R" ||
860 elementType == "H",
861 (std::string("Unknown 3D element type: ") + elementType).c_str());
862
863 std::string IsCompressed;
864 element->QueryStringAttribute("COMPRESSED", &IsCompressed);
865
866 ASSERTL0(
867 boost::iequals(IsCompressed,
869 "Compressed formats do not match. Expected :" +
871 IsCompressed);
872
873 // Extract the face body
874 TiXmlNode *child = element->FirstChild();
875 ASSERTL0(child, "Unable to extract the data from "
876 "the compressed face tag.");
877
878 std::string str;
879 if (child->Type() == TiXmlNode::TINYXML_TEXT)
880 {
881 str += child->ToText()->ValueStr();
882 }
883
884 int indx;
885 if (elementType == "A")
886 {
887 std::vector<SpatialDomains::MeshTet> data;
889 TriGeomSharedPtr tfaces[4];
890 for (int i = 0; i < data.size(); ++i)
891 {
892 indx = data[i].id;
893 for (int j = 0; j < 4; ++j)
894 {
896 m_meshGraph->GetGeometry2D(data[i].f[j]);
897 tfaces[j] = std::static_pointer_cast<TriGeom>(face);
898 }
899
900 TetGeomSharedPtr tetgeom(
902 tetGeoms[indx] = tetgeom;
903 m_meshGraph->PopulateFaceToElMap(tetgeom, 4);
904 }
905 }
906 else if (elementType == "P")
907 {
908 std::vector<SpatialDomains::MeshPyr> data;
910 Geometry2DSharedPtr faces[5];
911 for (int i = 0; i < data.size(); ++i)
912 {
913 indx = data[i].id;
914 int Ntfaces = 0;
915 int Nqfaces = 0;
916 for (int j = 0; j < 5; ++j)
917 {
919 m_meshGraph->GetGeometry2D(data[i].f[j]);
920
921 if (face == Geometry2DSharedPtr() ||
922 (face->GetShapeType() != LibUtilities::eTriangle &&
923 face->GetShapeType() != LibUtilities::eQuadrilateral))
924 {
925 std::stringstream errorstring;
926 errorstring << "Element " << indx
927 << " has invalid face: " << j;
928 ASSERTL0(false, errorstring.str().c_str());
929 }
930 else if (face->GetShapeType() == LibUtilities::eTriangle)
931 {
932 faces[j] = std::static_pointer_cast<TriGeom>(face);
933 Ntfaces++;
934 }
935 else if (face->GetShapeType() ==
937 {
938 faces[j] = std::static_pointer_cast<QuadGeom>(face);
939 Nqfaces++;
940 }
941 }
942 ASSERTL0((Ntfaces == 4) && (Nqfaces == 1),
943 "Did not identify the correct number of "
944 "triangular and quadrilateral faces for a "
945 "pyramid");
946
947 PyrGeomSharedPtr pyrgeom(
949
950 pyrGeoms[indx] = pyrgeom;
951 m_meshGraph->PopulateFaceToElMap(pyrgeom, 5);
952 }
953 }
954 else if (elementType == "R")
955 {
956 std::vector<SpatialDomains::MeshPrism> data;
958 Geometry2DSharedPtr faces[5];
959 for (int i = 0; i < data.size(); ++i)
960 {
961 indx = data[i].id;
962 int Ntfaces = 0;
963 int Nqfaces = 0;
964 for (int j = 0; j < 5; ++j)
965 {
967 m_meshGraph->GetGeometry2D(data[i].f[j]);
968 if (face == Geometry2DSharedPtr() ||
969 (face->GetShapeType() != LibUtilities::eTriangle &&
970 face->GetShapeType() != LibUtilities::eQuadrilateral))
971 {
972 std::stringstream errorstring;
973 errorstring << "Element " << indx
974 << " has invalid face: " << j;
975 ASSERTL0(false, errorstring.str().c_str());
976 }
977 else if (face->GetShapeType() == LibUtilities::eTriangle)
978 {
979 faces[j] = std::static_pointer_cast<TriGeom>(face);
980 Ntfaces++;
981 }
982 else if (face->GetShapeType() ==
984 {
985 faces[j] = std::static_pointer_cast<QuadGeom>(face);
986 Nqfaces++;
987 }
988 }
989 ASSERTL0((Ntfaces == 2) && (Nqfaces == 3),
990 "Did not identify the correct number of "
991 "triangular and quadrilateral faces for a "
992 "prism");
993
994 PrismGeomSharedPtr prismgeom(
996
997 prismGeoms[indx] = prismgeom;
998 m_meshGraph->PopulateFaceToElMap(prismgeom, 5);
999 }
1000 }
1001 else if (elementType == "H")
1002 {
1003 std::vector<SpatialDomains::MeshHex> data;
1005
1006 QuadGeomSharedPtr faces[6];
1007 for (int i = 0; i < data.size(); ++i)
1008 {
1009 indx = data[i].id;
1010 for (int j = 0; j < 6; ++j)
1011 {
1012 Geometry2DSharedPtr face =
1013 m_meshGraph->GetGeometry2D(data[i].f[j]);
1014 faces[j] = std::static_pointer_cast<QuadGeom>(face);
1015 }
1016
1017 HexGeomSharedPtr hexgeom(
1019 hexGeoms[indx] = hexgeom;
1020 m_meshGraph->PopulateFaceToElMap(hexgeom, 6);
1021 }
1022 }
1023 /// Keep looking
1024 element = element->NextSiblingElement();
1025 }
1026}
1027
1029 PointGeomMap &verts)
1030{
1031 if (verts.size() == 0)
1032 {
1033 return;
1034 }
1035
1036 TiXmlElement *vertTag = new TiXmlElement("VERTEX");
1037
1038 std::vector<MeshVertex> vertInfo;
1039
1040 for (auto &i : verts)
1041 {
1042 MeshVertex v;
1043 v.id = i.first;
1044 v.x = i.second->x();
1045 v.y = i.second->y();
1046 v.z = i.second->z();
1047 vertInfo.push_back(v);
1048 }
1049
1050 vertTag->SetAttribute("COMPRESSED",
1052 vertTag->SetAttribute("BITSIZE",
1054
1055 std::string vertStr;
1057
1058 vertTag->LinkEndChild(new TiXmlText(vertStr));
1059
1060 geomTag->LinkEndChild(vertTag);
1061}
1062
1064 SegGeomMap &edges)
1065{
1066 int meshDimension = m_meshGraph->GetMeshDimension();
1067
1068 if (edges.size() == 0)
1069 {
1070 return;
1071 }
1072
1073 TiXmlElement *edgeTag = new TiXmlElement(meshDimension == 1 ? "S" : "EDGE");
1074
1075 std::vector<MeshEdge> edgeInfo;
1076
1077 for (auto &i : edges)
1078 {
1079 MeshEdge e;
1080 e.id = i.first;
1081 e.v0 = i.second->GetVid(0);
1082 e.v1 = i.second->GetVid(1);
1083 edgeInfo.push_back(e);
1084 }
1085
1086 std::string edgeStr;
1088
1089 edgeTag->SetAttribute("COMPRESSED",
1091 edgeTag->SetAttribute("BITSIZE",
1093
1094 edgeTag->LinkEndChild(new TiXmlText(edgeStr));
1095
1096 if (meshDimension == 1)
1097 {
1098 TiXmlElement *tmp = new TiXmlElement("ELEMENT");
1099 tmp->LinkEndChild(edgeTag);
1100 geomTag->LinkEndChild(tmp);
1101 }
1102 else
1103 {
1104 geomTag->LinkEndChild(edgeTag);
1105 }
1106}
1107
1108void MeshGraphIOXmlCompressed::v_WriteTris(TiXmlElement *faceTag,
1109 TriGeomMap &tris)
1110{
1111 if (tris.size() == 0)
1112 {
1113 return;
1114 }
1115
1116 std::string tag = "T";
1117
1118 std::vector<MeshTri> triInfo;
1119
1120 for (auto &i : tris)
1121 {
1122 MeshTri t;
1123 t.id = i.first;
1124 t.e[0] = i.second->GetEid(0);
1125 t.e[1] = i.second->GetEid(1);
1126 t.e[2] = i.second->GetEid(2);
1127 triInfo.push_back(t);
1128 }
1129
1130 TiXmlElement *x = new TiXmlElement(tag);
1131 std::string triStr;
1133
1134 x->SetAttribute("COMPRESSED",
1136 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1137
1138 x->LinkEndChild(new TiXmlText(triStr));
1139
1140 faceTag->LinkEndChild(x);
1141}
1142
1144 QuadGeomMap &quads)
1145{
1146 if (quads.size() == 0)
1147 {
1148 return;
1149 }
1150
1151 std::string tag = "Q";
1152
1153 std::vector<MeshQuad> quadInfo;
1154
1155 for (auto &i : quads)
1156 {
1157 MeshQuad q;
1158 q.id = i.first;
1159 q.e[0] = i.second->GetEid(0);
1160 q.e[1] = i.second->GetEid(1);
1161 q.e[2] = i.second->GetEid(2);
1162 q.e[3] = i.second->GetEid(3);
1163 quadInfo.push_back(q);
1164 }
1165
1166 TiXmlElement *x = new TiXmlElement(tag);
1167 std::string quadStr;
1169
1170 x->SetAttribute("COMPRESSED",
1172 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1173
1174 x->LinkEndChild(new TiXmlText(quadStr));
1175
1176 faceTag->LinkEndChild(x);
1177}
1178
1179void MeshGraphIOXmlCompressed::v_WriteHexs(TiXmlElement *elmtTag,
1180 HexGeomMap &hexs)
1181{
1182 if (hexs.size() == 0)
1183 {
1184 return;
1185 }
1186
1187 std::string tag = "H";
1188
1189 std::vector<MeshHex> elementInfo;
1190
1191 for (auto &i : hexs)
1192 {
1193 MeshHex e;
1194 e.id = i.first;
1195 e.f[0] = i.second->GetFid(0);
1196 e.f[1] = i.second->GetFid(1);
1197 e.f[2] = i.second->GetFid(2);
1198 e.f[3] = i.second->GetFid(3);
1199 e.f[4] = i.second->GetFid(4);
1200 e.f[5] = i.second->GetFid(5);
1201 elementInfo.push_back(e);
1202 }
1203
1204 TiXmlElement *x = new TiXmlElement(tag);
1205 std::string elStr;
1207
1208 x->SetAttribute("COMPRESSED",
1210 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1211
1212 x->LinkEndChild(new TiXmlText(elStr));
1213
1214 elmtTag->LinkEndChild(x);
1215}
1216
1218 PrismGeomMap &pris)
1219{
1220 if (pris.size() == 0)
1221 {
1222 return;
1223 }
1224
1225 std::string tag = "R";
1226
1227 std::vector<MeshPrism> elementInfo;
1228
1229 for (auto &i : pris)
1230 {
1231 MeshPrism e;
1232 e.id = i.first;
1233 e.f[0] = i.second->GetFid(0);
1234 e.f[1] = i.second->GetFid(1);
1235 e.f[2] = i.second->GetFid(2);
1236 e.f[3] = i.second->GetFid(3);
1237 e.f[4] = i.second->GetFid(4);
1238 elementInfo.push_back(e);
1239 }
1240
1241 TiXmlElement *x = new TiXmlElement(tag);
1242 std::string elStr;
1244
1245 x->SetAttribute("COMPRESSED",
1247 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1248
1249 x->LinkEndChild(new TiXmlText(elStr));
1250
1251 elmtTag->LinkEndChild(x);
1252}
1253
1254void MeshGraphIOXmlCompressed::v_WritePyrs(TiXmlElement *elmtTag,
1255 PyrGeomMap &pyrs)
1256{
1257 if (pyrs.size() == 0)
1258 {
1259 return;
1260 }
1261
1262 std::string tag = "P";
1263
1264 std::vector<MeshPyr> elementInfo;
1265
1266 for (auto &i : pyrs)
1267 {
1268 MeshPyr e;
1269 e.id = i.first;
1270 e.f[0] = i.second->GetFid(0);
1271 e.f[1] = i.second->GetFid(1);
1272 e.f[2] = i.second->GetFid(2);
1273 e.f[3] = i.second->GetFid(3);
1274 e.f[4] = i.second->GetFid(4);
1275 elementInfo.push_back(e);
1276 }
1277
1278 TiXmlElement *x = new TiXmlElement(tag);
1279 std::string elStr;
1281
1282 x->SetAttribute("COMPRESSED",
1284 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1285
1286 x->LinkEndChild(new TiXmlText(elStr));
1287
1288 elmtTag->LinkEndChild(x);
1289}
1290
1291void MeshGraphIOXmlCompressed::v_WriteTets(TiXmlElement *elmtTag,
1292 TetGeomMap &tets)
1293{
1294 if (tets.size() == 0)
1295 {
1296 return;
1297 }
1298
1299 std::string tag = "A";
1300
1301 std::vector<MeshTet> elementInfo;
1302
1303 for (auto &i : tets)
1304 {
1305 MeshTet e;
1306 e.id = i.first;
1307 e.f[0] = i.second->GetFid(0);
1308 e.f[1] = i.second->GetFid(1);
1309 e.f[2] = i.second->GetFid(2);
1310 e.f[3] = i.second->GetFid(3);
1311 elementInfo.push_back(e);
1312 }
1313
1314 TiXmlElement *x = new TiXmlElement(tag);
1315 std::string elStr;
1317
1318 x->SetAttribute("COMPRESSED",
1320 x->SetAttribute("BITSIZE", LibUtilities::CompressData::GetBitSizeStr());
1321
1322 x->LinkEndChild(new TiXmlText(elStr));
1323
1324 elmtTag->LinkEndChild(x);
1325}
1326
1328 CurveMap &edges, CurveMap &faces)
1329{
1330 if (edges.size() == 0 && faces.size() == 0)
1331 {
1332 return;
1333 }
1334
1335 TiXmlElement *curveTag = new TiXmlElement("CURVED");
1336
1337 std::vector<MeshCurvedInfo> edgeInfo;
1338 std::vector<MeshCurvedInfo> faceInfo;
1339 MeshCurvedPts curvedPts;
1340 curvedPts.id = 0;
1341 int ptOffset = 0;
1342 int newIdx = 0;
1343 int edgeCnt = 0;
1344 int faceCnt = 0;
1345
1346 for (auto &i : edges)
1347 {
1348 MeshCurvedInfo cinfo;
1349 cinfo.id = edgeCnt++;
1350 cinfo.entityid = i.first;
1351 cinfo.npoints = i.second->m_points.size();
1352 cinfo.ptype = i.second->m_ptype;
1353 cinfo.ptid = 0;
1354 cinfo.ptoffset = ptOffset;
1355
1356 edgeInfo.push_back(cinfo);
1357
1358 for (int j = 0; j < i.second->m_points.size(); j++)
1359 {
1360 MeshVertex v;
1361 v.id = newIdx;
1362 v.x = i.second->m_points[j]->x();
1363 v.y = i.second->m_points[j]->y();
1364 v.z = i.second->m_points[j]->z();
1365 curvedPts.pts.push_back(v);
1366 curvedPts.index.push_back(newIdx);
1367 newIdx++;
1368 }
1369 ptOffset += cinfo.npoints;
1370 }
1371
1372 for (auto &i : faces)
1373 {
1374 MeshCurvedInfo cinfo;
1375 cinfo.id = faceCnt++;
1376 cinfo.entityid = i.first;
1377 cinfo.npoints = i.second->m_points.size();
1378 cinfo.ptype = i.second->m_ptype;
1379 cinfo.ptid = 0;
1380 cinfo.ptoffset = ptOffset;
1381
1382 faceInfo.push_back(cinfo);
1383
1384 for (int j = 0; j < i.second->m_points.size(); j++)
1385 {
1386 MeshVertex v;
1387 v.id = newIdx;
1388 v.x = i.second->m_points[j]->x();
1389 v.y = i.second->m_points[j]->y();
1390 v.z = i.second->m_points[j]->z();
1391 curvedPts.pts.push_back(v);
1392 curvedPts.index.push_back(newIdx);
1393 newIdx++;
1394 }
1395 ptOffset += cinfo.npoints;
1396 }
1397
1398 curveTag->SetAttribute("COMPRESSED",
1400 curveTag->SetAttribute("BITSIZE",
1402
1403 if (edgeInfo.size())
1404 {
1405 TiXmlElement *x = new TiXmlElement("E");
1406 std::string dataStr;
1408
1409 x->LinkEndChild(new TiXmlText(dataStr));
1410 curveTag->LinkEndChild(x);
1411 }
1412
1413 if (faceInfo.size())
1414 {
1415 TiXmlElement *x = new TiXmlElement("F");
1416 std::string dataStr;
1418
1419 x->LinkEndChild(new TiXmlText(dataStr));
1420 curveTag->LinkEndChild(x);
1421 }
1422
1423 if (edgeInfo.size() || faceInfo.size())
1424 {
1425 TiXmlElement *x = new TiXmlElement("DATAPOINTS");
1426 x->SetAttribute("ID", curvedPts.id);
1427 TiXmlElement *subx = new TiXmlElement("INDEX");
1428 std::string dataStr;
1430 dataStr);
1431 subx->LinkEndChild(new TiXmlText(dataStr));
1432 x->LinkEndChild(subx);
1433
1434 subx = new TiXmlElement("POINTS");
1436 dataStr);
1437 subx->LinkEndChild(new TiXmlText(dataStr));
1438 x->LinkEndChild(subx);
1439 curveTag->LinkEndChild(x);
1440 }
1441
1442 geomTag->LinkEndChild(curveTag);
1443}
1444} // namespace Nektar::SpatialDomains
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
Interpreter class for the evaluation of mathematical expressions.
Definition: Interpreter.h:76
int DefineFunction(const std::string &vlist, const std::string &expr)
Defines a function for the purposes of evaluation.
NekDouble Evaluate(const int id)
Evaluate a function which depends only on constants and/or parameters.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
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.
void v_WriteEdges(TiXmlElement *geomTag, SegGeomMap &edges) override
void v_WriteHexs(TiXmlElement *elmtTag, HexGeomMap &hexs) override
void v_WriteVertices(TiXmlElement *geomTag, PointGeomMap &verts) override
void v_WritePyrs(TiXmlElement *elmtTag, PyrGeomMap &pyrs) override
void v_WriteTets(TiXmlElement *elmtTag, TetGeomMap &tets) override
void v_WriteCurves(TiXmlElement *geomTag, CurveMap &edges, CurveMap &faces) override
void v_WritePrisms(TiXmlElement *elmtTag, PrismGeomMap &pris) override
void v_WriteTris(TiXmlElement *faceTag, TriGeomMap &tris) override
void v_WriteQuads(TiXmlElement *faceTag, QuadGeomMap &quads) override
static const int kNedges
Definition: QuadGeom.h:75
static const int kNedges
Get the orientation of face1.
Definition: TriGeom.h:70
int ZlibDecodeFromBase64Str(std::string &in64, std::vector< T > &out)
Definition: CompressData.h:216
int ZlibEncodeToBase64Str(std::vector< T > &in, std::string &out64)
Definition: CompressData.h:129
std::map< int, TriGeomSharedPtr > TriGeomMap
Definition: TriGeom.h:57
std::shared_ptr< QuadGeom > QuadGeomSharedPtr
Definition: HexGeom.h:45
std::map< int, PyrGeomSharedPtr > PyrGeomMap
Definition: PyrGeom.h:76
std::map< int, QuadGeomSharedPtr > QuadGeomMap
Definition: QuadGeom.h:53
std::shared_ptr< PrismGeom > PrismGeomSharedPtr
Definition: PrismGeom.h:82
std::map< int, SegGeomSharedPtr > SegGeomMap
Definition: SegGeom.h:50
std::shared_ptr< Curve > CurveSharedPtr
Definition: Curve.hpp:58
std::unordered_map< int, CurveSharedPtr > CurveMap
Definition: Curve.hpp:59
std::shared_ptr< HexGeom > HexGeomSharedPtr
Definition: HexGeom.h:84
std::shared_ptr< SegGeom > SegGeomSharedPtr
Definition: Geometry2D.h:58
std::map< int, TetGeomSharedPtr > TetGeomMap
Definition: TetGeom.h:84
std::shared_ptr< PyrGeom > PyrGeomSharedPtr
Definition: PyrGeom.h:75
MeshGraphIOFactory & GetMeshGraphIOFactory()
Definition: MeshGraphIO.cpp:45
std::shared_ptr< TetGeom > TetGeomSharedPtr
Definition: TetGeom.h:82
std::map< int, PrismGeomSharedPtr > PrismGeomMap
Definition: PrismGeom.h:83
std::shared_ptr< PointGeom > PointGeomSharedPtr
Definition: Geometry.h:57
std::shared_ptr< Geometry2D > Geometry2DSharedPtr
Definition: Geometry.h:62
std::shared_ptr< TriGeom > TriGeomSharedPtr
Definition: TriGeom.h:56
std::map< int, HexGeomSharedPtr > HexGeomMap
Definition: HexGeom.h:85
std::map< int, PointGeomSharedPtr > PointGeomMap
Definition: PointGeom.h:51
std::vector< double > q(NPUPPER *NPUPPER)
std::int32_t int32_t
double NekDouble
int64_t ptype
point offset of data entry for this curve
int64_t npoints
The entity id corresponding to the global edge/curve.
int64_t ptoffset
the id of point data map (currently always 0 since we are using just one set).
int64_t entityid
Id of this curved information.
int64_t ptid
The number of points in this curved entity.
std::vector< int64_t > index
Mapping to access the pts value. Given a 'ptoffset' value the npoints subsquent values provide the in...
std::vector< MeshVertex > pts
mapping to access pts value.