FreeFOAM The Cross-Platform CFD Toolkit
ensightMesh.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include <OpenFOAM/argList.H>
27 #include <OpenFOAM/Time.H>
28 #include "ensightMesh.H"
29 #include <finiteVolume/fvMesh.H>
33 #include <OpenFOAM/cellModeller.H>
34 #include <OpenFOAM/IOmanip.H>
35 #include "itoa.H"
36 #include "ensightWriteBinary.H"
37 #include <fstream>
38 
39 // * * * * * * * * * * * * * Private Functions * * * * * * * * * * * * * * //
40 
41 namespace Foam
42 {
43  //- Proxy-class to hold the patch processor list combination operator
44  class concatPatchProcs
45  {
46 
47  public:
48 
49  void operator()
50  (
51  List<labelList>& x,
52  const List<labelList>& y
53  ) const
54  {
55  forAll(y, i)
56  {
57  const labelList& yPatches = y[i];
58 
59  if (yPatches.size())
60  {
61  labelList& xPatches = x[i];
62 
63  label offset = xPatches.size();
64  xPatches.setSize(offset + yPatches.size());
65 
66  forAll(yPatches, i)
67  {
68  xPatches[i + offset] = yPatches[i];
69  }
70  }
71  }
72  }
73  };
74 } // End namespace Foam
75 
76 
77 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
78 
79 Foam::ensightMesh::ensightMesh
80 (
81  const fvMesh& mesh,
82  const argList& args,
83  const bool binary
84 )
85 :
86  mesh_(mesh),
87  binary_(binary),
88  patchPartOffset_(2),
89  meshCellSets_(mesh_.nCells()),
90  boundaryFaceSets_(mesh_.boundary().size()),
91  allPatchNames_(0),
92  allPatchProcs_(0),
93  patchNames_(0),
94  nPatchPrims_(0)
95 {
96  const cellShapeList& cellShapes = mesh.cellShapes();
97 
98  const cellModel& tet = *(cellModeller::lookup("tet"));
99  const cellModel& pyr = *(cellModeller::lookup("pyr"));
100  const cellModel& prism = *(cellModeller::lookup("prism"));
101  const cellModel& wedge = *(cellModeller::lookup("wedge"));
102  const cellModel& hex = *(cellModeller::lookup("hex"));
103 
104  if (!args.optionFound("noPatches"))
105  {
106  allPatchNames_ = wordList::subList
107  (
108  mesh_.boundaryMesh().names(), mesh_.boundary().size()
109  - mesh_.globalData().processorPatches().size()
110  );
111 
112  allPatchProcs_.setSize(allPatchNames_.size());
113 
114  forAll (allPatchProcs_, patchi)
115  {
116  if (mesh_.boundary()[patchi].size())
117  {
118  allPatchProcs_[patchi].setSize(1);
119  allPatchProcs_[patchi][0] = Pstream::myProcNo();
120  }
121  }
122 
123  combineReduce(allPatchProcs_, concatPatchProcs());
124 
125  if (args.optionFound("patches"))
126  {
127  wordList patchNameList(args.optionLookup("patches")());
128 
129  if (patchNameList.empty())
130  {
131  patchNameList = allPatchNames_;
132  }
133 
134  forAll (patchNameList, i)
135  {
136  patchNames_.insert(patchNameList[i]);
137  }
138  }
139  }
140 
141  if (patchNames_.size())
142  {
143  // no internalMesh
144  patchPartOffset_ = 1;
145  }
146  else
147  {
148  // Count the shapes
149  labelList& tets = meshCellSets_.tets;
150  labelList& pyrs = meshCellSets_.pyrs;
151  labelList& prisms = meshCellSets_.prisms;
152  labelList& wedges = meshCellSets_.wedges;
153  labelList& hexes = meshCellSets_.hexes;
154  labelList& polys = meshCellSets_.polys;
155 
156  label nTets = 0;
157  label nPyrs = 0;
158  label nPrisms = 0;
159  label nWedges = 0;
160  label nHexes = 0;
161  label nPolys = 0;
162 
163  forAll(cellShapes, cellI)
164  {
165  const cellShape& cellShape = cellShapes[cellI];
166  const cellModel& cellModel = cellShape.model();
167 
168  if (cellModel == tet)
169  {
170  tets[nTets++] = cellI;
171  }
172  else if (cellModel == pyr)
173  {
174  pyrs[nPyrs++] = cellI;
175  }
176  else if (cellModel == prism)
177  {
178  prisms[nPrisms++] = cellI;
179  }
180  else if (cellModel == wedge)
181  {
182  wedges[nWedges++] = cellI;
183  }
184  else if (cellModel == hex)
185  {
186  hexes[nHexes++] = cellI;
187  }
188  else
189  {
190  polys[nPolys++] = cellI;
191  }
192  }
193 
194  tets.setSize(nTets);
195  pyrs.setSize(nPyrs);
196  prisms.setSize(nPrisms);
197  wedges.setSize(nWedges);
198  hexes.setSize(nHexes);
199  polys.setSize(nPolys);
200 
201  meshCellSets_.nTets = nTets;
202  reduce(meshCellSets_.nTets, sumOp<label>());
203 
204  meshCellSets_.nPyrs = nPyrs;
205  reduce(meshCellSets_.nPyrs, sumOp<label>());
206 
207  meshCellSets_.nPrisms = nPrisms;
208  reduce(meshCellSets_.nPrisms, sumOp<label>());
209 
210  meshCellSets_.nHexesWedges = nHexes + nWedges;
211  reduce(meshCellSets_.nHexesWedges, sumOp<label>());
212 
213  meshCellSets_.nPolys = nPolys;
214  reduce(meshCellSets_.nPolys, sumOp<label>());
215  }
216 
217  if (!args.optionFound("noPatches"))
218  {
219  forAll (mesh.boundary(), patchi)
220  {
221  if (mesh.boundary()[patchi].size())
222  {
223  const polyPatch& p = mesh.boundaryMesh()[patchi];
224 
225  labelList& tris = boundaryFaceSets_[patchi].tris;
226  labelList& quads = boundaryFaceSets_[patchi].quads;
227  labelList& polys = boundaryFaceSets_[patchi].polys;
228 
229  tris.setSize(p.size());
230  quads.setSize(p.size());
231  polys.setSize(p.size());
232 
233  label nTris = 0;
234  label nQuads = 0;
235  label nPolys = 0;
236 
237  forAll(p, faceI)
238  {
239  const face& f = p[faceI];
240 
241  if (f.size() == 3)
242  {
243  tris[nTris++] = faceI;
244  }
245  else if (f.size() == 4)
246  {
247  quads[nQuads++] = faceI;
248  }
249  else
250  {
251  polys[nPolys++] = faceI;
252  }
253  }
254 
255  tris.setSize(nTris);
256  quads.setSize(nQuads);
257  polys.setSize(nPolys);
258  }
259  }
260  }
261 
262 
263  forAll(allPatchNames_, patchi)
264  {
265  const word& patchName = allPatchNames_[patchi];
266  nFacePrimitives nfp;
267 
268  if (patchNames_.empty() || patchNames_.found(patchName))
269  {
270  if (mesh.boundary()[patchi].size())
271  {
272  nfp.nPoints = mesh.boundaryMesh()[patchi].localPoints().size();
273  nfp.nTris = boundaryFaceSets_[patchi].tris.size();
274  nfp.nQuads = boundaryFaceSets_[patchi].quads.size();
275  nfp.nPolys = boundaryFaceSets_[patchi].polys.size();
276  }
277  }
278 
279  reduce(nfp.nPoints, sumOp<label>());
280  reduce(nfp.nTris, sumOp<label>());
281  reduce(nfp.nQuads, sumOp<label>());
282  reduce(nfp.nPolys, sumOp<label>());
283 
284  nPatchPrims_.insert(patchName, nfp);
285  }
286 }
287 
288 
289 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
290 
292 {}
293 
294 
295 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
296 
297 void Foam::ensightMesh::writePoints
298 (
299  const scalarField& pointsComponent,
300  OFstream& ensightGeometryFile
301 ) const
302 {
303  forAll(pointsComponent, pointI)
304  {
305  ensightGeometryFile<< setw(12) << float(pointsComponent[pointI]) << nl;
306  }
307 }
308 
309 
310 Foam::cellShapeList Foam::ensightMesh::map
311 (
312  const cellShapeList& cellShapes,
313  const labelList& prims
314 ) const
315 {
316  cellShapeList mcsl(prims.size());
317 
318  forAll(prims, i)
319  {
320  mcsl[i] = cellShapes[prims[i]];
321  }
322 
323  return mcsl;
324 }
325 
326 
327 Foam::cellShapeList Foam::ensightMesh::map
328 (
329  const cellShapeList& cellShapes,
330  const labelList& hexes,
331  const labelList& wedges
332 ) const
333 {
334  cellShapeList mcsl(hexes.size() + wedges.size());
335 
336  forAll(hexes, i)
337  {
338  mcsl[i] = cellShapes[hexes[i]];
339  }
340 
341  label offset = hexes.size();
342 
343  const cellModel& hex = *(cellModeller::lookup("hex"));
344  labelList hexLabels(8);
345 
346  forAll(wedges, i)
347  {
348  const cellShape& cellPoints = cellShapes[wedges[i]];
349 
350  hexLabels[0] = cellPoints[0];
351  hexLabels[1] = cellPoints[1];
352  hexLabels[2] = cellPoints[0];
353  hexLabels[3] = cellPoints[2];
354  hexLabels[4] = cellPoints[3];
355  hexLabels[5] = cellPoints[4];
356  hexLabels[6] = cellPoints[6];
357  hexLabels[7] = cellPoints[5];
358 
359  mcsl[i + offset] = cellShape(hex, hexLabels);
360  }
361 
362  return mcsl;
363 }
364 
365 
366 void Foam::ensightMesh::writePrims
367 (
368  const cellShapeList& cellShapes,
369  const label pointOffset,
370  OFstream& ensightGeometryFile
371 ) const
372 {
373  label po = pointOffset + 1;
374 
375  forAll(cellShapes, i)
376  {
377  const cellShape& cellPoints = cellShapes[i];
378 
379  forAll(cellPoints, pointI)
380  {
381  ensightGeometryFile<< setw(10) << cellPoints[pointI] + po;
382  }
383  ensightGeometryFile << nl;
384  }
385 }
386 
387 
388 void Foam::ensightMesh::writePrimsBinary
389 (
390  const cellShapeList& cellShapes,
391  const label pointOffset,
392  std::ofstream& ensightGeometryFile
393 ) const
394 {
395  label po = pointOffset + 1;
396 
397  // Create a temp int array
398  int numElem;
399 
400  numElem = cellShapes.size();
401 
402  if (cellShapes.size())
403  {
404  // All the cellShapes have the same number of elements!
405  int numIntElem = cellShapes.size()*cellShapes[0].size();
406  List<int> temp(numIntElem);
407 
408  int n = 0;
409 
410  forAll(cellShapes, i)
411  {
412  const cellShape& cellPoints = cellShapes[i];
413 
414  forAll(cellPoints, pointI)
415  {
416  temp[n] = cellPoints[pointI] + po;
417  n++;
418  }
419  }
420 
421  ensightGeometryFile.write
422  (
423  reinterpret_cast<char*>(temp.begin()),
424  numIntElem*sizeof(int)
425  );
426  }
427 }
428 
429 
430 void Foam::ensightMesh::writePolysNFaces
431 (
432  const labelList& polys,
433  const cellList& cellFaces,
434  OFstream& ensightGeometryFile
435 ) const
436 {
437  forAll(polys, i)
438  {
439  ensightGeometryFile
440  << setw(10) << cellFaces[polys[i]].size() << nl;
441  }
442 }
443 
444 
445 void Foam::ensightMesh::writePolysNPointsPerFace
446 (
447  const labelList& polys,
448  const cellList& cellFaces,
449  const faceList& faces,
450  OFstream& ensightGeometryFile
451 ) const
452 {
453  forAll(polys, i)
454  {
455  const labelList& cf = cellFaces[polys[i]];
456 
457  forAll(cf, faceI)
458  {
459  ensightGeometryFile
460  << setw(10) << faces[cf[faceI]].size() << nl;
461  }
462  }
463 }
464 
465 
466 void Foam::ensightMesh::writePolysPoints
467 (
468  const labelList& polys,
469  const cellList& cellFaces,
470  const faceList& faces,
471  const label pointOffset,
472  OFstream& ensightGeometryFile
473 ) const
474 {
475  label po = pointOffset + 1;
476 
477  forAll(polys, i)
478  {
479  const labelList& cf = cellFaces[polys[i]];
480 
481  forAll(cf, faceI)
482  {
483  const face& f = faces[cf[faceI]];
484 
485  forAll(f, pointI)
486  {
487  ensightGeometryFile << setw(10) << f[pointI] + po;
488  }
489  ensightGeometryFile << nl;
490  }
491  }
492 }
493 
494 
495 void Foam::ensightMesh::writeAllPolys
496 (
497  const labelList& pointOffsets,
498  OFstream& ensightGeometryFile
499 ) const
500 {
501  if (meshCellSets_.nPolys)
502  {
503  const cellList& cellFaces = mesh_.cells();
504  const faceList& faces = mesh_.faces();
505 
506  if (Pstream::master())
507  {
508  ensightGeometryFile
509  << "nfaced" << nl << setw(10) << meshCellSets_.nPolys << nl;
510  }
511 
512  // Number of faces for each poly cell
513  if (Pstream::master())
514  {
515  // Master
516  writePolysNFaces
517  (
518  meshCellSets_.polys,
519  cellFaces,
520  ensightGeometryFile
521  );
522  // Slaves
523  for (int slave=1; slave<Pstream::nProcs(); slave++)
524  {
525  IPstream fromSlave(Pstream::scheduled, slave);
526  labelList polys(fromSlave);
527  cellList cellFaces(fromSlave);
528 
529  writePolysNFaces
530  (
531  polys,
532  cellFaces,
533  ensightGeometryFile
534  );
535  }
536  }
537  else
538  {
540  toMaster<< meshCellSets_.polys << cellFaces;
541  }
542 
543  // Number of points for each face of the above list
544  if (Pstream::master())
545  {
546  // Master
547  writePolysNPointsPerFace
548  (
549  meshCellSets_.polys,
550  cellFaces,
551  faces,
552  ensightGeometryFile
553  );
554  // Slaves
555  for (int slave=1; slave<Pstream::nProcs(); slave++)
556  {
557  IPstream fromSlave(Pstream::scheduled, slave);
558  labelList polys(fromSlave);
559  cellList cellFaces(fromSlave);
560  faceList faces(fromSlave);
561 
562  writePolysNPointsPerFace
563  (
564  polys,
565  cellFaces,
566  faces,
567  ensightGeometryFile
568  );
569  }
570  }
571  else
572  {
574  toMaster<< meshCellSets_.polys << cellFaces << faces;
575  }
576 
577  // List of points id for each face of the above list
578  if (Pstream::master())
579  {
580  // Master
581  writePolysPoints
582  (
583  meshCellSets_.polys,
584  cellFaces,
585  faces,
586  0,
587  ensightGeometryFile
588  );
589  // Slaves
590  for (int slave=1; slave<Pstream::nProcs(); slave++)
591  {
592  IPstream fromSlave(Pstream::scheduled, slave);
593  labelList polys(fromSlave);
594  cellList cellFaces(fromSlave);
595  faceList faces(fromSlave);
596 
597  writePolysPoints
598  (
599  polys,
600  cellFaces,
601  faces,
602  pointOffsets[slave-1],
603  ensightGeometryFile
604  );
605  }
606  }
607  else
608  {
610  toMaster<< meshCellSets_.polys << cellFaces << faces;
611  }
612  }
613 }
614 
615 
616 void Foam::ensightMesh::writePolysNFacesBinary
617 (
618  const labelList& polys,
619  const cellList& cellFaces,
620  std::ofstream& ensightGeometryFile
621 ) const
622 {
623  forAll(polys, i)
624  {
626  (
627  cellFaces[polys[i]].size(),
628  ensightGeometryFile
629  );
630  }
631 }
632 
633 
634 void Foam::ensightMesh::writePolysNPointsPerFaceBinary
635 (
636  const labelList& polys,
637  const cellList& cellFaces,
638  const faceList& faces,
639  std::ofstream& ensightGeometryFile
640 ) const
641 {
642  forAll(polys, i)
643  {
644  const labelList& cf = cellFaces[polys[i]];
645 
646  forAll(cf, faceI)
647  {
649  (
650  faces[cf[faceI]].size(),
651  ensightGeometryFile
652  );
653  }
654  }
655 }
656 
657 
658 void Foam::ensightMesh::writePolysPointsBinary
659 (
660  const labelList& polys,
661  const cellList& cellFaces,
662  const faceList& faces,
663  const label pointOffset,
664  std::ofstream& ensightGeometryFile
665 ) const
666 {
667  label po = pointOffset + 1;
668 
669  forAll(polys, i)
670  {
671  const labelList& cf = cellFaces[polys[i]];
672 
673  forAll(cf, faceI)
674  {
675  const face& f = faces[cf[faceI]];
676 
677  forAll(f, pointI)
678  {
679  writeEnsDataBinary(f[pointI] + po,ensightGeometryFile);
680  }
681  }
682  }
683 }
684 
685 
686 void Foam::ensightMesh::writeAllPolysBinary
687 (
688  const labelList& pointOffsets,
689  std::ofstream& ensightGeometryFile
690 ) const
691 {
692  if (meshCellSets_.nPolys)
693  {
694  const cellList& cellFaces = mesh_.cells();
695  const faceList& faces = mesh_.faces();
696 
697  if (Pstream::master())
698  {
699  writeEnsDataBinary("nfaced",ensightGeometryFile);
700  writeEnsDataBinary(meshCellSets_.nPolys,ensightGeometryFile);
701  }
702 
703  // Number of faces for each poly cell
704  if (Pstream::master())
705  {
706  // Master
707  writePolysNFacesBinary
708  (
709  meshCellSets_.polys,
710  cellFaces,
711  ensightGeometryFile
712  );
713  // Slaves
714  for (int slave=1; slave<Pstream::nProcs(); slave++)
715  {
716  IPstream fromSlave(Pstream::scheduled, slave);
717  labelList polys(fromSlave);
718  cellList cellFaces(fromSlave);
719 
720  writePolysNFacesBinary
721  (
722  polys,
723  cellFaces,
724  ensightGeometryFile
725  );
726  }
727  }
728  else
729  {
731  toMaster<< meshCellSets_.polys << cellFaces;
732  }
733 
734  // Number of points for each face of the above list
735  if (Pstream::master())
736  {
737  // Master
738  writePolysNPointsPerFaceBinary
739  (
740  meshCellSets_.polys,
741  cellFaces,
742  faces,
743  ensightGeometryFile
744  );
745  // Slaves
746  for (int slave=1; slave<Pstream::nProcs(); slave++)
747  {
748  IPstream fromSlave(Pstream::scheduled, slave);
749  labelList polys(fromSlave);
750  cellList cellFaces(fromSlave);
751  faceList faces(fromSlave);
752 
753  writePolysNPointsPerFaceBinary
754  (
755  polys,
756  cellFaces,
757  faces,
758  ensightGeometryFile
759  );
760  }
761  }
762  else
763  {
765  toMaster<< meshCellSets_.polys << cellFaces << faces;
766  }
767 
768  // List of points id for each face of the above list
769  if (Pstream::master())
770  {
771  // Master
772  writePolysPointsBinary
773  (
774  meshCellSets_.polys,
775  cellFaces,
776  faces,
777  0,
778  ensightGeometryFile
779  );
780  // Slaves
781  for (int slave=1; slave<Pstream::nProcs(); slave++)
782  {
783  IPstream fromSlave(Pstream::scheduled, slave);
784  labelList polys(fromSlave);
785  cellList cellFaces(fromSlave);
786  faceList faces(fromSlave);
787 
788  writePolysPointsBinary
789  (
790  polys,
791  cellFaces,
792  faces,
793  pointOffsets[slave-1],
794  ensightGeometryFile
795  );
796  }
797  }
798  else
799  {
801  toMaster<< meshCellSets_.polys << cellFaces << faces;
802  }
803  }
804 }
805 
806 
807 void Foam::ensightMesh::writeAllPrims
808 (
809  const char* key,
810  const label nPrims,
811  const cellShapeList& cellShapes,
812  const labelList& pointOffsets,
813  OFstream& ensightGeometryFile
814 ) const
815 {
816  if (nPrims)
817  {
818  if (Pstream::master())
819  {
820  ensightGeometryFile << key << nl << setw(10) << nPrims << nl;
821 
822  writePrims(cellShapes, 0, ensightGeometryFile);
823 
824  for (int slave=1; slave<Pstream::nProcs(); slave++)
825  {
826  IPstream fromSlave(Pstream::scheduled, slave);
827  cellShapeList cellShapes(fromSlave);
828 
829  writePrims
830  (
831  cellShapes,
832  pointOffsets[slave-1],
833  ensightGeometryFile
834  );
835  }
836  }
837  else
838  {
840  toMaster<< cellShapes;
841  }
842  }
843 }
844 
845 
846 void Foam::ensightMesh::writeAllPrimsBinary
847 (
848  const char* key,
849  const label nPrims,
850  const cellShapeList& cellShapes,
851  const labelList& pointOffsets,
852  std::ofstream& ensightGeometryFile
853 ) const
854 {
855  if (nPrims)
856  {
857  if (Pstream::master())
858  {
859  writeEnsDataBinary(key,ensightGeometryFile);
860  writeEnsDataBinary(nPrims,ensightGeometryFile);
861 
862  writePrimsBinary(cellShapes, 0, ensightGeometryFile);
863 
864  for (int slave=1; slave<Pstream::nProcs(); slave++)
865  {
866  IPstream fromSlave(Pstream::scheduled, slave);
867  cellShapeList cellShapes(fromSlave);
868 
869  writePrimsBinary
870  (
871  cellShapes,
872  pointOffsets[slave-1],
873  ensightGeometryFile
874  );
875  }
876  }
877  else
878  {
880  toMaster<< cellShapes;
881  }
882  }
883 }
884 
885 
886 void Foam::ensightMesh::writeFacePrims
887 (
888  const faceList& patchFaces,
889  const label pointOffset,
890  OFstream& ensightGeometryFile
891 ) const
892 {
893  if (patchFaces.size())
894  {
895  label po = pointOffset + 1;
896 
897  forAll(patchFaces, i)
898  {
899  const face& patchFace = patchFaces[i];
900 
901  forAll(patchFace, pointI)
902  {
903  ensightGeometryFile << setw(10) << patchFace[pointI] + po;
904  }
905  ensightGeometryFile << nl;
906  }
907  }
908 }
909 
910 
911 void Foam::ensightMesh::writeFacePrimsBinary
912 (
913  const faceList& patchFaces,
914  const label pointOffset,
915  std::ofstream& ensightGeometryFile
916 ) const
917 {
918  if (patchFaces.size())
919  {
920  label po = pointOffset + 1;
921 
922  forAll(patchFaces, i)
923  {
924  const face& patchFace = patchFaces[i];
925 
926  forAll(patchFace, pointI)
927  {
929  (
930  patchFace[pointI] + po,
931  ensightGeometryFile
932  );
933  }
934  }
935  }
936 }
937 
938 
939 Foam::faceList Foam::ensightMesh::map
940 (
941  const faceList& patchFaces,
942  const labelList& prims
943 ) const
944 {
945  faceList ppf(prims.size());
946 
947  forAll (prims, i)
948  {
949  ppf[i] = patchFaces[prims[i]];
950  }
951 
952  return ppf;
953 }
954 
955 
956 void Foam::ensightMesh::writeAllFacePrims
957 (
958  const char* key,
959  const labelList& prims,
960  const label nPrims,
961  const faceList& patchFaces,
962  const labelList& pointOffsets,
963  const labelList& patchProcessors,
964  OFstream& ensightGeometryFile
965 ) const
966 {
967  if (nPrims)
968  {
969  if (Pstream::master())
970  {
971  ensightGeometryFile << key << nl << setw(10) << nPrims << nl;
972 
973  if (&prims != NULL)
974  {
975  writeFacePrims
976  (
977  map(patchFaces, prims),
978  0,
979  ensightGeometryFile
980  );
981  }
982 
983  forAll (patchProcessors, i)
984  {
985  if (patchProcessors[i] != 0)
986  {
987  label slave = patchProcessors[i];
988  IPstream fromSlave(Pstream::scheduled, slave);
989  faceList patchFaces(fromSlave);
990 
991  writeFacePrims
992  (
993  patchFaces,
994  pointOffsets[i],
995  ensightGeometryFile
996  );
997  }
998  }
999  }
1000  else if (&prims != NULL)
1001  {
1003  toMaster<< map(patchFaces, prims);
1004  }
1005  }
1006 }
1007 
1008 
1009 void Foam::ensightMesh::writeNSidedNPointsPerFace
1010 (
1011  const faceList& patchFaces,
1012  OFstream& ensightGeometryFile
1013 ) const
1014 {
1015  forAll(patchFaces, i)
1016  {
1017  ensightGeometryFile
1018  << setw(10) << patchFaces[i].size() << nl;
1019  }
1020 }
1021 
1022 
1023 void Foam::ensightMesh::writeNSidedPoints
1024 (
1025  const faceList& patchFaces,
1026  const label pointOffset,
1027  OFstream& ensightGeometryFile
1028 ) const
1029 {
1030  writeFacePrims
1031  (
1032  patchFaces,
1033  pointOffset,
1034  ensightGeometryFile
1035  );
1036 }
1037 
1038 
1039 void Foam::ensightMesh::writeAllNSided
1040 (
1041  const labelList& prims,
1042  const label nPrims,
1043  const faceList& patchFaces,
1044  const labelList& pointOffsets,
1045  const labelList& patchProcessors,
1046  OFstream& ensightGeometryFile
1047 ) const
1048 {
1049  if (nPrims)
1050  {
1051  if (Pstream::master())
1052  {
1053  ensightGeometryFile
1054  << "nsided" << nl << setw(10) << nPrims << nl;
1055  }
1056 
1057  // Number of points for each face
1058  if (Pstream::master())
1059  {
1060  if (&prims != NULL)
1061  {
1062  writeNSidedNPointsPerFace
1063  (
1064  map(patchFaces, prims),
1065  ensightGeometryFile
1066  );
1067  }
1068 
1069  forAll (patchProcessors, i)
1070  {
1071  if (patchProcessors[i] != 0)
1072  {
1073  label slave = patchProcessors[i];
1074  IPstream fromSlave(Pstream::scheduled, slave);
1075  faceList patchFaces(fromSlave);
1076 
1077  writeNSidedNPointsPerFace
1078  (
1079  patchFaces,
1080  ensightGeometryFile
1081  );
1082  }
1083  }
1084  }
1085  else if (&prims != NULL)
1086  {
1088  toMaster<< map(patchFaces, prims);
1089  }
1090 
1091  // List of points id for each face
1092  if (Pstream::master())
1093  {
1094  if (&prims != NULL)
1095  {
1096  writeNSidedPoints
1097  (
1098  map(patchFaces, prims),
1099  0,
1100  ensightGeometryFile
1101  );
1102  }
1103 
1104  forAll (patchProcessors, i)
1105  {
1106  if (patchProcessors[i] != 0)
1107  {
1108  label slave = patchProcessors[i];
1109  IPstream fromSlave(Pstream::scheduled, slave);
1110  faceList patchFaces(fromSlave);
1111 
1112  writeNSidedPoints
1113  (
1114  patchFaces,
1115  pointOffsets[i],
1116  ensightGeometryFile
1117  );
1118  }
1119  }
1120  }
1121  else if (&prims != NULL)
1122  {
1124  toMaster<< map(patchFaces, prims);
1125  }
1126  }
1127 }
1128 
1129 
1130 void Foam::ensightMesh::writeNSidedPointsBinary
1131 (
1132  const faceList& patchFaces,
1133  const label pointOffset,
1134  std::ofstream& ensightGeometryFile
1135 ) const
1136 {
1137  writeFacePrimsBinary
1138  (
1139  patchFaces,
1140  pointOffset,
1141  ensightGeometryFile
1142  );
1143 }
1144 
1145 
1146 void Foam::ensightMesh::writeNSidedNPointsPerFaceBinary
1147 (
1148  const faceList& patchFaces,
1149  std::ofstream& ensightGeometryFile
1150 ) const
1151 {
1152  forAll(patchFaces, i)
1153  {
1155  (
1156  patchFaces[i].size(),
1157  ensightGeometryFile
1158  );
1159  }
1160 }
1161 
1162 
1163 void Foam::ensightMesh::writeAllNSidedBinary
1164 (
1165  const labelList& prims,
1166  const label nPrims,
1167  const faceList& patchFaces,
1168  const labelList& pointOffsets,
1169  const labelList& patchProcessors,
1170  std::ofstream& ensightGeometryFile
1171 ) const
1172 {
1173  if (nPrims)
1174  {
1175  if (Pstream::master())
1176  {
1177  writeEnsDataBinary("nsided",ensightGeometryFile);
1178  writeEnsDataBinary(nPrims,ensightGeometryFile);
1179  }
1180 
1181  // Number of points for each face
1182  if (Pstream::master())
1183  {
1184  if (&prims != NULL)
1185  {
1186  writeNSidedNPointsPerFaceBinary
1187  (
1188  map(patchFaces, prims),
1189  ensightGeometryFile
1190  );
1191  }
1192 
1193  forAll (patchProcessors, i)
1194  {
1195  if (patchProcessors[i] != 0)
1196  {
1197  label slave = patchProcessors[i];
1198  IPstream fromSlave(Pstream::scheduled, slave);
1199  faceList patchFaces(fromSlave);
1200 
1201  writeNSidedNPointsPerFaceBinary
1202  (
1203  patchFaces,
1204  ensightGeometryFile
1205  );
1206  }
1207  }
1208  }
1209  else if (&prims != NULL)
1210  {
1212  toMaster<< map(patchFaces, prims);
1213  }
1214 
1215  // List of points id for each face
1216  if (Pstream::master())
1217  {
1218  if (&prims != NULL)
1219  {
1220  writeNSidedPointsBinary
1221  (
1222  map(patchFaces, prims),
1223  0,
1224  ensightGeometryFile
1225  );
1226  }
1227 
1228  forAll (patchProcessors, i)
1229  {
1230  if (patchProcessors[i] != 0)
1231  {
1232  label slave = patchProcessors[i];
1233  IPstream fromSlave(Pstream::scheduled, slave);
1234  faceList patchFaces(fromSlave);
1235 
1236  writeNSidedPointsBinary
1237  (
1238  patchFaces,
1239  pointOffsets[i],
1240  ensightGeometryFile
1241  );
1242  }
1243  }
1244  }
1245  else if (&prims != NULL)
1246  {
1248  toMaster<< map(patchFaces, prims);
1249  }
1250  }
1251 }
1252 
1253 
1254 void Foam::ensightMesh::writeAllFacePrimsBinary
1255 (
1256  const char* key,
1257  const labelList& prims,
1258  const label nPrims,
1259  const faceList& patchFaces,
1260  const labelList& pointOffsets,
1261  const labelList& patchProcessors,
1262  std::ofstream& ensightGeometryFile
1263 ) const
1264 {
1265  if (nPrims)
1266  {
1267  if (Pstream::master())
1268  {
1269  writeEnsDataBinary(key,ensightGeometryFile);
1270  writeEnsDataBinary(nPrims,ensightGeometryFile);
1271 
1272  if (&prims != NULL)
1273  {
1274  writeFacePrimsBinary
1275  (
1276  map(patchFaces, prims),
1277  0,
1278  ensightGeometryFile
1279  );
1280  }
1281 
1282  forAll (patchProcessors, i)
1283  {
1284  if (patchProcessors[i] != 0)
1285  {
1286  label slave = patchProcessors[i];
1287  IPstream fromSlave(Pstream::scheduled, slave);
1288  faceList patchFaces(fromSlave);
1289 
1290  writeFacePrimsBinary
1291  (
1292  patchFaces,
1293  pointOffsets[i],
1294  ensightGeometryFile
1295  );
1296  }
1297  }
1298  }
1299  else if (&prims != NULL)
1300  {
1302  toMaster<< map(patchFaces, prims);
1303  }
1304  }
1305 }
1306 
1307 
1309 (
1310  const fileName& postProcPath,
1311  const word& prepend,
1312  const label timeIndex,
1313  Ostream& ensightCaseFile
1314 ) const
1315 {
1316  if (binary_)
1317  {
1318  writeBinary(postProcPath, prepend, timeIndex, ensightCaseFile);
1319  }
1320  else
1321  {
1322  writeAscii(postProcPath, prepend, timeIndex, ensightCaseFile);
1323  }
1324 }
1325 
1326 
1327 void Foam::ensightMesh::writeAscii
1328 (
1329  const fileName& postProcPath,
1330  const word& prepend,
1331  const label timeIndex,
1332  Ostream& ensightCaseFile
1333 ) const
1334 {
1335  const Time& runTime = mesh_.time();
1336  const pointField& points = mesh_.points();
1337  const cellShapeList& cellShapes = mesh_.cellShapes();
1338 
1339  word timeFile = prepend;
1340 
1341  if (timeIndex == 0)
1342  {
1343  timeFile += "000.";
1344  }
1345  else if (mesh_.moving())
1346  {
1347  timeFile += itoa(timeIndex) + '.';
1348  }
1349 
1350  // set the filename of the ensight file
1351  fileName ensightGeometryFileName = timeFile + "mesh";
1352 
1353  OFstream *ensightGeometryFilePtr = NULL;
1354  if (Pstream::master())
1355  {
1356  ensightGeometryFilePtr = new OFstream
1357  (
1358  postProcPath/ensightGeometryFileName,
1359  runTime.writeFormat(),
1360  runTime.writeVersion(),
1362  );
1363  }
1364 
1365  OFstream& ensightGeometryFile = *ensightGeometryFilePtr;
1366 
1367  if (Pstream::master())
1368  {
1369  // Set Format
1370  ensightGeometryFile.setf
1371  (
1373  ios_base::floatfield
1374  );
1375  ensightGeometryFile.precision(5);
1376 
1377  ensightGeometryFile
1378  << "EnSight Geometry File" << nl
1379  << "written from OpenFOAM-" << Foam::FOAMfullVersion << nl
1380  << "node id assign" << nl
1381  << "element id assign" << nl;
1382  }
1383 
1384  labelList pointOffsets(Pstream::nProcs(), 0);
1385 
1386  if (patchNames_.empty())
1387  {
1388  label nPoints = points.size();
1389  Pstream::gather(nPoints, sumOp<label>());
1390 
1391  if (Pstream::master())
1392  {
1393  ensightGeometryFile
1394  << "part" << nl
1395  << setw(10) << 1 << nl
1396  << "internalMesh" << nl
1397  << "coordinates" << nl
1398  << setw(10) << nPoints
1399  << endl;
1400 
1401  for (direction d=0; d<vector::nComponents; d++)
1402  {
1403  writePoints(points.component(d), ensightGeometryFile);
1404  pointOffsets[0] = points.size();
1405 
1406  for (int slave=1; slave<Pstream::nProcs(); slave++)
1407  {
1408  IPstream fromSlave(Pstream::scheduled, slave);
1409  scalarField pointsComponent(fromSlave);
1410  writePoints(pointsComponent, ensightGeometryFile);
1411  pointOffsets[slave] =
1412  pointOffsets[slave-1]
1413  + pointsComponent.size();
1414  }
1415  }
1416  }
1417  else
1418  {
1419  for (direction d=0; d<vector::nComponents; d++)
1420  {
1422  toMaster<< points.component(d);
1423  }
1424  }
1425 
1426  writeAllPrims
1427  (
1428  "hexa8",
1429  meshCellSets_.nHexesWedges,
1430  map(cellShapes, meshCellSets_.hexes, meshCellSets_.wedges),
1431  pointOffsets,
1432  ensightGeometryFile
1433  );
1434 
1435  writeAllPrims
1436  (
1437  "penta6",
1438  meshCellSets_.nPrisms,
1439  map(cellShapes, meshCellSets_.prisms),
1440  pointOffsets,
1441  ensightGeometryFile
1442  );
1443 
1444  writeAllPrims
1445  (
1446  "pyramid5",
1447  meshCellSets_.nPyrs,
1448  map(cellShapes, meshCellSets_.pyrs),
1449  pointOffsets,
1450  ensightGeometryFile
1451  );
1452 
1453  writeAllPrims
1454  (
1455  "tetra4",
1456  meshCellSets_.nTets,
1457  map(cellShapes, meshCellSets_.tets),
1458  pointOffsets,
1459  ensightGeometryFile
1460  );
1461 
1462  writeAllPolys
1463  (
1464  pointOffsets,
1465  ensightGeometryFile
1466  );
1467  }
1468 
1469 
1470  label ensightPatchI = patchPartOffset_;
1471 
1472  forAll(allPatchNames_, patchi)
1473  {
1474  const word& patchName = allPatchNames_[patchi];
1475  const labelList& patchProcessors = allPatchProcs_[patchi];
1476 
1477  if (patchNames_.empty() || patchNames_.found(patchName))
1478  {
1479  const nFacePrimitives& nfp = nPatchPrims_.find(patchName)();
1480 
1481  const labelList *trisPtr = NULL;
1482  const labelList *quadsPtr = NULL;
1483  const labelList *polysPtr = NULL;
1484 
1485  const pointField *patchPointsPtr = NULL;
1486  const faceList *patchFacesPtr = NULL;
1487 
1488  if (mesh_.boundary()[patchi].size())
1489  {
1490  const polyPatch& p = mesh_.boundaryMesh()[patchi];
1491 
1492  trisPtr = &boundaryFaceSets_[patchi].tris;
1493  quadsPtr = &boundaryFaceSets_[patchi].quads;
1494  polysPtr = &boundaryFaceSets_[patchi].polys;
1495 
1496  patchPointsPtr = &(p.localPoints());
1497  patchFacesPtr = &(p.localFaces());
1498  }
1499 
1500  const labelList& tris = *trisPtr;
1501  const labelList& quads = *quadsPtr;
1502  const labelList& polys = *polysPtr;
1503  const pointField& patchPoints = *patchPointsPtr;
1504  const faceList& patchFaces = *patchFacesPtr;
1505 
1506  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
1507  {
1508  labelList patchPointOffsets(Pstream::nProcs(), 0);
1509 
1510  if (Pstream::master())
1511  {
1512  ensightGeometryFile
1513  << "part" << nl
1514  << setw(10) << ensightPatchI++ << nl
1515  << patchName << nl
1516  << "coordinates" << nl
1517  << setw(10) << nfp.nPoints
1518  << endl;
1519 
1520  for (direction d=0; d<vector::nComponents; d++)
1521  {
1522  if (patchPointsPtr)
1523  {
1524  writePoints
1525  (
1526  patchPoints.component(d),
1527  ensightGeometryFile
1528  );
1529  }
1530 
1531  patchPointOffsets = 0;
1532 
1533  forAll (patchProcessors, i)
1534  {
1535  if (patchProcessors[i] != 0)
1536  {
1537  label slave = patchProcessors[i];
1538  IPstream fromSlave(Pstream::scheduled, slave);
1539  scalarField patchPointsComponent(fromSlave);
1540 
1541  writePoints
1542  (
1543  patchPointsComponent,
1544  ensightGeometryFile
1545  );
1546 
1547  if (i < Pstream::nProcs()-1)
1548  {
1549  patchPointOffsets[i+1] =
1550  patchPointOffsets[i]
1551  + patchPointsComponent.size();
1552  }
1553  }
1554  else
1555  {
1556  if (i < Pstream::nProcs()-1)
1557  {
1558  patchPointOffsets[i+1] =
1559  patchPointOffsets[i]
1560  + patchPoints.size();
1561  }
1562  }
1563  }
1564  }
1565  }
1566  else if (patchPointsPtr)
1567  {
1568  for (direction d=0; d<vector::nComponents; d++)
1569  {
1570  OPstream toMaster
1571  (
1574  );
1575  toMaster<< patchPoints.component(d);
1576  }
1577  }
1578 
1579  writeAllFacePrims
1580  (
1581  "tria3",
1582  tris,
1583  nfp.nTris,
1584  patchFaces,
1585  patchPointOffsets,
1586  patchProcessors,
1587  ensightGeometryFile
1588  );
1589 
1590  writeAllFacePrims
1591  (
1592  "quad4",
1593  quads,
1594  nfp.nQuads,
1595  patchFaces,
1596  patchPointOffsets,
1597  patchProcessors,
1598  ensightGeometryFile
1599  );
1600 
1601  writeAllNSided
1602  (
1603  polys,
1604  nfp.nPolys,
1605  patchFaces,
1606  patchPointOffsets,
1607  patchProcessors,
1608  ensightGeometryFile
1609  );
1610  }
1611  }
1612  }
1613 
1614  if (Pstream::master())
1615  {
1616  delete ensightGeometryFilePtr;
1617  }
1618 }
1619 
1620 
1621 void Foam::ensightMesh::writeBinary
1622 (
1623  const fileName& postProcPath,
1624  const word& prepend,
1625  const label timeIndex,
1626  Ostream& ensightCaseFile
1627 ) const
1628 {
1629  //const Time& runTime = mesh.time();
1630  const pointField& points = mesh_.points();
1631  const cellShapeList& cellShapes = mesh_.cellShapes();
1632 
1633  word timeFile = prepend;
1634 
1635  if (timeIndex == 0)
1636  {
1637  timeFile += "000.";
1638  }
1639  else if (mesh_.moving())
1640  {
1641  timeFile += itoa(timeIndex) + '.';
1642  }
1643 
1644  // set the filename of the ensight file
1645  fileName ensightGeometryFileName = timeFile + "mesh";
1646 
1647  std::ofstream *ensightGeometryFilePtr = NULL;
1648 
1649  if (Pstream::master())
1650  {
1651  ensightGeometryFilePtr = new std::ofstream
1652  (
1653  (postProcPath/ensightGeometryFileName).c_str(),
1654  ios_base::out | ios_base::binary | ios_base::trunc
1655  );
1656  // Check on file opened?
1657  }
1658 
1659  std::ofstream& ensightGeometryFile = *ensightGeometryFilePtr;
1660 
1661  if (Pstream::master())
1662  {
1663  writeEnsDataBinary("C binary", ensightGeometryFile);
1664  writeEnsDataBinary("EnSight Geometry File", ensightGeometryFile);
1665  writeEnsDataBinary("written from OpenFOAM", ensightGeometryFile);
1666  writeEnsDataBinary("node id assign", ensightGeometryFile);
1667  writeEnsDataBinary("element id assign", ensightGeometryFile);
1668  }
1669 
1670  labelList pointOffsets(Pstream::nProcs(), 0);
1671 
1672  if (patchNames_.empty())
1673  {
1674  label nPoints = points.size();
1675  Pstream::gather(nPoints, sumOp<label>());
1676 
1677  if (Pstream::master())
1678  {
1679  writeEnsDataBinary("part",ensightGeometryFile);
1680  writeEnsDataBinary(1,ensightGeometryFile);
1681  writeEnsDataBinary("internalMesh",ensightGeometryFile);
1682  writeEnsDataBinary("coordinates",ensightGeometryFile);
1683  writeEnsDataBinary(nPoints,ensightGeometryFile);
1684 
1685  for (direction d=0; d<vector::nComponents; d++)
1686  {
1687  //writePointsBinary(points.component(d), ensightGeometryFile);
1688  writeEnsDataBinary(points.component(d), ensightGeometryFile);
1689  pointOffsets[0] = points.size();
1690 
1691  for (int slave=1; slave<Pstream::nProcs(); slave++)
1692  {
1693  IPstream fromSlave(Pstream::scheduled, slave);
1694  scalarField pointsComponent(fromSlave);
1695  //writePointsBinary(pointsComponent, ensightGeometryFile);
1696  writeEnsDataBinary(pointsComponent, ensightGeometryFile);
1697  pointOffsets[slave] =
1698  pointOffsets[slave-1]
1699  + pointsComponent.size();
1700  }
1701  }
1702  }
1703  else
1704  {
1705  for (direction d=0; d<vector::nComponents; d++)
1706  {
1708  toMaster<< points.component(d);
1709  }
1710  }
1711 
1712  writeAllPrimsBinary
1713  (
1714  "hexa8",
1715  meshCellSets_.nHexesWedges,
1716  map(cellShapes, meshCellSets_.hexes, meshCellSets_.wedges),
1717  pointOffsets,
1718  ensightGeometryFile
1719  );
1720 
1721  writeAllPrimsBinary
1722  (
1723  "penta6",
1724  meshCellSets_.nPrisms,
1725  map(cellShapes, meshCellSets_.prisms),
1726  pointOffsets,
1727  ensightGeometryFile
1728  );
1729 
1730  writeAllPrimsBinary
1731  (
1732  "pyramid5",
1733  meshCellSets_.nPyrs,
1734  map(cellShapes, meshCellSets_.pyrs),
1735  pointOffsets,
1736  ensightGeometryFile
1737  );
1738 
1739  writeAllPrimsBinary
1740  (
1741  "tetra4",
1742  meshCellSets_.nTets,
1743  map(cellShapes, meshCellSets_.tets),
1744  pointOffsets,
1745  ensightGeometryFile
1746  );
1747 
1748  writeAllPolysBinary
1749  (
1750  pointOffsets,
1751  ensightGeometryFile
1752  );
1753 
1754  }
1755 
1756  label ensightPatchI = patchPartOffset_;
1757  label iCount = 0;
1758 
1759  forAll(allPatchNames_, patchi)
1760  {
1761  iCount ++;
1762  const word& patchName = allPatchNames_[patchi];
1763  const labelList& patchProcessors = allPatchProcs_[patchi];
1764 
1765  if (patchNames_.empty() || patchNames_.found(patchName))
1766  {
1767  const nFacePrimitives& nfp = nPatchPrims_.find(patchName)();
1768 
1769  const labelList *trisPtr = NULL;
1770  const labelList *quadsPtr = NULL;
1771  const labelList *polysPtr = NULL;
1772 
1773  const pointField *patchPointsPtr = NULL;
1774  const faceList *patchFacesPtr = NULL;
1775 
1776  if (mesh_.boundary()[patchi].size())
1777  {
1778  const polyPatch& p = mesh_.boundaryMesh()[patchi];
1779 
1780  trisPtr = &boundaryFaceSets_[patchi].tris;
1781  quadsPtr = &boundaryFaceSets_[patchi].quads;
1782  polysPtr = &boundaryFaceSets_[patchi].polys;
1783 
1784  patchPointsPtr = &(p.localPoints());
1785  patchFacesPtr = &(p.localFaces());
1786  }
1787 
1788  const labelList& tris = *trisPtr;
1789  const labelList& quads = *quadsPtr;
1790  const labelList& polys = *polysPtr;
1791  const pointField& patchPoints = *patchPointsPtr;
1792  const faceList& patchFaces = *patchFacesPtr;
1793 
1794  if (nfp.nTris || nfp.nQuads || nfp.nPolys)
1795  {
1796  labelList patchPointOffsets(Pstream::nProcs(), 0);
1797 
1798  if (Pstream::master())
1799  {
1800  writeEnsDataBinary("part",ensightGeometryFile);
1801  writeEnsDataBinary(ensightPatchI++,ensightGeometryFile);
1802  //writeEnsDataBinary(patchName.c_str(),ensightGeometryFile);
1803  writeEnsDataBinary(patchName.c_str(),ensightGeometryFile);
1804  writeEnsDataBinary("coordinates",ensightGeometryFile);
1805  writeEnsDataBinary(nfp.nPoints,ensightGeometryFile);
1806 
1807  for (direction d=0; d<vector::nComponents; d++)
1808  {
1809  if (patchPointsPtr)
1810  {
1811  //writePointsBinary
1813  (
1814  patchPoints.component(d),
1815  ensightGeometryFile
1816  );
1817  }
1818 
1819  patchPointOffsets = 0;
1820 
1821 
1822  forAll (patchProcessors, i)
1823  {
1824  if (patchProcessors[i] != 0)
1825  {
1826  label slave = patchProcessors[i];
1827  IPstream fromSlave(Pstream::scheduled, slave);
1828  scalarField patchPointsComponent(fromSlave);
1829 
1830  //writePointsBinary
1832  (
1833  patchPointsComponent,
1834  ensightGeometryFile
1835  );
1836 
1837  if (i < Pstream::nProcs()-1)
1838  {
1839  patchPointOffsets[i+1] =
1840  patchPointOffsets[i]
1841  + patchPointsComponent.size();
1842  }
1843  }
1844  else
1845  {
1846  if (i < Pstream::nProcs()-1)
1847  {
1848  patchPointOffsets[i+1] =
1849  patchPointOffsets[i]
1850  + patchPoints.size();
1851  }
1852  }
1853  }
1854  }
1855  }
1856  else if (patchPointsPtr)
1857  {
1858  for (direction d=0; d<vector::nComponents; d++)
1859  {
1860  OPstream toMaster
1861  (
1864  );
1865  toMaster<< patchPoints.component(d);
1866  }
1867  }
1868 
1869  writeAllFacePrimsBinary
1870  (
1871  "tria3",
1872  tris,
1873  nfp.nTris,
1874  patchFaces,
1875  patchPointOffsets,
1876  patchProcessors,
1877  ensightGeometryFile
1878  );
1879 
1880  writeAllFacePrimsBinary
1881  (
1882  "quad4",
1883  quads,
1884  nfp.nQuads,
1885  patchFaces,
1886  patchPointOffsets,
1887  patchProcessors,
1888  ensightGeometryFile
1889  );
1890 
1891  writeAllNSidedBinary
1892  (
1893  polys,
1894  nfp.nPolys,
1895  patchFaces,
1896  patchPointOffsets,
1897  patchProcessors,
1898  ensightGeometryFile
1899  );
1900  }
1901  }
1902  }
1903 
1904 
1905  if (Pstream::master())
1906  {
1907  delete ensightGeometryFilePtr;
1908  }
1909 }
1910 
1911 
1912 // ************************ vim: set sw=4 sts=4 et: ************************ //