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