preserveBafflesConstraint.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) 2015-2016 OpenFOAM Foundation
9  Copyright (C) 2018 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
31 #include "syncTools.H"
32 #include "localPointRegion.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace decompositionConstraints
39 {
40  defineTypeName(preserveBaffles);
41 
43  (
44  decompositionConstraint,
45  preserveBaffles,
46  dictionary
47  );
48 }
49 }
50 
51 
52 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
53 
55 (
56  const dictionary& dict
57 )
58 :
59  decompositionConstraint(dict, typeName)
60 {
62  {
63  Info<< type()
64  << " : setting constraints to preserve baffles"
65  //<< returnReduce(bafflePairs.size(), sumOp<label>())
66  << endl;
67  }
68 }
69 
70 
72 :
74 {
76  {
77  Info<< type()
78  << " : setting constraints to preserve baffles"
79  //<< returnReduce(bafflePairs.size(), sumOp<label>())
80  << endl;
81  }
82 }
83 
84 
85 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
86 
88 (
89  const polyMesh& mesh,
90  boolList& blockedFace,
91  PtrList<labelList>& specifiedProcessorFaces,
92  labelList& specifiedProcessor,
93  List<labelPair>& explicitConnections
94 ) const
95 {
96  const labelPairList bafflePairs
97  (
99  );
100 
102  {
103  Info<< type() << " : setting constraints to preserve "
104  << returnReduce(bafflePairs.size(), sumOp<label>())
105  << " baffles" << endl;
106  }
107 
108 
109  // Merge into explicitConnections
110  {
111  // Convert into face-to-face addressing
112  labelList faceToFace(mesh.nFaces(), -1);
113  for (const labelPair& p : explicitConnections)
114  {
115  faceToFace[p[0]] = p[1];
116  faceToFace[p[1]] = p[0];
117  }
118 
119  // Merge in bafflePairs
120  for (const labelPair& p : bafflePairs)
121  {
122  if (faceToFace[p[0]] == -1 && faceToFace[p[1]] == -1)
123  {
124  faceToFace[p[0]] = p[1];
125  faceToFace[p[1]] = p[0];
126  }
127  else if (labelPair::compare(p, labelPair(p[0], faceToFace[p[0]])))
128  {
129  // Connection already present
130  }
131  else
132  {
133  const label p0Slave = faceToFace[p[0]];
134  const label p1Slave = faceToFace[p[1]];
135  IOWarningInFunction(coeffDict_)
136  << "When adding baffle between faces "
137  << p[0] << " at " << mesh.faceCentres()[p[0]]
138  << " and "
139  << p[1] << " at " << mesh.faceCentres()[p[1]]
140  << " : face " << p[0] << " already is connected to face "
141  << p0Slave << " at " << mesh.faceCentres()[p0Slave]
142  << " and face " << p[1] << " already is connected to face "
143  << p1Slave << " at " << mesh.faceCentres()[p1Slave]
144  << endl;
145  }
146  }
147 
148  // Convert back into labelPairList
149  label n = 0;
150  forAll(faceToFace, faceI)
151  {
152  label otherFaceI = faceToFace[faceI];
153  if (otherFaceI != -1 && faceI < otherFaceI)
154  {
155  // I am master of slave
156  n++;
157  }
158  }
159  explicitConnections.setSize(n);
160  n = 0;
161  forAll(faceToFace, faceI)
162  {
163  label otherFaceI = faceToFace[faceI];
164  if (otherFaceI != -1 && faceI < otherFaceI)
165  {
166  explicitConnections[n++] = labelPair(faceI, otherFaceI);
167  }
168  }
169  }
170 
171  // Make sure blockedFace is uptodate
172  blockedFace.setSize(mesh.nFaces(), true);
173  for (const labelPair& p : explicitConnections)
174  {
175  blockedFace[p.first()] = false;
176  blockedFace[p.second()] = false;
177  }
178  syncTools::syncFaceList(mesh, blockedFace, andEqOp<bool>());
179 }
180 
181 
183 (
184  const polyMesh& mesh,
185  const boolList& blockedFace,
186  const PtrList<labelList>& specifiedProcessorFaces,
187  const labelList& specifiedProcessor,
188  const List<labelPair>& explicitConnections,
189  labelList& decomposition
190 ) const
191 {
192  const labelPairList bafflePairs
193  (
195  );
196 
197  label nChanged = 0;
198 
199  for (const labelPair& baffle : bafflePairs)
200  {
201  const label f0 = baffle.first();
202  const label f1 = baffle.second();
203 
204  const label procI = decomposition[mesh.faceOwner()[f0]];
205 
206  if (mesh.isInternalFace(f0))
207  {
208  const label nei0 = mesh.faceNeighbour()[f0];
209  if (decomposition[nei0] != procI)
210  {
211  decomposition[nei0] = procI;
212  ++nChanged;
213  }
214  }
215 
216  const label own1 = mesh.faceOwner()[f1];
217  if (decomposition[own1] != procI)
218  {
219  decomposition[own1] = procI;
220  ++nChanged;
221  }
222  if (mesh.isInternalFace(f1))
223  {
224  const label nei1 = mesh.faceNeighbour()[f1];
225  if (decomposition[nei1] != procI)
226  {
227  decomposition[nei1] = procI;
228  }
229  }
230  }
231 
233  {
234  reduce(nChanged, sumOp<label>());
235  Info<< type() << " : changed decomposition on " << nChanged
236  << " cells" << endl;
237  }
238 }
239 
240 
241 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::decompositionConstraint
Abstract class for handling decomposition constraints.
Definition: decompositionConstraint.H:58
Foam::decompositionConstraints::defineTypeName
defineTypeName(geometric)
Foam::decompositionConstraints::preserveBaffles::apply
virtual void apply(const polyMesh &mesh, const boolList &blockedFace, const PtrList< labelList > &specifiedProcessorFaces, const labelList &specifiedProcessor, const List< labelPair > &explicitConnections, labelList &decomposition) const
Apply any additional post-decomposition constraints.
Definition: preserveBafflesConstraint.C:183
Foam::returnReduce
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
Definition: PstreamReduceOps.H:94
Foam::andEqOp
Definition: ops.H:85
Foam::decompositionConstraints::preserveBaffles::preserveBaffles
preserveBaffles()
Construct from components.
Definition: preserveBafflesConstraint.C:71
localPointRegion.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
syncTools.H
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::sumOp
Definition: ops.H:213
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::labelPair
Pair< label > labelPair
A pair of labels.
Definition: Pair.H:54
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::reduce
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
Definition: PstreamReduceOps.H:51
Foam::syncTools::syncFaceList
static void syncFaceList(const polyMesh &mesh, UList< T > &faceValues, const CombineOp &cop)
Synchronize values on all mesh faces.
Definition: syncTools.H:396
Foam::Pair< label >::compare
static int compare(const Pair< label > &a, const Pair< label > &b)
Compare Pairs.
Definition: PairI.H:33
Foam::faceToFace
A topoSetFaceSource to select faces based on usage in another faceSet.
Definition: faceToFace.H:159
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:59
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
preserveBafflesConstraint.H
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::decompositionConstraints::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionConstraint, geometric, dictionary)
Foam::Pair< label >
Foam::decompositionConstraints::preserveBaffles::add
virtual void add(const polyMesh &mesh, boolList &blockedFace, PtrList< labelList > &specifiedProcessorFaces, labelList &specifiedProcessor, List< labelPair > &explicitConnections) const
Add my constraints to list of constraints.
Definition: preserveBafflesConstraint.C:88
Foam::List< bool >
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::roots::type
type
Types of root.
Definition: Roots.H:54
IOWarningInFunction
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
Definition: messageStream.H:340
Foam::localPointRegion::findDuplicateFacePairs
static labelPairList findDuplicateFacePairs(const polyMesh &)
Helper routine to find all baffles (two boundary faces.
Definition: localPointRegion.C:625