ccmWriterMesh.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2016-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*----------------------------------------------------------------------------*/
27
28#include "ccmWriter.H"
29#include "ccmInternal.H" // include last to avoid any strange interactions
30
31
32// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33
34Foam::label Foam::ccm::writer::prostarCellFaceId
35(
36 const label& cellId,
37 const label& faceI
38) const
39{
40 const faceList& faces = mesh_.faces();
41 const cellShape& shape = mesh_.cellShapes()[cellId];
42 const labelList& cFaces = mesh_.cells()[cellId];
43
44 label cellFaceId = cFaces.find(faceI);
45 label mapIndex = shape.model().index();
46
47 if (ccm::debug > 1)
48 {
49 Info<< " face[" << faceI << "]: " << faces[faceI] << nl
50 << " owner: " << cellId
51 << " cellFace: " << cellFaceId
52 << " cellFaces: " << cFaces
53 << " shape [index " << mapIndex << "] " << shape
54 << endl;
55
56 Info << "cellFaces" << nl;
57 forAll(cFaces, cFaceI)
58 {
59 Info<< " face [" << cFaces[cFaceI] << "] = "
60 << faces[cFaces[cFaceI]] << nl;
61 }
62
63 Info << "shapeFaces" << nl;
64 {
65 const faceList sFaces = shape.faces();
66 forAll(sFaces, sFaceI)
67 {
68 Info<< " sFace[" << sFaceI << "] = "
69 << sFaces[sFaceI] << nl;
70 }
71 }
72 }
73
74 // The face order returned by primitiveMesh::cells() is not
75 // necessarily the same as defined by
76 // primitiveMesh::cellShapes().
77 // Thus do the lookup for registered primitive types.
78 // Finally, remap the cellModel face number to the ProstarFaceId
79 if (prostarShapeLookup_.found(mapIndex))
80 {
81 const faceList sFaces = shape.faces();
82 forAll(sFaces, sFaceI)
83 {
84 if (faces[faceI] == sFaces[sFaceI])
85 {
86 cellFaceId = sFaceI;
87
88 if (ccm::debug > 1)
89 {
90 Info << " FoamCellFace: " << sFaceI;
91 }
92
93 break;
94 }
95 }
96
97 mapIndex = prostarShapeLookup_[mapIndex];
98 cellFaceId = STARCDCore::foamToStarFaceAddr[mapIndex][cellFaceId];
99
100 if (ccm::debug > 1)
101 {
102 Info<< " ProstarShape: " << mapIndex
103 << " ProstarFaceId: " << cellFaceId + 1
104 << endl;
105 }
106 }
107
108 // PROSTAR used 1-based indices
109 return cellFaceId + 1;
110}
111
112
113// write the faces, the number of vertices appears before each entry
114void Foam::ccm::writer::writeFaces
115(
116 const ccmID& nodeId,
117 const ccmID& mapId,
118 bool isBoundary,
119 label size,
120 label start
121) const
122{
123 if (globalState_->hasError() || size <= 0)
124 {
125 return;
126 }
127
128 CCMIOEntity nodeType;
129 if (isBoundary)
130 {
131 nodeType = kCCMIOBoundaryFaces;
132 }
133 else
134 {
135 nodeType = kCCMIOInternalFaces;
136 }
137
138 const faceList& faces = mesh_.faces();
139 const labelList& owner = mesh_.faceOwner();
140 const labelList& neigh = mesh_.faceNeighbour();
141
142 // 1. write vertices to define faces
143
144 // Determine allocation size
145 label streamSize = size;
146 for (label faceI = 0; faceI < size; ++faceI)
147 {
148 streamSize += faces[faceI + start].size();
149 }
150
151 // Build a vertex stream
152 List<int> ccmStream(streamSize);
153 streamSize = 0;
154 for (label faceI = 0; faceI < size; ++faceI)
155 {
156 const labelList& vrtList = faces[faceI + start];
157
158 ccmStream[streamSize++] = vrtList.size();
159 forAll(vrtList, i)
160 {
161 ccmStream[streamSize++] = vrtList[i] + 1;
162 }
163 }
164
165 if (ccm::debug)
166 {
167 Info<<"CCMIOWriteFaces() size:" << size << " streamSize:"
168 << streamSize << endl;
169 }
170
171 CCMIOWriteFaces
172 (
173 &(globalState_->error),
174 nodeId,
175 nodeType,
176 mapId,
177 streamSize,
178 ccmStream.cdata(),
179 kCCMIOStart,
180 kCCMIOEnd
181 );
182
183 if (globalState_->hasError())
184 {
185 return;
186 }
187
188 // 2. write face owner/neighbours
189 if (isBoundary)
190 {
191 streamSize = size;
192 if (ccmStream.size() < streamSize)
193 {
194 ccmStream.setSize(streamSize);
195 }
196
197 for (int faceI = 0; faceI < size; ++faceI)
198 {
199 ccmStream[faceI] = owner[faceI + start] + 1;
200 }
201 }
202 else
203 {
204 // kCCMIOInternalFaces
205 streamSize = 2 * size;
206 if (ccmStream.size() < streamSize)
207 {
208 ccmStream.setSize(streamSize);
209 }
210
211 for (int faceI = 0; faceI < size; ++faceI)
212 {
213 ccmStream[faceI*2] = owner[faceI + start] + 1;
214 ccmStream[faceI*2+1] = neigh[faceI + start] + 1;
215 }
216 }
217
218 if (ccm::debug)
219 {
220 Info<<"CCMIOWriteFaceCells() size:" << size;
221 if (isBoundary)
222 {
223 Info<< " boundary";
224 }
225 else
226 {
227 Info<<" internal";
228 }
229 Info<<"Faces"<<endl;
230 }
231
232 CCMIOWriteFaceCells
233 (
234 &(globalState_->error),
235 nodeId,
236 nodeType,
237 mapId,
238 ccmStream.cdata(),
239 kCCMIOStart,
240 kCCMIOEnd
241 );
242
243 // determine corresponding ProstarFaceId
244 if (isBoundary)
245 {
246 streamSize = size;
247 if (ccmStream.size() < streamSize)
248 {
249 ccmStream.setSize(streamSize);
250 }
251
252 for (int faceIdx = 0; faceIdx < size; ++faceIdx)
253 {
254 label faceI = faceIdx + start;
255 // boundary face - owner only
256 ccmStream[faceIdx] = prostarCellFaceId(owner[faceI], faceI);
257 }
258
259 CCMIOWriteOpt1i
260 (
261 &(globalState_->error),
262 nodeId,
263 "ProstarFaceId",
264 size,
265 ccmStream.cdata(),
266 kCCMIOStart,
267 kCCMIOEnd
268 );
269 }
270 else
271 {
272 // kCCMIOInternalFaces
273 streamSize = 2 * size;
274 if (ccmStream.size() < streamSize)
275 {
276 ccmStream.setSize(streamSize);
277 }
278
279 for (int faceIdx = 0; faceIdx < size; ++faceIdx)
280 {
281 label faceI = faceIdx + start;
282 ccmStream[faceIdx*2] = prostarCellFaceId(owner[faceI], faceI);
283 ccmStream[faceIdx*2+1] = prostarCellFaceId(neigh[faceI], faceI);
284 }
285
286 CCMIOWriteOpt2i
287 (
288 &(globalState_->error),
289 nodeId,
290 "ProstarFaceId",
291 size,
292 2,
293 ccmStream.cdata(),
294 kCCMIOStart,
295 kCCMIOEnd
296 );
297 }
298}
299
300
301// writeVertices
302// 1) write the vertex map (starting with 1)
303// 2) write the vertex data
304//
305void Foam::ccm::writer::writeVertices
306(
307 const ccmID& verticesNode
308) const
309{
310 const pointField& pts = mesh_.points();
311
312 Info<< "writing points: " << pts.size() << endl;
313
314 // 1. mapping data array index to the vertex Id (starting with 1)
315 ccmID vertexMap;
316 addLinearMap
317 (
318 "vertex Map",
319 vertexMap,
320 pts.size()
321 );
322
323 // 2. write vertices - scale [m] -> [mm] for consistency with PROSTAR
324 // but is probably immaterial
325 const float scaling(0.001); // to recover meters
326 List<float> vrts(3*pts.size());
327 forAll(pts, i)
328 {
329 vrts[3*i] = 1000.0 * pts[i].x();
330 vrts[3*i+1] = 1000.0 * pts[i].y();
331 vrts[3*i+2] = 1000.0 * pts[i].z();
332 }
333
334 CCMIOWriteVerticesf
335 (
336 &(globalState_->error),
337 verticesNode,
338 3, scaling,
339 vertexMap,
340 vrts.cdata(),
341 kCCMIOStart,
342 kCCMIOEnd
343 );
344 assertNoError("writing 'Vertices' node");
345}
346
347
348// writeInternalFaces
349// 1) write the face map (starting with 1)
350// 2) write owner/neighbour
351//
352void Foam::ccm::writer::writeInternalFaces
353(
354 const ccmID& topoNode
355) const
356{
357 label nFaces = mesh_.nInternalFaces();
358
359 Info<< "writing internalFaces: " << nFaces << endl;
360 if (nFaces > 0)
361 {
362 ccmID nodeId;
363 CCMIONewEntity
364 (
365 &(globalState_->error),
366 topoNode,
367 kCCMIOInternalFaces,
368 nullptr,
369 &nodeId
370 );
371 assertNoError("creating internalFaces node");
372
373 writeFaces
374 (
375 nodeId,
376 maps_->internalFaces,
377 false,
378 nFaces
379 );
380 assertNoError("writing internalFaces");
381 }
382}
383
384
385// writeBoundaryFaces:
386// - write faces with owner
387void Foam::ccm::writer::writeBoundaryFaces
388(
389 const ccmID& topoNode
390) const
391{
392 const polyBoundaryMesh& patches = mesh_.boundaryMesh();
393
394 // 4. write boundary faces
395 Info<< "writing boundaryFaces:" << flush;
396
397 label defaultId = findDefaultBoundary();
398
399 // write Default_Boundary_Region as BoundaryFaces-0
400 if (defaultId >= 0)
401 {
402 Info<< " " << patches[defaultId].size() << flush;
403
404 ccmID nodeId;
405 CCMIONewIndexedEntity
406 (
407 &(globalState_->error),
408 topoNode,
409 kCCMIOBoundaryFaces,
410 0,
411 nullptr,
412 &nodeId
413 );
414 assertNoError
415 (
416 "creating boundaryFaces node patch: " + patches[defaultId].name()
417 );
418
419 writeFaces
420 (
421 nodeId,
422 maps_->boundary[defaultId],
423 true,
424 patches[defaultId].size(),
425 patches[defaultId].start()
426 );
427 assertNoError
428 (
429 "writing boundaryFaces patch: " + patches[defaultId].name()
430 );
431 }
432
433 //
434 // write boundary faces - skip Default_Boundary_Region
435 //
436 forAll(patches, patchI)
437 {
438 label regionId = patchI;
439 if (regionId == defaultId)
440 {
441 continue; // skip - already written
442 }
443 else if (defaultId == -1 || regionId < defaultId)
444 {
445 regionId++;
446 }
447
448 Info<< " " << patches[patchI].size() << flush;
449
450 if (patches[patchI].size() > 0)
451 {
452 ccmID nodeId;
453 CCMIONewIndexedEntity
454 (
455 &(globalState_->error),
456 topoNode,
457 kCCMIOBoundaryFaces,
458 regionId,
459 nullptr,
460 &nodeId
461 );
462 assertNoError
463 (
464 "creating boundaryFaces node patch: " + patches[patchI].name()
465 );
466
467 writeFaces
468 (
469 nodeId,
470 maps_->boundary[patchI],
471 true,
472 patches[patchI].size(),
473 patches[patchI].start()
474 );
475 assertNoError
476 (
477 "writing boundaryFaces patch: " + patches[patchI].name()
478 );
479 }
480 }
481
482 Info<< endl;
483}
484
485
486// writeCells:
487// - write faces with owner/neighbour
488// - write interfaces
489void Foam::ccm::writer::writeCells
490(
491 const ccmID& topoNode
492)
493{
494 Info<< "writing cells: " << mesh_.nCells() << endl;
495
496 // 1. cellTableId
497 // - if possible, read from constant/polyMesh/cellTableId
498 // - otherwise use cellZone information
499 List<int> mapData(mesh_.nCells(), -1);
500 bool useCellZones = false;
501
502 IOList<label> ioList
503 (
504 IOobject
505 (
506 "cellTableId",
507 "constant",
509 mesh_,
512 false
513 )
514 );
515
516
517 if (ioList.headerOk())
518 {
519 if (ioList.size() == mesh_.nCells())
520 {
521 mapData.transfer(ioList);
522 }
523 else
524 {
526 << ioList.objectPath() << endl
527 << " Has incorrect number of cells: "
528 << ioList.size() << " instead of " << mesh_.nCells()
529 << " - use cellZone information instead"
530 << endl;
531
532 ioList.clear();
533 useCellZones = true;
534 }
535 }
536 else
537 {
538 useCellZones = true;
539 }
540
541
542 if (useCellZones)
543 {
544 if (!cellTable_.size())
545 {
546 Info<< "created cellTable from cellZones" << endl;
547 cellTable_ = mesh_;
548 }
549
550 // track if there are unzoned cells
551 label nUnzoned = mesh_.nCells();
552
553 // get the cellZone <-> cellTable correspondence
554 Info<< "matching cellZones to cellTable" << endl;
555
556 forAll(mesh_.cellZones(), zoneI)
557 {
558 const cellZone& cZone = mesh_.cellZones()[zoneI];
559 if (cZone.size())
560 {
561 nUnzoned -= cZone.size();
562
563 label tableId = cellTable_.findIndex(cZone.name());
564 if (tableId < 0)
565 {
566 dictionary dict;
567
568 dict.add("Label", cZone.name());
569 dict.add("MaterialType", "fluid");
570 tableId = cellTable_.append(dict);
571 }
572
573 for (auto id : cZone)
574 {
575 mapData[id] = tableId;
576 }
577 }
578 }
579
580 if (nUnzoned)
581 {
582 dictionary dict;
583
584 dict.add("Label", "__unzonedCells__");
585 dict.add("MaterialType", "fluid");
586 label tableId = cellTable_.append(dict);
587
588 for (auto& id : mapData)
589 {
590 if (id < 0)
591 {
592 id = tableId;
593 }
594 }
595 }
596 }
597
598 // 2. mapping data array index to the cell Id (starting with 1)
599 ccmID cellsNode;
600 CCMIONewEntity
601 (
602 &(globalState_->error),
603 topoNode,
604 kCCMIOCells,
605 nullptr,
606 &cellsNode
607 );
608 assertNoError("creating 'Cells' node");
609
610 CCMIOWriteCells
611 (
612 &(globalState_->error),
613 cellsNode,
614 maps_->cells,
615 mapData.data(),
616 kCCMIOStart, kCCMIOEnd
617 );
618 assertNoError("writing 'Cells' node");
619
620 // 3. add cell topology information, if possible
621 const cellShapeList& shapes = mesh_.cellShapes();
622 forAll(shapes, cellI)
623 {
624 label mapIndex = shapes[cellI].model().index();
625
626 // A registered primitive type
627 if (prostarShapeLookup_.found(mapIndex))
628 {
629 mapData[cellI] = prostarShapeLookup_[mapIndex];
630 }
631 else
632 {
633 mapData[cellI] = STARCDCore::starcdPoly; // Treat as polyhedral
634 }
635 }
636
637 CCMIOWriteOpt1i
638 (
639 &(globalState_->error),
640 cellsNode,
641 "CellTopologyType",
642 mesh_.nCells(),
643 mapData.cdata(),
644 kCCMIOStart, kCCMIOEnd
645 );
646
647 // 4. write interfaces
648 writeInterfaces(cellsNode);
649}
650
651
652// write interfaces
653// 1) PROSTAR baffles
654//
655void Foam::ccm::writer::writeInterfaces
656(
657 const ccmID& cellsNode
658) const
659{
660 IOList<labelList> interfaces
661 (
662 IOobject
663 (
664 "interfaces",
665 "constant",
667 mesh_,
670 false
671 )
672 );
673
674 if (interfaces.headerOk() && interfaces.size() > 0)
675 {
676 List<int> mapData(2 * interfaces.size());
677
678 forAll(interfaces, i)
679 {
680 mapData[i*2] = interfaces[i][0];
681 mapData[i*2+1] = interfaces[i][1];
682 }
683
684 ccmID nodeId;
685 if
686 (
687 CCMIONewEntity
688 (
689 &(globalState_->error),
690 cellsNode,
691 kCCMIOInterfaces,
692 0,
693 &nodeId
694 )
695 != kCCMIONoErr
696
697 || CCMIOWriteOpt2i
698 (
699 &(globalState_->error),
700 nodeId,
701 "FaceIds",
702 interfaces.size(),
703 2,
704 mapData.cdata(), kCCMIOStart, kCCMIOEnd
705 )
706 != kCCMIONoErr
707 )
708 {
709 assertNoError("writing interfaces 'FaceIds'");
710 }
711 }
712}
713
714
715// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
716
718{
719 // use first ("default") state node
720 ccmID stateNode;
721 // int stateI = 0;
722
723 // create "default" State
724 CCMIONewState
725 (
726 &(globalState_->error),
727 (globalState_->root),
728 "default",
729 nullptr,
730 nullptr,
731 &stateNode
732 );
733 assertNoError("could not create default state");
734
735 // use first processor - create a new one
736 ccmID processorNode;
737 int procI = 0;
738
739 if
740 (
741 CCMIONextEntity
742 (
743 nullptr,
744 stateNode,
745 kCCMIOProcessor,
746 &procI,
747 &processorNode
748 )
749 != kCCMIONoErr
750 )
751 {
752 CCMIONewEntity
753 (
754 &(globalState_->error),
755 stateNode,
756 kCCMIOProcessor,
757 nullptr,
758 &processorNode
759 );
760 assertNoError("could not create processor node");
761 }
762
763 // remove old data (if it exists)
764 CCMIOClearProcessor
765 (
766 nullptr, stateNode, processorNode,
767 true, // Clear vertices
768 true, // Clear topology
769 true, // Clear initial field
770 true, // Clear solution
771 true // Clear lagrangian
772 );
773
774#if 0
775 CCMIOWriteOptstr
776 (
777 nullptr,
778 processorNode,
779 "CreatingProgram",
780 "ccm::writer"
781 );
782#endif
783
784 //
785 // create vertices and topology nodes
786 //
787 ccmID verticesNode, topoNode;
788 if
789 (
790 CCMIONewEntity
791 (
792 &(globalState_->error),
793 (globalState_->root),
794 kCCMIOVertices,
795 nullptr,
796 &verticesNode
797 )
798 == kCCMIONoErr
799
800 && CCMIONewEntity
801 (
802 &(globalState_->error),
803 (globalState_->root),
804 kCCMIOTopology,
805 nullptr,
806 &topoNode
807 )
808 == kCCMIONoErr
809 )
810 {
811 writeVertices(verticesNode);
812 writeInternalFaces(topoNode);
813 writeBoundaryFaces(topoNode);
814 writeCells(topoNode);
815 writeProblem(stateNode);
816 }
817
818 // Now we have the mesh (vertices and topology),
819 // we can write out the processor information.
820 CCMIOWriteProcessor
821 (
822 &(globalState_->error),
823 processorNode,
824 nullptr, &verticesNode, // no verticesFile, write verticesNode
825 nullptr, &topoNode, // no topologyFile, write topoNode
826 nullptr, nullptr, // initialField unchanged
827 nullptr, nullptr // no solutionFile, solution unchanged
828 );
829 assertNoError("Error after writing geometry processor");
830
831}
832
833
834// ************************************************************************* //
Internal bits for wrapping libccmio - do not use directly.
bool found(const Key &key) const
Return true if hashed entry is found in table.
Definition: HashTableI.H:100
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
Definition: UList.C:212
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
void writeGeometry()
Write the mesh.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1108
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:324
const cellShapeList & cellShapes() const
Return cell shapes.
const cellList & cells() const
const polyBoundaryMesh & patches
const labelList nFaces(UPstream::listGatherValues< label >(aMesh.nFaces()))
label cellId
#define WarningInFunction
Report a warning using Foam::Warning.
List< cellShape > cellShapeList
List of cellShapes and PtrList of List of cellShape.
Definition: cellShapeList.H:45
List< label > labelList
A List of labels.
Definition: List.H:66
messageStream Info
Information stream (stdout output on master, null elsewhere)
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:44
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
List< face > faceList
A List of faces.
Definition: faceListFwd.H:47
Ostream & flush(Ostream &os)
Flush stream.
Definition: Ostream.H:364
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333