domainDecompositionDryRun.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) 2018-2023 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
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 
29 #include "volFields.H"
30 #include "globalMeshData.H"
31 #include "decompositionModel.H"
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
37 (
38  const IOobject& io,
39  const fileName& decompDictFile,
40  const label nDomains,
41  const word& methodName
42 )
43 :
44  mesh_(io),
45  decompDictFile_(decompDictFile),
46  nDomainsOverride_(nDomains),
47  methodNameOverride_(methodName)
48 {}
49 
50 
51 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
52 
54 (
55  const bool writeCellDist,
56  const bool verbose
57 )
58 {
59  cpuTime decompositionTime;
60 
61  Info<< "\nCalculating distribution of cells. nCells = "
62  << mesh_.nCells() << endl;
63 
64  const decompositionModel& model = decompositionModel::New
65  (
66  mesh_,
67  decompDictFile_
68  );
69 
70  // Allow overrides for testing
71 
72  dictionary& modelDict = const_cast<decompositionModel&>(model);
73 
74  if (nDomainsOverride_ > 0)
75  {
76  modelDict.add
77  (
78  word("numberOfSubdomains"),
79  nDomainsOverride_,
80  true
81  );
82  }
83 
84  if (!methodNameOverride_.empty())
85  {
86  modelDict.add
87  (
88  word("method"),
89  methodNameOverride_,
90  true
91  );
92  }
93 
94  scalarField cellWeights;
95  word weightName;
96  if (model.readIfPresent("weightField", weightName))
97  {
98  volScalarField weights
99  (
100  IOobject
101  (
102  weightName,
103  mesh_.time().timeName(),
104  mesh_,
107  ),
108  mesh_
109  );
110  cellWeights = weights.primitiveField();
111  }
112 
113  decompositionMethod& method = model.decomposer();
114 
115  CompactListList<label> cellCells;
117  (
118  mesh_,
119  identity(mesh_.nCells()),
120  mesh_.nCells(),
121  false, // false = local only
122  cellCells
123  );
124 
125  labelList cellToProc = method.decompose(mesh_, cellWeights);
126 
127  Info<< "\nFinished decomposition into "
128  << method.nDomains() << " domains in "
129  << decompositionTime.elapsedCpuTime() << " s" << nl << nl;
130 
131  decompositionInformation info
132  (
133  cellCells,
134  cellToProc,
135  method.nDomains()
136  );
137 
138  if (writeCellDist)
139  {
140  // Write decomposition for visualization
141  // - write as VTU to avoid any impact
142  writeVTK("cellDist", cellToProc);
143 
144 // Less likely that this is actually required, but may be useful...
145 //
146 // // Write decomposition as labelList for use with 'manual'
147 // // decomposition method.
148 // labelIOList cellDecomposition
149 // (
150 // IOobject
151 // (
152 // "cellDecomposition",
153 // mesh_.facesInstance(),
154 // mesh_,
155 // IOobject::NO_READ,
156 // IOobject::NO_WRITE,
157 // IOobject::NO_REGISTER
158 // ),
159 // std::move(cellToProc)
160 // );
161 // cellDecomposition.write();
162 //
163 // Info<< nl << "Wrote decomposition to "
164 // << cellDecomposition.objectRelPath()
165 // << " for use in manual decomposition." << endl;
166 
167  Info<< nl;
168  }
169 
170  if (verbose)
171  {
172  info.printDetails(Info);
173  Info<< nl;
174  }
175  info.printSummary(Info);
176 }
177 
178 
179 // ************************************************************************* //
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
static void calcCellCells(const polyMesh &mesh, const labelList &agglom, const label nLocalCoarse, const bool parallel, CompactListList< label > &cellCells)
Determine (local or global) cellCells from mesh agglomeration.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:531
Ignore writing from objectRegistry::writeObject()
static const decompositionModel & New(const polyMesh &mesh, const fileName &decompDictFile="", const dictionary *fallback=nullptr)
Read and register on mesh, optionally with alternative decomposeParDict path/name or with fallback co...
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:81
void execute(const bool writeCellDist, const bool verbose=false)
Perform dry-run.
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i), works like std::iota() but returning a...
Definition: labelLists.C:44
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
cpuTimePosix cpuTime
Selection of preferred clock mechanism for the elapsed cpu time.
Definition: cpuTimeFwd.H:36
messageStream Info
Information stream (stdout output on master, null elsewhere)
void writeVTK(OFstream &os, const Type &value)
List< label > labelList
A List of labels.
Definition: List.H:62
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, IOobject::NO_REGISTER)
domainDecompositionDryRun(const IOobject &io, const fileName &decompDictFile="", const label nDomains=0, const word &methodName="")
Construct from components.