procFacesGAMGProcAgglomeration.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2013-2019 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 
28 #include "GAMGAgglomeration.H"
29 #include "Random.H"
30 #include "lduMesh.H"
31 #include "processorLduInterface.H"
32 #include "processorGAMGInterface.H"
33 #include "pairGAMGAgglomeration.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(procFacesGAMGProcAgglomeration, 0);
40 
42  (
43  GAMGProcAgglomeration,
44  procFacesGAMGProcAgglomeration,
45  GAMGAgglomeration
46  );
47 }
48 
49 
50 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
51 
53 Foam::procFacesGAMGProcAgglomeration::singleCellMesh
54 (
55  const label singleCellMeshComm,
56  const lduMesh& mesh,
57  scalarField& faceWeights
58 ) const
59 {
60  // Count number of faces per processor
61  List<Map<label>> procFaces(UPstream::nProcs(mesh.comm()));
62  Map<label>& myNeighbours = procFaces[UPstream::myProcNo(mesh.comm())];
63 
64  {
65  const lduInterfacePtrsList interfaces(mesh.interfaces());
66  forAll(interfaces, intI)
67  {
68  if (interfaces.set(intI))
69  {
70  const processorLduInterface& pp =
71  refCast<const processorLduInterface>
72  (
73  interfaces[intI]
74  );
75  label size = interfaces[intI].faceCells().size();
76  myNeighbours.insert(pp.neighbProcNo(), size);
77  }
78  }
79  }
80 
81  Pstream::gatherList(procFaces, Pstream::msgType(), mesh.comm());
82  Pstream::scatterList(procFaces, Pstream::msgType(), mesh.comm());
83 
84  autoPtr<lduPrimitiveMesh> singleCellMeshPtr;
85 
86  if (Pstream::master(mesh.comm()))
87  {
88  // I am master
89  label nCells = Pstream::nProcs(mesh.comm());
90 
91  DynamicList<label> l(3*nCells);
92  DynamicList<label> u(3*nCells);
93  DynamicList<scalar> weight(3*nCells);
94 
95  DynamicList<label> nbrs;
96  DynamicList<scalar> weights;
97 
98  forAll(procFaces, proci)
99  {
100  const Map<label>& neighbours = procFaces[proci];
101 
102  // Add all the higher processors
103  nbrs.clear();
104  weights.clear();
105  forAllConstIter(Map<label>, neighbours, iter)
106  {
107  if (iter.key() > proci)
108  {
109  nbrs.append(iter.key());
110  weights.append(iter());
111  }
112  sort(nbrs);
113  forAll(nbrs, i)
114  {
115  l.append(proci);
116  u.append(nbrs[i]);
117  weight.append(weights[i]);
118  }
119  }
120  }
121 
122  faceWeights.transfer(weight);
123 
124  PtrList<const lduInterface> primitiveInterfaces(0);
125  const lduSchedule ps(0);
126 
127  singleCellMeshPtr.reset
128  (
129  new lduPrimitiveMesh
130  (
131  nCells,
132  l,
133  u,
134  primitiveInterfaces,
135  ps,
136  singleCellMeshComm
137  )
138  );
139  }
140  return singleCellMeshPtr;
141 }
142 
143 
145 Foam::procFacesGAMGProcAgglomeration::processorAgglomeration
146 (
147  const lduMesh& mesh
148 ) const
149 {
150  label singleCellMeshComm = UPstream::allocateCommunicator
151  (
152  mesh.comm(),
153  labelList(1, label(0)) // only processor 0
154  );
155 
156  scalarField faceWeights;
157  autoPtr<lduPrimitiveMesh> singleCellMeshPtr
158  (
159  singleCellMesh
160  (
161  singleCellMeshComm,
162  mesh,
163  faceWeights
164  )
165  );
166 
167  tmp<labelField> tfineToCoarse(new labelField(0));
168  labelField& fineToCoarse = tfineToCoarse.ref();
169 
170  if (singleCellMeshPtr.valid())
171  {
172  // On master call the agglomerator
173  const lduPrimitiveMesh& singleCellMesh = singleCellMeshPtr();
174 
175  label nCoarseProcs;
177  (
178  nCoarseProcs,
179  singleCellMesh,
180  faceWeights
181  );
182 
183  labelList coarseToMaster(nCoarseProcs, labelMax);
184  forAll(fineToCoarse, celli)
185  {
186  label coarseI = fineToCoarse[celli];
187  coarseToMaster[coarseI] = min(coarseToMaster[coarseI], celli);
188  }
189 
190  // Sort according to master and redo restriction
191  labelList newToOld;
192  sortedOrder(coarseToMaster, newToOld);
193  labelList oldToNew(invert(newToOld.size(), newToOld));
194 
195  fineToCoarse = UIndirectList<label>(oldToNew, fineToCoarse)();
196  }
197 
198  Pstream::scatter(fineToCoarse, Pstream::msgType(), mesh.comm());
199  UPstream::freeCommunicator(singleCellMeshComm);
200 
201  return tfineToCoarse;
202 }
203 
204 
205 bool Foam::procFacesGAMGProcAgglomeration::doProcessorAgglomeration
206 (
207  const lduMesh& mesh
208 ) const
209 {
210  // Check the need for further agglomeration on all processors
211  bool doAgg = mesh.lduAddr().size() < nAgglomeratingCells_;
212  mesh.reduce(doAgg, orOp<bool>());
213  return doAgg;
214 }
215 
216 
217 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
218 
220 (
221  GAMGAgglomeration& agglom,
222  const dictionary& controlDict
223 )
224 :
225  GAMGProcAgglomeration(agglom, controlDict),
226  nAgglomeratingCells_(controlDict.lookup<label>("nAgglomeratingCells"))
227 {}
228 
229 
230 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
231 
233 {
234  forAllReverse(comms_, i)
235  {
236  if (comms_[i] != -1)
237  {
238  UPstream::freeCommunicator(comms_[i]);
239  }
240  }
241 }
242 
243 
244 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
245 
247 {
248  if (debug)
249  {
250  Pout<< nl << "Starting mesh overview" << endl;
251  printStats(Pout, agglom_);
252  }
253 
254  if (agglom_.size() >= 1)
255  {
256  Random rndGen(0);
257 
258  for
259  (
260  label fineLevelIndex = 2;
261  fineLevelIndex < agglom_.size();
262  fineLevelIndex++
263  )
264  {
265  if (agglom_.hasMeshLevel(fineLevelIndex))
266  {
267  // Get the fine mesh
268  const lduMesh& levelMesh = agglom_.meshLevel(fineLevelIndex);
269 
270  label levelComm = levelMesh.comm();
271  label nProcs = UPstream::nProcs(levelComm);
272 
273  if (nProcs > 1 && doProcessorAgglomeration(levelMesh))
274  {
275  tmp<labelField> tprocAgglomMap
276  (
277  processorAgglomeration(levelMesh)
278  );
279  const labelField& procAgglomMap = tprocAgglomMap();
280 
281  // Master processor
282  labelList masterProcs;
283  // Local processors that agglomerate. agglomProcIDs[0] is in
284  // masterProc.
285  List<label> agglomProcIDs;
287  (
288  levelComm,
289  procAgglomMap,
290  masterProcs,
291  agglomProcIDs
292  );
293 
294  // Allocate a communicator for the processor-agglomerated
295  // matrix
296  comms_.append
297  (
299  (
300  levelComm,
301  masterProcs
302  )
303  );
304 
305 
306  // Use processor agglomeration maps to do the actual
307  // collecting.
309  (
310  fineLevelIndex,
311  procAgglomMap,
312  masterProcs,
313  agglomProcIDs,
314  comms_.last()
315  );
316  }
317  }
318  }
319  }
320 
321  // Print a bit
322  if (debug)
323  {
324  Pout<< nl << "Agglomerated mesh overview" << endl;
325  printStats(Pout, agglom_);
326  }
327 
328  return true;
329 }
330 
331 
332 // ************************************************************************* //
static void scatterList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Scatter data. Reverse of gatherList.
Field< label > labelField
Specialisation of Field<T> for label.
Definition: labelField.H:49
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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 sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
procFacesGAMGProcAgglomeration(GAMGAgglomeration &agglom, const dictionary &controlDict)
Construct given agglomerator and controls.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:59
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: UList.H:446
List< lduScheduleEntry > lduSchedule
Definition: lduSchedule.H:74
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
virtual label comm() const =0
Return communicator used for parallel communication.
Macros for easy insertion into run-time selection tables.
Random rndGen(label(0))
void agglomerate(const lduMesh &mesh, const scalarField &faceWeights)
Agglomerate all levels starting from the given face weights.
labelList invert(const label len, const labelUList &)
Invert one-to-one map. Unmapped elements will be -1.
Definition: ListOps.C:37
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
void sort(UList< T > &)
Definition: UList.C:115
UPtrList< const lduInterface > lduInterfacePtrsList
List of coupled interface fields to be used in coupling.
static const label labelMax
Definition: label.H:62
List< label > labelList
A List of labels.
Definition: labelList.H:56
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Random number generator.
Definition: Random.H:57
static void scatter(const List< commsStruct > &comms, T &Value, const int tag, const label comm)
Scatter data. Distribute without modification. Reverse of gather.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool agglomerate()
Modify agglomeration. Return true if modified.
virtual bool agglomerate()=0
Modify agglomeration. Return true if modified.
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
Processor agglomeration of GAMGAgglomerations.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
A class for managing temporary objects.
Definition: PtrList.H:53
static label allocateCommunicator(const label parent, const labelList &subRanks, const bool doPstream=true)
Allocate a new communicator.
Definition: UPstream.C:248
T & last()
Return the last element of the list.
Definition: UListI.H:128
Geometric agglomerated algebraic multigrid agglomeration class.
static void gatherList(const List< commsStruct > &comms, List< T > &Values, const int tag, const label comm)
Gather data but keep individual values separate.
static void calculateRegionMaster(const label comm, const labelList &procAgglomMap, labelList &masterProcs, List< label > &agglomProcIDs)
Given fine to coarse processor map determine:
static void freeCommunicator(const label communicator, const bool doPstream=true)
Free a previously allocated communicator.
Definition: UPstream.C:314
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844