GAMGProcAgglomeration.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) 2013-2014 OpenFOAM Foundation
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 "GAMGProcAgglomeration.H"
27 #include "GAMGAgglomeration.H"
28 #include "lduMesh.H"
29 
30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34  defineTypeNameAndDebug(GAMGProcAgglomeration, 0);
35  defineRunTimeSelectionTable(GAMGProcAgglomeration, GAMGAgglomeration);
36 }
37 
38 
39 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
40 
42 (
43  Ostream& os,
44  GAMGAgglomeration& agglom
45 ) const
46 {
47  for (label levelI = 0; levelI <= agglom.size(); levelI++)
48  {
49  if (agglom.hasMeshLevel(levelI))
50  {
51  os << agglom.meshLevel(levelI).info() << endl;
52  }
53  else
54  {
55  os << "Level " << levelI << " has no fine mesh:" << endl;
56  }
57 
58  if
59  (
60  levelI < agglom.restrictAddressing_.size()
61  && agglom.restrictAddressing_.set(levelI)
62  )
63  {
64  const labelList& cellRestrict =
65  agglom.restrictAddressing(levelI);
66  const labelList& faceRestrict =
67  agglom.faceRestrictAddressing(levelI);
68 
69  os << "Level " << levelI << " agglomeration:" << nl
70  << " nCoarseCells:" << agglom.nCells(levelI) << nl
71  << " nCoarseFaces:" << agglom.nFaces(levelI) << nl
72  << " cellRestriction:"
73  << " size:" << cellRestrict.size()
74  << " max:" << max(cellRestrict)
75  << nl
76  << " faceRestriction:"
77  << " size:" << faceRestrict.size()
78  << " max:" << max(faceRestrict)
79  << nl;
80 
81  const labelListList& patchFaceRestrict =
82  agglom.patchFaceRestrictAddressing(levelI);
83  forAll(patchFaceRestrict, i)
84  {
85  if (patchFaceRestrict[i].size())
86  {
87  const labelList& faceRestrict =
88  patchFaceRestrict[i];
89  os << " " << i
90  << " size:" << faceRestrict.size()
91  << " max:" << max(faceRestrict)
92  << nl;
93  }
94  }
95  }
96  if
97  (
98  levelI < agglom.procCellOffsets_.size()
99  && agglom.procCellOffsets_.set(levelI)
100  )
101  {
102  os << " procCellOffsets:" << agglom.procCellOffsets_[levelI]
103  << nl
104  << " procAgglomMap:" << agglom.procAgglomMap_[levelI]
105  << nl
106  << " procIDs:" << agglom.agglomProcIDs_[levelI]
107  << nl
108  << " comm:" << agglom.procCommunicator_[levelI]
109  << endl;
110  }
111 
112  os << endl;
113  }
114  os << endl;
115 }
116 
117 
119 (
120  const lduMesh& mesh
121 )
122 {
123  const lduAddressing& addr = mesh.lduAddr();
124  lduInterfacePtrsList interfaces = mesh.interfaces();
125 
126  const label myProcID = Pstream::myProcNo(mesh.comm());
127 
128  globalIndex globalNumbering
129  (
130  addr.size(),
132  mesh.comm(),
134  );
135 
136  labelList globalIndices(addr.size());
137  forAll(globalIndices, cellI)
138  {
139  globalIndices[cellI] = globalNumbering.toGlobal(myProcID, cellI);
140  }
141 
142 
143  // Get the interface cells
144  PtrList<labelList> nbrGlobalCells(interfaces.size());
145  {
146  // Initialise transfer of restrict addressing on the interface
147  forAll(interfaces, inti)
148  {
149  if (interfaces.set(inti))
150  {
151  interfaces[inti].initInternalFieldTransfer
152  (
154  globalIndices
155  );
156  }
157  }
158 
159  if (Pstream::parRun())
160  {
162  }
163 
164  forAll(interfaces, inti)
165  {
166  if (interfaces.set(inti))
167  {
168  nbrGlobalCells.set
169  (
170  inti,
171  new labelList
172  (
173  interfaces[inti].internalFieldTransfer
174  (
176  globalIndices
177  )
178  )
179  );
180  }
181  }
182  }
183 
184 
185  // Scan the neighbour list to find out how many times the cell
186  // appears as a neighbour of the face. Done this way to avoid guessing
187  // and resizing list
188  labelList nNbrs(addr.size(), 1);
189 
190  const labelUList& nbr = addr.upperAddr();
191  const labelUList& own = addr.lowerAddr();
192 
193  {
194  forAll(nbr, faceI)
195  {
196  nNbrs[nbr[faceI]]++;
197  nNbrs[own[faceI]]++;
198  }
199 
200  forAll(interfaces, inti)
201  {
202  if (interfaces.set(inti))
203  {
204  const labelUList& faceCells = interfaces[inti].faceCells();
205 
206  forAll(faceCells, i)
207  {
208  nNbrs[faceCells[i]]++;
209  }
210  }
211  }
212  }
213 
214 
215  // Create cell-cells addressing
216  labelListList cellCells(addr.size());
217 
218  forAll(cellCells, cellI)
219  {
220  cellCells[cellI].setSize(nNbrs[cellI], -1);
221  }
222 
223  // Reset the list of number of neighbours to zero
224  nNbrs = 0;
225 
226  // Scatter the neighbour faces
227  forAll(nbr, faceI)
228  {
229  label c0 = own[faceI];
230  label c1 = nbr[faceI];
231 
232  cellCells[c0][nNbrs[c0]++] = globalIndices[c1];
233  cellCells[c1][nNbrs[c1]++] = globalIndices[c0];
234  }
235  forAll(interfaces, inti)
236  {
237  if (interfaces.set(inti))
238  {
239  const labelUList& faceCells = interfaces[inti].faceCells();
240 
241  forAll(faceCells, i)
242  {
243  label c0 = faceCells[i];
244  cellCells[c0][nNbrs[c0]++] = nbrGlobalCells[inti][i];
245  }
246  }
247  }
248 
249  forAll(cellCells, cellI)
250  {
251  Foam::stableSort(cellCells[cellI]);
252  }
253 
254  // Replace the initial element (always -1) with the local cell
255  forAll(cellCells, cellI)
256  {
257  cellCells[cellI][0] = globalIndices[cellI];
258  }
259 
260  return cellCells;
261 }
262 
263 
265 (
266  const label fineLevelIndex,
267  const labelList& procAgglomMap,
268  const labelList& masterProcs,
269  const List<label>& agglomProcIDs,
270  const label procAgglomComm
271 )
272 {
273  const lduMesh& levelMesh = agglom_.meshLevels_[fineLevelIndex];
274  label levelComm = levelMesh.comm();
275 
276  if (Pstream::myProcNo(levelComm) != -1)
277  {
278  // Collect meshes and restrictAddressing onto master
279  // Overwrites the fine mesh (meshLevels_[index-1]) and addressing
280  // from fine mesh to coarse mesh (restrictAddressing_[index]).
281  agglom_.procAgglomerateLduAddressing
282  (
283  levelComm,
284  procAgglomMap,
285  agglomProcIDs,
286  procAgglomComm,
287 
288  fineLevelIndex //fine level index
289  );
290 
291  // Combine restrict addressing only onto master
292  for
293  (
294  label levelI = fineLevelIndex+1;
295  levelI < agglom_.meshLevels_.size();
296  levelI++
297  )
298  {
299  agglom_.procAgglomerateRestrictAddressing
300  (
301  levelComm,
302  agglomProcIDs,
303  levelI
304  );
305  }
306 
307  if (Pstream::myProcNo(levelComm) == agglomProcIDs[0])
308  {
309  // On master. Recreate coarse meshes from restrict addressing
310  for
311  (
312  label levelI = fineLevelIndex;
313  levelI < agglom_.meshLevels_.size();
314  levelI++
315  )
316  {
317  agglom_.agglomerateLduAddressing(levelI);
318  }
319  }
320  else
321  {
322  // Agglomerated away. Clear mesh storage.
323  for
324  (
325  label levelI = fineLevelIndex+1;
326  levelI <= agglom_.size();
327  levelI++
328  )
329  {
330  agglom_.clearLevel(levelI);
331  }
332  }
333  }
334 
335  // Should check!
336  return true;
337 }
338 
339 
340 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
341 
342 Foam::GAMGProcAgglomeration::GAMGProcAgglomeration
343 (
344  GAMGAgglomeration& agglom,
345  const dictionary& controlDict
346 )
347 :
348  agglom_(agglom)
349 {}
350 
351 
353 (
354  const word& type,
355  GAMGAgglomeration& agglom,
356  const dictionary& controlDict
357 )
358 {
359  if (debug)
360  {
361  Info<< "GAMGProcAgglomeration::New(const word&, GAMGAgglomeration&"
362  ", const dictionary&) : "
363  "constructing GAMGProcAgglomeration"
364  << endl;
365  }
366 
367  GAMGAgglomerationConstructorTable::iterator cstrIter =
368  GAMGAgglomerationConstructorTablePtr_->find(type);
369 
370  if (cstrIter == GAMGAgglomerationConstructorTablePtr_->end())
371  {
373  (
374  "GAMGProcAgglomeration::New(const word&, GAMGAgglomeration&"
375  ", const dictionary&) "
376  ) << "Unknown GAMGProcAgglomeration type "
377  << type << " for GAMGAgglomeration " << agglom.type() << nl << nl
378  << "Valid GAMGProcAgglomeration types are :" << endl
379  << GAMGAgglomerationConstructorTablePtr_->sortedToc()
380  << exit(FatalError);
381  }
382 
383  return autoPtr<GAMGProcAgglomeration>(cstrIter()(agglom, controlDict));
384 }
385 
386 
387 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
388 
390 {}
391 
392 
393 // ************************************************************************* //
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:380
label size() const
Return number of equations.
virtual const labelUList & upperAddr() const =0
Return upper addressing.
void printStats(Ostream &os, GAMGAgglomeration &agglom) const
Debug: write agglomeration info.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
virtual label comm() const =0
Return communicator used for parallel communication.
Geometric agglomerated algebraic multigrid agglomeration class.
static autoPtr< GAMGProcAgglomeration > New(const word &type, GAMGAgglomeration &agglom, const dictionary &controlDict)
Return the selected agglomerator.
virtual lduInterfacePtrsList interfaces() const =0
Return a list of pointers for each patch.
labelList procCommunicator_
Communicator for given level.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const labelList & faceRestrictAddressing(const label leveli) const
Return face restrict addressing of given level.
messageStream Info
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:51
PtrList< labelList > procCellOffsets_
Mapping from processor to procMeshLevel cells.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
InfoProxy< lduMesh > info() const
Return info proxy.
Definition: lduMesh.H:96
PtrList< labelList > procAgglomMap_
Per level, per processor the processor it agglomerates into.
virtual bool agglomerate()=0
Modify agglomeration. Return true if modified.
virtual ~GAMGProcAgglomeration()
Destructor.
Namespace for OpenFOAM.
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
static void waitRequests(const label start=0)
Wait until all requests (from start onwards) have finished.
Definition: UPstream.C:106
const lduMesh & meshLevel(const label leveli) const
Return LDU mesh of given level.
static const char nl
Definition: Ostream.H:260
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool hasMeshLevel(const label leveli) const
Do we have mesh for given level?
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
static labelListList globalCellCells(const lduMesh &)
Debug: calculate global cell-cells.
PtrList< labelList > agglomProcIDs_
Per level the set of processors to agglomerate. Element 0 is.
#define forAll(list, i)
Definition: UList.H:421
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:451
virtual const labelUList & lowerAddr() const =0
Return lower addressing.
label nFaces(const label leveli) const
Return number of coarse faces (before processor agglomeration)
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:404
bool set(const label) const
Is element set.
Definition: UPtrListI.H:80
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
const labelField & restrictAddressing(const label leveli) const
Return cell restrict addressing of given level.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:31
error FatalError
runTime controlDict().lookup("adjustTimeStep") >> adjustTimeStep
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
label nCells(const label leveli) const
Return number of coarse cells (before processor agglomeration)
defineRunTimeSelectionTable(reactionRateFlameArea, dictionary)
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
void stableSort(UList< T > &)
Definition: UList.C:121
const labelListList & patchFaceRestrictAddressing(const label leveli) const
The class contains the addressing required by the lduMatrix: upper, lower and losort.
PtrList< labelField > restrictAddressing_
Cell restriction addressing array.
defineTypeNameAndDebug(combustionModel, 0)