zoltanRenumber.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) 2011-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 Class
25  Foam::zoltanRenumber
26 
27 Description
28  Renumber using Zoltan
29 
30  Zoltan install:
31  - in your ~/.bashrc:
32  export ZOLTAN_ARCH_DIR=\
33  $WM_THIRD_PARTY_DIR/platforms/linux64Gcc/Zoltan_XXX
34  - unpack into $WM_THIRD_PARTY_DIR
35  - cd Zoltan_XXX
36  - mkdir build
37  - cd build
38  - export CCFLAGS="-fPIC"
39  - export CXXFLAGS="-fPIC"
40  - export CFLAGS="-fPIC"
41  - export LDFLAGS="-shared"
42  - ../configure \
43  --prefix=$ZOLTAN_ARCH_DIR \
44  --with-ccflags=-fPIC --with-cxxflags=-fPIC --with-ldflags=-shared
45 
46 SourceFiles
47  zoltanRenumber.C
48 
49 \*---------------------------------------------------------------------------*/
50 
51 #include "zoltanRenumber.H"
53 #include "IFstream.H"
54 #include "labelIOList.H"
55 #include "polyMesh.H"
56 #include "globalMeshData.H"
57 #include "globalIndex.H"
58 #include "uint.H"
59 #include "PstreamGlobals.H"
60 
61 #include "zoltan.h"
62 #include <mpi.h>
63 
64 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
65 
66 namespace Foam
67 {
68  defineTypeNameAndDebug(zoltanRenumber, 0);
69 
71  (
72  renumberMethod,
73  zoltanRenumber,
74  dictionary
75  );
76 }
77 
78 
79 static int get_number_of_vertices(void *data, int *ierr)
80 {
81  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
82  *ierr = ZOLTAN_OK;
83  return mesh.nCells();
84 }
85 
86 
87 static void get_vertex_list(void *data, int sizeGID, int sizeLID,
88  ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID,
89  int wgt_dim, float *obj_wgts, int *ierr)
90 {
91  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
92 
93  *ierr = ZOLTAN_OK;
94 
95  /* In this example, return the IDs of our vertices, but no weights.
96  * Zoltan will assume equally weighted vertices.
97  */
98 
99  wgt_dim = 0;
100  obj_wgts = nullptr;
101 
102  for (Foam::label i=0; i<mesh.nCells(); i++)
103  {
104  globalID[i] = i; // should be global
105  localID[i] = i;
106  }
107 }
108 
110 static void get_num_edges_list(void *data, int sizeGID, int sizeLID,
111  int num_obj,
112  ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID,
113  int *numEdges, int *ierr)
114 {
115  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
116 
117  if ((sizeGID != 1) || (sizeLID != 1) || (num_obj != mesh.nCells()))
118  {
119  *ierr = ZOLTAN_FATAL;
120  return;
121  }
122 
123  for (Foam::label i=0; i < num_obj ;i++)
124  {
125  Foam::label celli = localID[i];
126  const Foam::cell& cFaces = mesh.cells()[celli];
127  forAll(cFaces, cFacei)
128  {
129  Foam::label n = 0;
130  if (mesh.isInternalFace(cFaces[cFacei]))
131  {
132  n++;
133  }
134  numEdges[i] = n;
135  }
136  }
137 
138  *ierr = ZOLTAN_OK;
139 }
140 
142 static void get_edge_list(void *data, int sizeGID, int sizeLID,
143  int num_obj, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID,
144  int *num_edges,
145  ZOLTAN_ID_PTR nborGID, int *nborProc,
146  int wgt_dim, float *ewgts, int *ierr)
147 {
148  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
149 
150  if
151  (
152  (sizeGID != 1)
153  || (sizeLID != 1)
154  || (num_obj != mesh.nCells())
155  || (wgt_dim != 1)
156  )
157  {
158  *ierr = ZOLTAN_FATAL;
159  return;
160  }
161 
162  ZOLTAN_ID_TYPE* nextNbor = nborGID;
163  int* nextProc = nborProc;
164  float* nextWgt = ewgts;
165 
166  for (Foam::label i=0; i < num_obj; i++)
167  {
168  Foam::label celli = localID[i];
169 
170  const Foam::cell& cFaces = mesh.cells()[celli];
171  forAll(cFaces, cFacei)
172  {
173  Foam::label n = 0;
174 
175  Foam::label facei = cFaces[cFacei];
176  if (mesh.isInternalFace(facei))
177  {
178  Foam::label nbr = mesh.faceOwner()[facei];
179  if (nbr == celli)
180  {
181  nbr = mesh.faceNeighbour()[facei];
182  }
183 
184  // Note: global index
185  *nextNbor++ = nbr;
186  *nextProc++ = 0;
187  *nextWgt++ = 1.0;
188 
189  n++;
190  }
191  if (n != num_edges[i])
192  {
193  *ierr = ZOLTAN_FATAL;
194  return;
195  }
196  }
197  }
198  *ierr = ZOLTAN_OK;
199 }
200 
202 static int get_mesh_dim(void *data, int *ierr)
203 {
204  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
205 
206  return mesh.nSolutionD();
207 }
208 
209 
210 static void get_geom_list
211 (
212  void *data,
213  int num_gid_entries,
214  int num_lid_entries,
215  int num_obj,
216  ZOLTAN_ID_PTR global_ids,
217  ZOLTAN_ID_PTR local_ids,
218  int num_dim,
219  double *geom_vec,
220  int *ierr
221 )
222 {
223  const Foam::polyMesh& mesh = *static_cast<const Foam::polyMesh*>(data);
224 
225  if
226  (
227  (num_gid_entries != 1)
228  || (num_lid_entries != 1)
229  || (num_obj != mesh.nCells())
230  || (num_dim != mesh.nSolutionD())
231  )
232  {
233  *ierr = ZOLTAN_FATAL;
234  return;
235  }
236 
237  double* p = geom_vec;
238 
239 
240  const Foam::Vector<Foam::label>& sol = mesh.solutionD();
241 
242  const Foam::pointField& cc = mesh.cellCentres();
243 
244  for (Foam::label celli = 0; celli < num_obj; celli++)
245  {
246  const Foam::point& pt = cc[celli];
247 
248  for (Foam::direction cmpt = 0; cmpt < Foam::vector::nComponents; cmpt++)
249  {
250  if (sol[cmpt] == 1)
251  {
252  *p++ = pt[cmpt];
253  }
254  }
255  }
256  *ierr = ZOLTAN_OK;
257 }
258 
259 
260 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
263 :
264  renumberMethod(renumberDict),
265  coeffsDict_(renumberDict.optionalSubDict(typeName+"Coeffs"))
266 {}
267 
268 
269 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
270 
272 (
273  const polyMesh& pMesh,
274  const pointField& points
275 ) const
276 {
277  stringList args(1);
278  args[0] = "zoltanRenumber";
279 
280  int argc = args.size();
281  char* argv[argc];
282  for (label i = 0; i < argc; i++)
283  {
284  argv[i] = strdup(args[i].c_str());
285  }
286 
287  float ver;
288  int rc = Zoltan_Initialize(argc, argv, &ver);
289 
290  Foam::Pout<< "Initialised to " << ver << Foam::endl;
291 
292  if (rc != ZOLTAN_OK)
293  {
295  << "Failed initialising Zoltan" << exit(FatalError);
296  }
297 
298  struct Zoltan_Struct *zz = Zoltan_Create(PstreamGlobals::MPI_COMM_FOAM);
299 
300  polyMesh& mesh = const_cast<polyMesh&>(pMesh);
301 
302 
303  forAllConstIter(IDLList<entry>, coeffsDict_, iter)
304  {
305  if (!iter().isDict())
306  {
307  const word& key = iter().keyword();
308  const word value(iter().stream());
309 
310  Info<< typeName << " : setting parameter " << key
311  << " to " << value << endl;
312 
313  Zoltan_Set_Param(zz, key.c_str(), value.c_str());
314  }
315  }
316 
317 
318  Zoltan_Set_Num_Obj_Fn(zz, get_number_of_vertices, &mesh);
319  Zoltan_Set_Obj_List_Fn(zz, get_vertex_list, &mesh);
320 
321  // Callbacks for geometry
322  Zoltan_Set_Num_Geom_Fn(zz, get_mesh_dim, &mesh);
323  Zoltan_Set_Geom_Multi_Fn(zz, get_geom_list, &mesh);
324 
325  // Callbacks for connectivity
326  Zoltan_Set_Num_Edges_Multi_Fn(zz, get_num_edges_list, &mesh);
327  Zoltan_Set_Edge_List_Multi_Fn(zz, get_edge_list, &mesh);
328 
329 
330 
331  // Note: !global indices
332  List<ZOLTAN_ID_TYPE> wantedCells(mesh.nCells());
333 
334  globalIndex globalCells(mesh.nCells());
335  forAll(wantedCells, i)
336  {
337  // wantedCells[i] = i;
338  wantedCells[i] = globalCells.toGlobal(i);
339  }
340 
341  List<ZOLTAN_ID_TYPE> oldToNew(mesh.nCells());
342 
343  int err = Zoltan_Order
344  (
345  zz,
346  1, // int num_gid_entries,
347  mesh.globalData().nTotalCells(), // int num_obj,
348  wantedCells.begin(),
349  oldToNew.begin()
350  );
351 
352  if (err != ZOLTAN_OK)
353  {
355  << "Failed Zoltan_Order" << exit(FatalError);
356  }
357 
358 
359  for (label i = 0; i < argc; i++)
360  {
361  free(argv[i]);
362  }
363 
364 
365  labelList order(oldToNew.size());
366  forAll(order, i)
367  {
368  order[i] = oldToNew[i];
369  }
370  return order;
371 }
372 
373 
374 // ************************************************************************* //
Template class for intrusive linked lists.
Definition: ILList.H:50
static int get_mesh_dim(void *data, int *ierr)
#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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1175
uint8_t direction
Definition: direction.H:45
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
Definition: polyMesh.C:831
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
label nCells() const
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
static void get_geom_list(void *data, int num_gid_entries, int num_lid_entries, int num_obj, ZOLTAN_ID_PTR global_ids, ZOLTAN_ID_PTR local_ids, int num_dim, double *geom_vec, int *ierr)
label nTotalCells() const
Return total number of cells in decomposed mesh.
const cellList & cells() const
virtual labelList renumber(const pointField &) const
Return the order in which cells need to be visited, i.e.
Macros for easy insertion into run-time selection tables.
Abstract base class for renumbering.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
static void get_edge_list(void *data, int sizeGID, int sizeLID, int num_obj, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, int *num_edges, ZOLTAN_ID_PTR nborGID, int *nborProc, int wgt_dim, float *ewgts, int *ierr)
static const direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:99
dynamicFvMesh & mesh
System uinteger.
A class for handling words, derived from string.
Definition: word.H:59
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1169
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
const vectorField & cellCentres() const
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
bool isInternalFace(const label faceIndex) const
Return true if given face label is internal to the mesh.
Templated 3D Vector derived from VectorSpace adding construction from 3 components, element access using x(), y() and z() member functions and the inner-product (dot-product) and cross product operators.
Definition: Vector.H:57
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
defineTypeNameAndDebug(combustionModel, 0)
static int get_number_of_vertices(void *data, int *ierr)
label nSolutionD() const
Return the number of valid solved-for dimensions in the mesh.
Definition: polyMesh.C:842
static void get_vertex_list(void *data, int sizeGID, int sizeLID, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, int wgt_dim, float *obj_wgts, int *ierr)
zoltanRenumber(const dictionary &renumberDict)
Construct given the renumber dictionary.
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
messageStream Info
label n
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
volScalarField & p
Foam::argList args(argc, argv)
Namespace for OpenFOAM.
static void get_num_edges_list(void *data, int sizeGID, int sizeLID, int num_obj, ZOLTAN_ID_PTR globalID, ZOLTAN_ID_PTR localID, int *numEdges, int *ierr)