1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | foam-extend: Open Source CFD
4 \\ / O peration | Version: 3.2
5 \\ / A nd | Web: http://www.foam-extend.org
6 \\/ M anipulation | For copyright notice see file Copyright
7 -------------------------------------------------------------------------------
9 This file is part of foam-extend.
11 foam-extend is free software: you can redistribute it and/or modify it
12 under the terms of the GNU General Public License as published by the
13 Free Software Foundation, either version 3 of the License, or (at your
14 option) any later version.
16 foam-extend is distributed in the hope that it will be useful, but
17 WITHOUT ANY WARRANTY; without even the implied warranty of
18 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19 General Public License for more details.
21 You should have received a copy of the GNU General Public License
22 along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "domainDecomposition.H"
28 #include "dictionary.H"
29 #include "labelIOList.H"
30 #include "processorPolyPatch.H"
32 #include "OSspecific.H"
34 #include "globalMeshData.H"
35 #include "DynamicList.H"
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 domainDecomposition::domainDecomposition(const IOobject& io)
53 nProcs_(readInt(decompositionDict_.lookup("numberOfSubdomains"))),
55 cellToProc_(nCells()),
56 procPointAddressing_(nProcs_),
57 procFaceAddressing_(nProcs_),
58 nInternalProcFaces_(nProcs_),
59 nLiveProcFaces_(nProcs_),
60 procCellAddressing_(nProcs_),
61 procBoundaryAddressing_(nProcs_),
62 procPatchSize_(nProcs_),
63 procPatchStartIndex_(nProcs_),
64 procNeighbourProcessors_(nProcs_),
65 procProcessorPatchSize_(nProcs_),
66 procProcessorPatchStartIndex_(nProcs_),
67 globallySharedPoints_(0),
68 cyclicParallel_(false)
70 if (decompositionDict_.found("distributed"))
72 distributed_ = Switch(decompositionDict_.lookup("distributed"));
77 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
79 domainDecomposition::~domainDecomposition()
83 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 bool domainDecomposition::writeDecomposition()
87 Info<< "\nConstructing processor meshes" << endl;
89 // Make a lookup map for globally shared points
90 Map<label> sharedPointLookup(2*globallySharedPoints_.size());
92 forAll (globallySharedPoints_, pointi)
94 sharedPointLookup.insert(globallySharedPoints_[pointi], pointi);
98 // Mark point/faces/cells that are in zones. Bad coding - removed
101 label totProcFaces = 0;
102 label maxProcPatches = 0;
103 label maxProcFaces = 0;
106 // Write out the meshes
107 for (label procI = 0; procI < nProcs_; procI++)
109 // Create processor points
110 const labelList& curPointLabels = procPointAddressing_[procI];
112 // Access list of all points in the mesh. HJ, 27/Mar/2009
113 const pointField& meshPoints = allPoints();
115 labelList pointLookup(meshPoints.size(), -1);
117 pointField procPoints(curPointLabels.size());
119 forAll (curPointLabels, pointi)
121 procPoints[pointi] = meshPoints[curPointLabels[pointi]];
123 pointLookup[curPointLabels[pointi]] = pointi;
126 // Create processor faces
127 const labelList& curFaceLabels = procFaceAddressing_[procI];
129 // Access list of all faces in the mesh. HJ, 27/Mar/2009
130 const faceList& meshFaces = allFaces();
132 labelList faceLookup(meshFaces.size(), -1);
134 faceList procFaces(curFaceLabels.size());
136 forAll (curFaceLabels, facei)
138 // Mark the original face as used
139 // Remember to decrement the index by one (turning index)
141 label curF = mag(curFaceLabels[facei]) - 1;
143 faceLookup[curF] = facei;
145 // get the original face
146 labelList origFaceLabels;
148 if (curFaceLabels[facei] >= 0)
151 origFaceLabels = meshFaces[curF];
155 origFaceLabels = meshFaces[curF].reverseFace();
158 // translate face labels into local point list
159 face& procFaceLabels = procFaces[facei];
161 procFaceLabels.setSize(origFaceLabels.size());
163 forAll (origFaceLabels, pointi)
165 procFaceLabels[pointi] = pointLookup[origFaceLabels[pointi]];
169 // Create cell lookup
170 labelList cellLookup(nCells(), -1);
171 const labelList& curCellLabels = procCellAddressing_[procI];
173 forAll (curCellLabels, cellI)
175 cellLookup[curCellLabels[cellI]] = cellI;
178 // Get complete owner-neighour addressing in the mesh
179 const labelList& own = faceOwner();
180 const labelList& nei = faceNeighbour();
182 // Calculate owner and neighbour list
183 // Owner list is sized to number of live faces
184 // Neighbour list is sized to number of internal faces
186 labelList procOwner(nLiveProcFaces_[procI]);
188 // Note: loop over owner, not all faces: sizes are different
189 forAll (procOwner, faceI)
191 // Remember to decrement the index by one (turning index)
193 label curF = mag(curFaceLabels[faceI]) - 1;
195 if (curFaceLabels[faceI] >= 0)
197 procOwner[faceI] = cellLookup[own[curF]];
201 procOwner[faceI] = cellLookup[nei[curF]];
205 labelList procNeighbour(nInternalProcFaces_[procI]);
207 // Note: loop over neighbour, not all faces: sizes are different
208 forAll (procNeighbour, faceI)
210 // Remember to decrement the index by one (turning index)
212 label curF = mag(curFaceLabels[faceI]) - 1;
214 if (curFaceLabels[faceI] >= 0)
216 procNeighbour[faceI] = cellLookup[nei[curF]];
220 procNeighbour[faceI] = cellLookup[own[curF]];
224 // Create processor cells. No longer needed: using owner and neighbour
226 // const cellList& meshCells = cells();
228 // cellList procCells(curCellLabels.size());
230 // forAll (curCellLabels, cellI)
232 // const labelList& origCellLabels = meshCells[curCellLabels[cellI]];
234 // cell& curCell = procCells[cellI];
236 // curCell.setSize(origCellLabels.size());
238 // forAll (origCellLabels, cellFaceI)
240 // curCell[cellFaceI] = faceLookup[origCellLabels[cellFaceI]];
244 // Create processor mesh without a boundary
246 fileName processorCasePath
248 time().caseName()/fileName(word("processor") + Foam::name(procI))
251 // make the processor directory
252 mkDir(time().rootPath()/processorCasePath);
257 Time::controlDictName,
270 this->polyMesh::name(), // region name of undecomposed mesh
274 xferMove(procPoints),
277 xferMove(procNeighbour),
278 false // Do not sync par
279 // xferMove(procCells) // Old-fashioned mesh creation using cells.
280 // Deprecated: using face owner/neighbour
284 // Create processor boundary patches
285 const labelList& curBoundaryAddressing = procBoundaryAddressing_[procI];
287 const labelList& curPatchSizes = procPatchSize_[procI];
289 const labelList& curPatchStarts = procPatchStartIndex_[procI];
291 const labelList& curNeighbourProcessors =
292 procNeighbourProcessors_[procI];
294 const labelList& curProcessorPatchSizes =
295 procProcessorPatchSize_[procI];
297 const labelList& curProcessorPatchStarts =
298 procProcessorPatchStartIndex_[procI];
300 const polyPatchList& meshPatches = boundaryMesh();
302 List<polyPatch*> procPatches
305 + curProcessorPatchSizes.size(),
306 reinterpret_cast<polyPatch*>(0)
311 forAll (curPatchSizes, patchi)
313 procPatches[nPatches] =
314 meshPatches[curBoundaryAddressing[patchi]].clone
316 procMesh.boundaryMesh(),
318 curPatchSizes[patchi],
319 curPatchStarts[patchi]
325 forAll (curProcessorPatchSizes, procPatchI)
327 procPatches[nPatches] =
328 new processorPolyPatch
330 word("procBoundary") + Foam::name(procI)
332 + Foam::name(curNeighbourProcessors[procPatchI]),
333 curProcessorPatchSizes[procPatchI],
334 curProcessorPatchStarts[procPatchI],
336 procMesh.boundaryMesh(),
338 curNeighbourProcessors[procPatchI]
344 // Add boundary patches
345 procMesh.addPatches(procPatches);
347 // Create and add zones
350 // This coding was all wrong, as each point/face/cell may only
351 // belong to a single zone.
352 // Additionally, ordering of points/faces/cells in the processor mesh
353 // needs to match the ordering in global mesh zones. Full rewrite.
356 // Create zones if needed
359 pointZones().size() > 0
360 || faceZones().size() > 0
361 || cellZones().size() > 0
365 List<pointZone*> procPz(pointZones().size());
368 const pointZoneMesh& pz = pointZones();
370 // Go through all the zoned points and find out if they
371 // belong to a processor. If so, add it to the zone as
375 const labelList& zonePoints = pz[zoneI];
377 labelList procZonePoints(zonePoints.size());
378 label nZonePoints = 0;
380 forAll (zonePoints, pointI)
382 const label localIndex =
383 pointLookup[zonePoints[pointI]];
387 // Point live on processor: add to zone
388 procZonePoints[nZonePoints] = localIndex;
394 procZonePoints.setSize(nZonePoints);
396 procPz[zoneI] = new pointZone
401 procMesh.pointZones()
408 List<faceZone*> procFz(faceZones().size());
411 const faceZoneMesh& fz = faceZones();
415 const labelList& zoneFaces = fz[zoneI];
416 const boolList& flipMap = fz[zoneI].flipMap();
418 // Go through all the zoned faces and find out if they
419 // belong to a processor. If so, add it to the zone as
422 labelList procZoneFaces(zoneFaces.size());
423 boolList procZoneFaceFlips(zoneFaces.size());
424 label nZoneFaces = 0;
426 forAll (zoneFaces, faceI)
428 const label localIndex = faceLookup[zoneFaces[faceI]];
432 // Face is present on the processor
434 // Add the face to the zone
435 procZoneFaces[nZoneFaces] = localIndex;
438 bool flip = flipMap[faceI];
440 if (curFaceLabels[localIndex] < 0)
445 procZoneFaceFlips[nZoneFaces] = flip;
451 procZoneFaces.setSize(nZoneFaces);
452 procZoneFaceFlips.setSize(nZoneFaces);
454 procFz[zoneI] = new faceZone
466 List<cellZone*> procCz(cellZones().size());
469 const cellZoneMesh& cz = cellZones();
471 // Go through all the zoned cells and find out if they
472 // belong to a processor. If so, add it to the zone as
477 const labelList& zoneCells = cz[zoneI];
479 labelList procZoneCells(zoneCells.size());
480 label nZoneCells = 0;
482 forAll (zoneCells, cellI)
484 const label localIndex = cellLookup[zoneCells[cellI]];
488 procZoneCells[nZoneCells] = localIndex;
494 procZoneCells.setSize(nZoneCells);
496 procCz[zoneI] = new cellZone
507 procMesh.addZones(procPz, procFz, procCz);
510 // Set the precision of the points data to 10
511 IOstream::defaultPrecision(10);
516 << "Processor " << procI << nl
517 << " Number of cells = " << procMesh.nCells()
520 label nBoundaryFaces = 0;
521 label nProcPatches = 0;
522 label nProcFaces = 0;
524 forAll (procMesh.boundaryMesh(), patchi)
528 procMesh.boundaryMesh()[patchi].type()
529 == processorPolyPatch::typeName
532 const processorPolyPatch& ppp =
533 refCast<const processorPolyPatch>
535 procMesh.boundaryMesh()[patchi]
538 Info<< " Number of faces shared with processor "
539 << ppp.neighbProcNo() << " = " << ppp.size() << endl;
542 nProcFaces += ppp.size();
546 nBoundaryFaces += procMesh.boundaryMesh()[patchi].size();
550 Info<< " Number of processor patches = " << nProcPatches << nl
551 << " Number of processor faces = " << nProcFaces << nl
552 << " Number of boundary faces = " << nBoundaryFaces << endl;
554 totProcFaces += nProcFaces;
555 maxProcPatches = max(maxProcPatches, nProcPatches);
556 maxProcFaces = max(maxProcFaces, nProcFaces);
558 // create and write the addressing information
559 labelIOList pointProcAddressing
563 "pointProcAddressing",
564 procMesh.facesInstance(),
570 procPointAddressing_[procI]
572 pointProcAddressing.write();
574 labelIOList faceProcAddressing
578 "faceProcAddressing",
579 procMesh.facesInstance(),
585 procFaceAddressing_[procI]
587 faceProcAddressing.write();
589 labelIOList cellProcAddressing
593 "cellProcAddressing",
594 procMesh.facesInstance(),
600 procCellAddressing_[procI]
602 cellProcAddressing.write();
604 labelIOList boundaryProcAddressing
608 "boundaryProcAddressing",
609 procMesh.facesInstance(),
615 procBoundaryAddressing_[procI]
617 boundaryProcAddressing.write();
621 << "Number of processor faces = " << totProcFaces/2 << nl
622 << "Max number of processor patches = " << maxProcPatches << nl
623 << "Max number of faces between processors = " << maxProcFaces
630 // ************************************************************************* //