FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 <finiteVolume/fvMesh.H>
31 #include <OpenFOAM/syncTools.H>
32 
33 //extern "C"
34 //{
35 //# include <mgridgen.h>
36 //}
37 
38 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
39 
40 void Foam::MGridGenGAMGAgglomeration::
41 makeCompactCellFaceAddressingAndFaceWeights
42 (
43  const lduAddressing& fineAddressing,
44  List<idxtype>& cellCells,
45  List<idxtype>& cellCellOffsets,
46  const vectorField& Si,
47  List<scalar>& faceWeights
48 )
49 {
50  const label nFineCells = fineAddressing.size();
51  const label nFineFaces = fineAddressing.upperAddr().size();
52 
53  const unallocLabelList& upperAddr = fineAddressing.upperAddr();
54  const unallocLabelList& lowerAddr = fineAddressing.lowerAddr();
55 
56  // Number of neighbours for each cell
57  labelList nNbrs(nFineCells, 0);
58 
59  forAll (upperAddr, facei)
60  {
61  nNbrs[upperAddr[facei]]++;
62  }
63 
64  forAll (lowerAddr, facei)
65  {
66  nNbrs[lowerAddr[facei]]++;
67  }
68 
69  // Set the sizes of the addressing and faceWeights arrays
70  cellCellOffsets.setSize(nFineCells + 1);
71  cellCells.setSize(2*nFineFaces);
72  faceWeights.setSize(2*nFineFaces);
73 
74 
75  cellCellOffsets[0] = 0;
76  forAll (nNbrs, celli)
77  {
78  cellCellOffsets[celli+1] = cellCellOffsets[celli] + nNbrs[celli];
79  }
80 
81  // reset the whole list to use as counter
82  nNbrs = 0;
83 
84  forAll (upperAddr, facei)
85  {
86  label own = upperAddr[facei];
87  label nei = lowerAddr[facei];
88 
89  label l1 = cellCellOffsets[own] + nNbrs[own]++;
90  label l2 = cellCellOffsets[nei] + nNbrs[nei]++;
91 
92  cellCells[l1] = nei;
93  cellCells[l2] = own;
94 
95  faceWeights[l1] = mag(Si[facei]);
96  faceWeights[l2] = mag(Si[facei]);
97  }
98 }
99 
100 
101 Foam::tmp<Foam::labelField> Foam::MGridGenGAMGAgglomeration::agglomerate
102 (
103  label& nCoarseCells,
104  const label minSize,
105  const label maxSize,
106  const lduAddressing& fineAddressing,
107  const scalarField& V,
108  const vectorField& Sf,
109  const scalarField& Sb
110 )
111 {
112  const label nFineCells = fineAddressing.size();
113 
114  // Compact addressing for cellCells
115  List<idxtype> cellCells;
116  List<idxtype> cellCellOffsets;
117 
118  // Face weights = face areas of the internal faces
119  List<scalar> faceWeights;
120 
121  // Create the compact addressing for cellCells and faceWeights
122  makeCompactCellFaceAddressingAndFaceWeights
123  (
124  fineAddressing,
125  cellCells,
126  cellCellOffsets,
127  Sf,
128  faceWeights
129  );
130 
131  // agglomeration options.
132  List<int> options(4, 0);
133  options[0] = 4; // globular agglom
134  options[1] = 6; // objective F3 and F2
135  options[2] = 128; // debugging output level
136  options[3] = fvMesh_.nGeometricD(); // Dimensionality of the grid
137 
138 
139  // output: cell -> processor addressing
140  List<int> finalAgglom(nFineCells);
141  int nMoves = -1;
142 
143  MGridGen
144  (
145  nFineCells,
146  cellCellOffsets.begin(),
147  const_cast<scalar*>(V.begin()),
148  const_cast<scalar*>(Sb.begin()),
149  cellCells.begin(),
150  faceWeights.begin(),
151  minSize,
152  maxSize,
153  options.begin(),
154  &nMoves,
155  &nCoarseCells,
156  finalAgglom.begin()
157  );
158 
159  return tmp<labelField>(new labelField(finalAgglom));
160 }
161 
162 
163 // ************************ vim: set sw=4 sts=4 et: ************************ //