faceAgglomerate.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-2018 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 Application
25  faceAgglomerate
26 
27 Description
28 
29  Agglomerate boundary faces using the pairPatchAgglomeration algorithm.
30  It writes a map from the fine to coarse grid.
31 
32 \*---------------------------------------------------------------------------*/
33 
34 #include "argList.H"
35 #include "fvMesh.H"
36 #include "Time.H"
37 #include "volFields.H"
38 #include "CompactListList.H"
39 #include "unitConversion.H"
40 #include "pairPatchAgglomeration.H"
41 #include "labelListIOList.H"
42 #include "syncTools.H"
43 #include "globalIndex.H"
44 
45 using namespace Foam;
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 int main(int argc, char *argv[])
50 {
51  #include "addRegionOption.H"
52  #include "addDictOption.H"
53  #include "setRootCase.H"
54  #include "createTime.H"
55  #include "createNamedMesh.H"
56 
57  const word dictName("viewFactorsDict");
58 
60 
61  // Read control dictionary
62  const IOdictionary agglomDict(dictIO);
63 
64  bool writeAgglom = readBool(agglomDict.lookup("writeFacesAgglomeration"));
65 
66 
67  const polyBoundaryMesh& boundary = mesh.boundaryMesh();
68 
69  labelListIOList finalAgglom
70  (
71  IOobject
72  (
73  "finalAgglom",
75  mesh,
78  false
79  ),
80  boundary.size()
81  );
82 
83  label nCoarseFaces = 0;
84 
85  forAllConstIter(dictionary, agglomDict, iter)
86  {
87  labelList patchids = boundary.findIndices(iter().keyword());
88  forAll(patchids, i)
89  {
90  label patchi = patchids[i];
91  const polyPatch& pp = boundary[patchi];
92  if (!pp.coupled())
93  {
94  Info << "\nAgglomerating patch : " << pp.name() << endl;
95  pairPatchAgglomeration agglomObject
96  (
97  pp,
98  agglomDict.subDict(pp.name())
99  );
100  agglomObject.agglomerate();
101  finalAgglom[patchi] =
102  agglomObject.restrictTopBottomAddressing();
103 
104  if (finalAgglom[patchi].size())
105  {
106  nCoarseFaces += max(finalAgglom[patchi] + 1);
107  }
108  }
109  }
110  }
111 
112 
113  // - All patches which are not agglomarated are identity for finalAgglom
114  forAll(boundary, patchid)
115  {
116  if (finalAgglom[patchid].size() == 0)
117  {
118  finalAgglom[patchid] = identity(boundary[patchid].size());
119  }
120  }
121 
122  // Sync agglomeration across coupled patches
123  labelList nbrAgglom(mesh.nFaces() - mesh.nInternalFaces(), -1);
124 
125  forAll(boundary, patchid)
126  {
127  const polyPatch& pp = boundary[patchid];
128  if (pp.coupled())
129  {
130  finalAgglom[patchid] = identity(pp.size());
131  forAll(pp, i)
132  {
133  nbrAgglom[pp.start() - mesh.nInternalFaces() + i] =
134  finalAgglom[patchid][i];
135  }
136  }
137  }
138 
140  forAll(boundary, patchid)
141  {
142  const polyPatch& pp = boundary[patchid];
143  if (pp.coupled() && !refCast<const coupledPolyPatch>(pp).owner())
144  {
145  forAll(pp, i)
146  {
147  finalAgglom[patchid][i] =
148  nbrAgglom[pp.start() - mesh.nInternalFaces() + i];
149  }
150  }
151  }
152 
153  finalAgglom.write();
154 
155  if (writeAgglom)
156  {
157  globalIndex index(nCoarseFaces);
158  volScalarField facesAgglomeration
159  (
160  IOobject
161  (
162  "facesAgglomeration",
163  mesh.time().timeName(),
164  mesh,
167  ),
168  mesh,
169  dimensionedScalar("facesAgglomeration", dimless, 0)
170  );
171 
172  volScalarField::Boundary& facesAgglomerationBf =
173  facesAgglomeration.boundaryFieldRef();
174 
175  label coarsePatchIndex = 0;
176  forAll(boundary, patchid)
177  {
178  const polyPatch& pp = boundary[patchid];
179  if (pp.size() > 0)
180  {
181  fvPatchScalarField& bFacesAgglomeration =
182  facesAgglomerationBf[patchid];
183 
184  forAll(bFacesAgglomeration, j)
185  {
186  bFacesAgglomeration[j] =
187  index.toGlobal
188  (
190  finalAgglom[patchid][j] + coarsePatchIndex
191  );
192  }
193 
194  coarsePatchIndex += max(finalAgglom[patchid]) + 1;
195  }
196  }
197 
198  Info<< "\nWriting facesAgglomeration" << endl;
199  facesAgglomeration.write();
200  }
201 
202  Info<< "End\n" << endl;
203  return 0;
204 }
205 
206 
207 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:424
#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
const word & name() const
Return name.
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:801
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
label nInternalFaces() const
label nFaces() const
Unit conversion functions.
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:427
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
labelList identity(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
bool readBool(Istream &)
Definition: boolIO.C:60
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:52
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
dynamicFvMesh & mesh
virtual bool coupled() const
Return true if this patch is geometrically coupled (i.e. faces and.
Definition: polyPatch.H:310
A class for handling words, derived from string.
Definition: word.H:59
const word dictName("particleTrackDict")
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Foam::polyBoundaryMesh.
IOobject dictIO(dictName, runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE)
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Primitive patch pair agglomerate method.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
messageStream Info
virtual Ostream & write(const token &)=0
Write next token to stream.
static void swapBoundaryFaceList(const polyMesh &mesh, UList< T > &l)
Swap coupled boundary face values.
Definition: syncTools.H:430
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Namespace for OpenFOAM.