assemblyFaceAreaPairGAMGAgglomeration.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2018-2019 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
29 #include "lduMatrix.H"
32 #include "surfaceFields.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38  defineTypeNameAndDebug(assemblyFaceAreaPairGAMGAgglomeration, 0);
39 
41  (
42  GAMGAgglomeration,
43  assemblyFaceAreaPairGAMGAgglomeration,
44  lduMatrix
45  );
46 
48  (
49  GAMGAgglomeration,
50  assemblyFaceAreaPairGAMGAgglomeration,
51  geometry
52  );
53 }
54 
55 
56 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * * //
57 
60 {}
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
66 (
67  const lduMatrix& matrix,
68  const dictionary& controlDict
69 )
70 :
72 {
73  const lduMesh& ldumesh = matrix.mesh();
74 
75  if (isA<lduPrimitiveMeshAssembly>(ldumesh))
76  {
78  refCast<const lduPrimitiveMeshAssembly>(ldumesh);
79 
80  vectorField faceAreas(mesh.lduAddr().upperAddr().size(), Zero);
81 
82  const labelListList& faceMap = mesh.faceMap();
83 
84  for (label i=0; i < mesh.meshes().size(); ++i)
85  {
86  const fvMesh& m = refCast<const fvMesh>(mesh.meshes()[i]);
87  const labelList& subFaceMap = faceMap[i];
88  const vectorField& areas = m.Sf().internalField();
89 
90  forAll(subFaceMap, facei)
91  {
92  faceAreas[subFaceMap[facei]] = areas[facei];
93  }
94 
95  const polyBoundaryMesh& patches = m.boundaryMesh();
96 
97  // Fill faceAreas for new faces
98  forAll(patches, patchI)
99  {
100  const polyPatch& pp = patches[patchI];
101  label globalPatchID = mesh.patchMap()[i][patchI];
102 
103  if (globalPatchID == -1)
104  {
105  if (pp.masterImplicit())
106  {
107  const vectorField& sf = m.boundary()[patchI].Sf();
108 
109  if (isA<cyclicAMIPolyPatch>(pp))
110  {
111  const cyclicAMIPolyPatch& mpp =
112  refCast<const cyclicAMIPolyPatch>(pp);
113 
114  const scalarListList& srcWeight =
115  mpp.AMI().srcWeights();
116 
117  label subFaceI = 0;
118  forAll(pp.faceCells(), faceI)
119  {
120  const scalarList& w = srcWeight[faceI];
121 
122  for(label j=0; j<w.size(); j++)
123  {
124  const label globalFaceI =
125  mesh.faceBoundMap()[i][patchI][subFaceI];
126 
127  if (globalFaceI != -1)
128  {
129  faceAreas[globalFaceI] = w[j]*sf[faceI];
130  }
131  subFaceI++;
132  }
133  }
134  }
135  else if (isA<cyclicACMIPolyPatch>(pp))
136  {
137  const cyclicACMIPolyPatch& mpp =
138  refCast<const cyclicACMIPolyPatch>(pp);
139 
140  const scalarListList& srcWeight =
141  mpp.AMI().srcWeights();
142 
143  const scalarField& mask = mpp.mask();
144  const scalar tol = mpp.tolerance();
145  label subFaceI = 0;
146  forAll(mask, faceI)
147  {
148  const scalarList& w = srcWeight[faceI];
149 
150  for(label j=0; j<w.size(); j++)
151  {
152  if (mask[faceI] > tol)
153  {
154  const label globalFaceI =
155  mesh.faceBoundMap()[i]
156  [patchI][subFaceI];
157 
158  faceAreas[globalFaceI] = w[j]*sf[faceI];
159  }
160  subFaceI++;
161  }
162  }
163  }
164  else
165  {
166  forAll(pp.faceCells(), faceI)
167  {
168  const label globalFaceI =
169  mesh.faceBoundMap()[i][patchI][faceI];
170 
171  if (globalFaceI != -1)
172  {
173  faceAreas[globalFaceI] = sf[faceI];
174  }
175  }
176  }
177  }
178  }
179  }
180  }
181 
183  (
184  mesh,
185  mag
186  (
188  (
189  faceAreas/sqrt(mag(faceAreas)),
190  vector(1, 1.01, 1.02)
191  )
192  )
193  );
194  }
195  else
196  {
197  // Revert to faceAreaPairGAMGAgglomeration
198  const fvMesh& fvmesh = refCast<const fvMesh>(matrix.mesh());
199 
201  (
202  matrix.mesh(),
203  mag
204  (
206  (
207  fvmesh.Sf().primitiveField()
208  /sqrt(fvmesh.magSf().primitiveField()),
209  vector(1, 1.01, 1.02)
210  //vector::one
211  )
212  )
213  );
214  }
215 }
216 
217 
220 (
221  const lduMatrix& matrix,
222  const scalarField& cellVolumes,
223  const vectorField& faceAreas,
224  const dictionary& controlDict
225 )
226 :
227  pairGAMGAgglomeration(matrix.mesh(), controlDict)
228 {
230  (
231  matrix.mesh(),
232  mag
233  (
235  (
236  faceAreas/sqrt(mag(faceAreas)),
237  vector(1, 1.01, 1.02)
238  )
239  )
240  );
241 }
242 
243 
244 // ************************************************************************* //
Foam::surfaceFields.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:584
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:118
const surfaceVectorField & Sf() const
Return cell face area vectors.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:120
Agglomerate using the pair algorithm.
const Internal & internalField() const noexcept
Return a const-reference to the dimensioned internal field.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:53
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:413
void agglomerate(const lduMesh &mesh, const scalarField &faceWeights)
Agglomerate all levels starting from the given face weights.
dynamicFvMesh & mesh
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
virtual const labelUList & upperAddr() const =0
Return upper addressing.
List< scalar > scalarList
A List of scalars.
Definition: scalarList.H:61
virtual const lduAddressing & lduAddr() const =0
Return ldu addressing.
Vector< scalar > vector
Definition: vector.H:57
runTime controlDict().readEntry("adjustTimeStep"
The central control dictionary, the contents of which are either taken directly from the FOAM_CONTROL...
Definition: debug.C:142
dimensioned< Type > cmptMultiply(const dimensioned< Type > &, const dimensioned< Type > &)
assemblyFaceAreaPairGAMGAgglomeration(const lduMatrix &matrix, const dictionary &controlDict)
Construct given mesh and controls.
defineTypeNameAndDebug(combustionModel, 0)
An assembly of lduMatrix that is specific inter-region coupling through mapped patches.
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
lduMatrix is a general matrix class in which the coefficients are stored as three arrays...
Definition: lduMatrix.H:79
List< scalarList > scalarListList
A List of scalarList.
Definition: scalarList.H:63
const lduMesh & mesh() const noexcept
Return the LDU mesh from which the addressing is obtained.
Definition: lduMatrix.H:698
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:79
const polyBoundaryMesh & patches
Field< vector > vectorField
Specialisation of Field<T> for vector.
Namespace for OpenFOAM.
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:705
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:157
const labelListList & faceMap(const label fineLeveli) const
Mapping from processor to procMesh face.