applications/utilities/parallelProcessing/decomposePar/distributeCells.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2008 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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "domainDecomposition.H"
  28 #include "decompositionMethod.H"
  29 #include "cpuTime.H"
  30 #include "cyclicPolyPatch.H"
  31 #include "cellSet.H"
  32 #include "regionSplit.H"
  33
  34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  35
  36 void domainDecomposition::distributeCells()
  37 {
  38     Info<< "\nCalculating distribution of cells" << endl;
  39
  40     cpuTime decompositionTime;
  41
  42
  43     // See if any faces need to have owner and neighbour on same processor
  44     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
  45
  46     labelHashSet sameProcFaces;
  47
  48     if (decompositionDict_.found("preservePatches"))
  49     {
  50         wordList pNames(decompositionDict_.lookup("preservePatches"));
  51
  52         Info<< "Keeping owner and neighbour of faces in patches " << pNames
  53             << " on same processor" << endl;
  54
  55         const polyBoundaryMesh& patches = boundaryMesh();
  56
  57         forAll(pNames, i)
  58         {
  59             label patchI = patches.findPatchID(pNames[i]);
  60
  61             if (patchI == -1)
  62             {
  63                 FatalErrorIn("domainDecomposition::distributeCells()")
  64                     << "Unknown preservePatch " << pNames[i]
  65                     << endl << "Valid patches are " << patches.names()
  66                     << exit(FatalError);
  67             }
  68
  69             const polyPatch& pp = patches[patchI];
  70
  71             forAll(pp, i)
  72             {
  73                 sameProcFaces.insert(pp.start() + i);
  74             }
  75         }
  76     }
  77     if (decompositionDict_.found("preserveFaceZones"))
  78     {
  79         wordList zNames(decompositionDict_.lookup("preserveFaceZones"));
  80
  81         Info<< "Keeping owner and neighbour of faces in zones " << zNames
  82             << " on same processor" << endl;
  83
  84         const faceZoneMesh& fZones = faceZones();
  85
  86         forAll(zNames, i)
  87         {
  88             label zoneI = fZones.findZoneID(zNames[i]);
  89
  90             if (zoneI == -1)
  91             {
  92                 FatalErrorIn("domainDecomposition::distributeCells()")
  93                     << "Unknown preserveFaceZone " << zNames[i]
  94                     << endl << "Valid faceZones are " << fZones.names()
  95                     << exit(FatalError);
  96             }
  97
  98             const faceZone& fz = fZones[zoneI];
  99
 100             forAll(fz, i)
 101             {
 102                 sameProcFaces.insert(fz[i]);
 103             }
 104         }
 105     }
 106
 107     if (sameProcFaces.size() > 0)
 108     {
 109         Info<< "Selected " << sameProcFaces.size()
 110             << " faces whose owner and neighbour cell should be kept on the"
 111             << " same processor" << endl;
 112     }
 113
 114
 115
 116     // Construct decomposition method and either do decomposition on
 117     // cell centres or on agglomeration
 118
 119
 120     autoPtr<decompositionMethod> decomposePtr = decompositionMethod::New
 121     (
 122         decompositionDict_,
 123         *this
 124     );
 125
 126     if (sameProcFaces.size() == 0)
 127     {
 128         cellToProc_ = decomposePtr().decompose(cellCentres());
 129     }
 130     else
 131     {
 132         // Faces where owner and neighbour are not 'connected' (= all except
 133         // sameProcFaces)
 134         boolList blockedFace(nFaces(), true);
 135
 136         forAllConstIter(labelHashSet, sameProcFaces, iter)
 137         {
 138             blockedFace[iter.key()] = false;
 139         }
 140
 141         // Connect coupled boundary faces
 142         const polyBoundaryMesh& patches =  boundaryMesh();
 143
 144         forAll(patches, patchI)
 145         {
 146             const polyPatch& pp = patches[patchI];
 147
 148             if (pp.coupled())
 149             {
 150                 forAll(pp, i)
 151                 {
 152                     blockedFace[pp.start()+i] = false;
 153                 }
 154             }
 155         }
 156
 157         // Determine global regions, separated by blockedFaces
 158         regionSplit globalRegion(*this, blockedFace);
 159
 160
 161         // Determine region cell centres
 162         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 163
 164         // This just takes the first cell in the region. Otherwise the problem
 165         // is with cyclics - if we'd average the region centre might be
 166         // somewhere in the middle of the domain which might not be anywhere
 167         // near any of the cells.
 168
 169         const point greatPoint(GREAT, GREAT, GREAT);
 170
 171         pointField regionCentres(globalRegion.nRegions(), greatPoint);
 172
 173         forAll(globalRegion, cellI)
 174         {
 175             label regionI = globalRegion[cellI];
 176
 177             if (regionCentres[regionI] == greatPoint)
 178             {
 179                 regionCentres[regionI] = cellCentres()[cellI];
 180             }
 181         }
 182
 183         // Do decomposition on agglomeration
 184         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
 185         cellToProc_ = decomposePtr().decompose(globalRegion, regionCentres);
 186     }
 187
 188     Info<< "\nFinished decomposition in "
 189         << decompositionTime.elapsedCpuTime()
 190         << " s" << endl;
 191 }
 192
 193
 194 // ************************************************************************* //