MGridGenGAMGAgglomerate.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) 2011-2013 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 Description
25  Agglomerate one level using the MGridGen algorithm.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 #include "fvMesh.H"
31 #include "syncTools.H"
32 
33 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
34 
35 void Foam::MGridGenGAMGAgglomeration::
36 makeCompactCellFaceAddressingAndFaceWeights
37 (
38  const lduAddressing& fineAddressing,
39  List<idxtype>& cellCells,
40  List<idxtype>& cellCellOffsets,
41  const scalarField& magSi,
42  List<scalar>& faceWeights
43 )
44 {
45  const label nFineCells = fineAddressing.size();
46  const label nFineFaces = fineAddressing.upperAddr().size();
47 
48  const labelUList& upperAddr = fineAddressing.upperAddr();
49  const labelUList& lowerAddr = fineAddressing.lowerAddr();
50 
51  // Number of neighbours for each cell
52  labelList nNbrs(nFineCells, 0);
53 
54  forAll(upperAddr, facei)
55  {
56  nNbrs[upperAddr[facei]]++;
57  }
58 
59  forAll(lowerAddr, facei)
60  {
61  nNbrs[lowerAddr[facei]]++;
62  }
63 
64  // Set the sizes of the addressing and faceWeights arrays
65  cellCellOffsets.setSize(nFineCells + 1);
66  cellCells.setSize(2*nFineFaces);
67  faceWeights.setSize(2*nFineFaces);
68 
69 
70  cellCellOffsets[0] = 0;
71  forAll(nNbrs, celli)
72  {
73  cellCellOffsets[celli+1] = cellCellOffsets[celli] + nNbrs[celli];
74  }
75 
76  // reset the whole list to use as counter
77  nNbrs = 0;
78 
79  forAll(upperAddr, facei)
80  {
81  label own = upperAddr[facei];
82  label nei = lowerAddr[facei];
83 
84  label l1 = cellCellOffsets[own] + nNbrs[own]++;
85  label l2 = cellCellOffsets[nei] + nNbrs[nei]++;
86 
87  cellCells[l1] = nei;
88  cellCells[l2] = own;
89 
90  faceWeights[l1] = magSi[facei];
91  faceWeights[l2] = magSi[facei];
92  }
93 }
94 
95 
96 Foam::tmp<Foam::labelField> Foam::MGridGenGAMGAgglomeration::agglomerate
97 (
98  label& nCoarseCells,
99  const label minSize,
100  const label maxSize,
101  const lduAddressing& fineAddressing,
102  const scalarField& V,
103  const scalarField& magSf,
104  const scalarField& magSb
105 )
106 {
107  const label nFineCells = fineAddressing.size();
108 
109  // Compact addressing for cellCells
110  List<idxtype> cellCells;
111  List<idxtype> cellCellOffsets;
112 
113  // Face weights = face areas of the internal faces
114  List<scalar> faceWeights;
115 
116  // Create the compact addressing for cellCells and faceWeights
117  makeCompactCellFaceAddressingAndFaceWeights
118  (
119  fineAddressing,
120  cellCells,
121  cellCellOffsets,
122  magSf,
123  faceWeights
124  );
125 
126  // agglomeration options.
127  List<int> options(4, 0);
128  options[0] = 4; // globular agglom
129  options[1] = 6; // objective F3 and F2
130  options[2] = 128; // debugging output level
131  options[3] = fvMesh_.nGeometricD(); // Dimensionality of the grid
132 
133 
134  // output: cell -> processor addressing
135  List<int> finalAgglom(nFineCells);
136  int nMoves = -1;
137 
138  MGridGen
139  (
140  nFineCells,
141  cellCellOffsets.begin(),
142  const_cast<scalar*>(V.begin()),
143  const_cast<scalar*>(magSb.begin()),
144  cellCells.begin(),
145  faceWeights.begin(),
146  minSize,
147  maxSize,
148  options.begin(),
149  &nMoves,
150  &nCoarseCells,
151  finalAgglom.begin()
152  );
153 
154  {
155  label nNewCoarseCells = 0;
156  labelList newRestrictAddr;
157  bool ok = checkRestriction
158  (
159  newRestrictAddr,
160  nNewCoarseCells,
161  fineAddressing,
162  finalAgglom,
163  nCoarseCells
164  );
165 
166  if (!ok)
167  {
168  nCoarseCells = nNewCoarseCells;
169  finalAgglom.transfer(newRestrictAddr);
170  }
171  }
172 
173  return tmp<labelField>(new labelField(finalAgglom));
174 }
175 
176 
177 // ************************************************************************* //
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:428
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 MGridGen(int, int *, Foam::scalar *, Foam::scalar *, int *, Foam::scalar *, int, int, int *, int *, int *, int *)
Definition: dummyMGridGen.C:49
static bool checkRestriction(labelList &newRestrict, label &nNewCoarse, const lduAddressing &fineAddressing, const labelUList &restrict, const label nCoarse)
Given restriction determines if coarse cells are connected.
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
Definition: polyMesh.C:818
UList< label > labelUList
Definition: UList.H:64
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
List< label > labelList
A List of labels.
Definition: labelList.H:56
volScalarField scalarField(fieldObject, mesh)
A class for managing temporary objects.
Definition: PtrList.H:53