pointConstraints.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) 2013-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 \*---------------------------------------------------------------------------*/
25 
26 #include "pointConstraints.H"
27 #include "emptyPointPatch.H"
28 #include "polyMesh.H"
29 #include "pointMesh.H"
30 #include "globalMeshData.H"
31 #include "twoDPointCorrector.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37  defineTypeNameAndDebug(pointConstraints, 0);
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 void Foam::pointConstraints::makePatchPatchAddressing()
44 {
45  if (debug)
46  {
47  Pout<< "pointConstraints::makePatchPatchAddressing() : "
48  << "constructing boundary addressing"
49  << endl << incrIndent;
50  }
51 
52  const pointMesh& pMesh = mesh();
53  const polyMesh& mesh = pMesh();
54 
55  const pointBoundaryMesh& pbm = pMesh.boundary();
56  const polyBoundaryMesh& bm = mesh.boundaryMesh();
57 
58 
59  // first count the total number of patch-patch points
60 
61  label nPatchPatchPoints = 0;
62 
63  forAll(pbm, patchi)
64  {
65  if (!isA<emptyPointPatch>(pbm[patchi]) && !pbm[patchi].coupled())
66  {
67  const labelList& bp = bm[patchi].boundaryPoints();
68 
69  nPatchPatchPoints += bp.size();
70 
71  if (debug)
72  {
73  Pout<< indent << "On patch:" << pbm[patchi].name()
74  << " nBoundaryPoints:" << bp.size() << endl;
75  }
76  }
77  }
78 
79  if (debug)
80  {
81  Pout<< indent << "Found nPatchPatchPoints:" << nPatchPatchPoints
82  << endl;
83  }
84 
85 
86  // Go through all patches and mark up the external edge points
87  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
88 
89  // From meshpoint to index in patchPatchPointConstraints_.
90  Map<label> patchPatchPointSet(2*nPatchPatchPoints);
91 
92  // Constraints (initialised to unconstrained)
93  patchPatchPointConstraints_.setSize(nPatchPatchPoints);
94  patchPatchPointConstraints_ = pointConstraint();
95 
96  // From constraint index to mesh point
97  labelList patchPatchPoints(nPatchPatchPoints);
98 
99  label pppi = 0;
100 
101  forAll(pbm, patchi)
102  {
103  if (!isA<emptyPointPatch>(pbm[patchi]) && !pbm[patchi].coupled())
104  {
105  const labelList& bp = bm[patchi].boundaryPoints();
106  const labelList& meshPoints = pbm[patchi].meshPoints();
107 
108  forAll(bp, pointi)
109  {
110  label ppp = meshPoints[bp[pointi]];
111 
112  Map<label>::iterator iter = patchPatchPointSet.find(ppp);
113 
114  label constraintI = -1;
115 
116  if (iter == patchPatchPointSet.end())
117  {
118  patchPatchPointSet.insert(ppp, pppi);
119  patchPatchPoints[pppi] = ppp;
120  constraintI = pppi++;
121  }
122  else
123  {
124  constraintI = iter();
125  }
126 
127  // Apply to patch constraints
128  pbm[patchi].applyConstraint
129  (
130  bp[pointi],
131  patchPatchPointConstraints_[constraintI]
132  );
133  }
134  }
135  }
136 
137  if (debug)
138  {
139  Pout<< indent << "Have (local) constrained points:"
140  << nPatchPatchPoints << endl;
141  }
142 
143 
144  // Extend set with constraints across coupled points
145  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
146 
147  {
148  const globalMeshData& gd = mesh.globalData();
149  const labelListList& globalPointSlaves = gd.globalPointSlaves();
150  const mapDistribute& globalPointSlavesMap = gd.globalPointSlavesMap();
151  const Map<label>& cpPointMap = gd.coupledPatch().meshPointMap();
152  const labelList& cpMeshPoints = gd.coupledPatch().meshPoints();
153 
154  // Constraints on coupled points
155  List<pointConstraint> constraints
156  (
157  globalPointSlavesMap.constructSize()
158  );
159 
160  // Copy from patchPatch constraints into coupledConstraints.
161  forAll(pbm, patchi)
162  {
163  if (!isA<emptyPointPatch>(pbm[patchi]) && !pbm[patchi].coupled())
164  {
165  const labelList& bp = bm[patchi].boundaryPoints();
166  const labelList& meshPoints = pbm[patchi].meshPoints();
167 
168  forAll(bp, pointi)
169  {
170  label ppp = meshPoints[bp[pointi]];
171 
172  Map<label>::const_iterator fnd = cpPointMap.find(ppp);
173  if (fnd != cpPointMap.end())
174  {
175  // Can just copy (instead of apply) constraint
176  // will already be consistent across multiple patches.
177  constraints[fnd()] = patchPatchPointConstraints_
178  [
179  patchPatchPointSet[ppp]
180  ];
181  }
182  }
183  }
184  }
185 
186  // Exchange data
187  globalPointSlavesMap.distribute(constraints);
188 
189  // Combine master with slave constraints
190  forAll(globalPointSlaves, pointi)
191  {
192  const labelList& slaves = globalPointSlaves[pointi];
193 
194  // Combine master constraint with slave constraints
195  forAll(slaves, i)
196  {
197  constraints[pointi].combine(constraints[slaves[i]]);
198  }
199  // Duplicate master constraint into slave slots
200  forAll(slaves, i)
201  {
202  constraints[slaves[i]] = constraints[pointi];
203  }
204  }
205 
206  // Send back
207  globalPointSlavesMap.reverseDistribute
208  (
209  cpMeshPoints.size(),
210  constraints
211  );
212 
213  // Add back into patchPatch constraints
214  forAll(constraints, coupledPointi)
215  {
216  if (constraints[coupledPointi].first() != 0)
217  {
218  label meshPointi = cpMeshPoints[coupledPointi];
219 
220  Map<label>::iterator iter = patchPatchPointSet.find(meshPointi);
221 
222  label constraintI = -1;
223 
224  if (iter == patchPatchPointSet.end())
225  {
226  // Pout<< indent << "on meshpoint:" << meshPointi
227  // << " coupled:" << coupledPointi
228  // << " at:" << mesh.points()[meshPointi]
229  // << " have new constraint:"
230  // << constraints[coupledPointi]
231  // << endl;
232 
233  // Allocate new constraint
234  if (patchPatchPoints.size() <= pppi)
235  {
236  patchPatchPoints.setSize(pppi+100);
237  }
238  patchPatchPointSet.insert(meshPointi, pppi);
239  patchPatchPoints[pppi] = meshPointi;
240  constraintI = pppi++;
241  }
242  else
243  {
244  // Pout<< indent << "on meshpoint:" << meshPointi
245  // << " coupled:" << coupledPointi
246  // << " at:" << mesh.points()[meshPointi]
247  // << " have possibly extended constraint:"
248  // << constraints[coupledPointi]
249  // << endl;
250 
251  constraintI = iter();
252  }
253 
254  // Combine (new or existing) constraint with one
255  // on coupled.
256  patchPatchPointConstraints_[constraintI].combine
257  (
258  constraints[coupledPointi]
259  );
260  }
261  }
262  }
263 
264 
265 
266  nPatchPatchPoints = pppi;
267  patchPatchPoints.setSize(nPatchPatchPoints);
268  patchPatchPointConstraints_.setSize(nPatchPatchPoints);
269 
270 
271  if (debug)
272  {
273  Pout<< indent << "Have (global) constrained points:"
274  << nPatchPatchPoints << endl;
275  }
276 
277 
278  // Copy out all non-trivial constraints
279  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
280 
281  patchPatchPointConstraintPoints_.setSize(nPatchPatchPoints);
282  patchPatchPointConstraintTensors_.setSize(nPatchPatchPoints);
283 
284  label nConstraints = 0;
285 
286  forAll(patchPatchPointConstraints_, i)
287  {
288  // Note: check for more than zero constraints. (could check for
289  // more than one constraint but what about coupled points that
290  // inherit the constraintness)
291  if (patchPatchPointConstraints_[i].first() != 0)
292  {
293  patchPatchPointConstraintPoints_[nConstraints] =
294  patchPatchPoints[i];
295 
296  patchPatchPointConstraintTensors_[nConstraints] =
297  patchPatchPointConstraints_[i].constraintTransformation();
298 
299  nConstraints++;
300  }
301  }
302 
303  if (debug)
304  {
305  Pout<< indent << "Have non-trivial constrained points:"
306  << nConstraints << endl;
307  }
308 
309  patchPatchPointConstraints_.setSize(nConstraints);
310  patchPatchPointConstraintPoints_.setSize(nConstraints);
311  patchPatchPointConstraintTensors_.setSize(nConstraints);
312 
313 
314  if (debug)
315  {
316  Pout<< decrIndent
317  << "pointConstraints::makePatchPatchAddressing() : "
318  << "finished constructing boundary addressing"
319  << endl;
320  }
321 }
322 
323 
324 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
325 
327 :
329 {
330  if (debug)
331  {
332  Pout<< "pointConstraints::pointConstraints(const pointMesh&): "
333  << "Constructing from pointMesh " << pm.name()
334  << endl;
335  }
336 
337  makePatchPatchAddressing();
338 }
339 
340 
341 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * //
342 
344 {
345  if (debug)
346  {
347  Pout<< "pointConstraints::~pointConstraints()" << endl;
348  }
349 }
350 
351 
352 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
353 
355 {
356  makePatchPatchAddressing();
357 }
358 
359 
361 {
362  return true;
363 }
364 
365 
367 (
368  pointVectorField& pf,
369  const bool overrideFixedValue
370 ) const
371 {
372  // Override constrained pointPatchField types with the constraint value.
373  // This relies on only constrained pointPatchField implementing the evaluate
374  // function
376 
377  // Sync any dangling points
379  (
380  pf.mesh()(),
381  pf.primitiveFieldRef(),
383  );
384 
385  // Apply multiple constraints on edge/corner points
386  constrainCorners(pf);
387 
388  // Apply any 2D motion constraints (or should they go before
389  // corner constraints?)
391  (
392  mesh()().points(),
393  pf.primitiveFieldRef()
394  );
395 
396  if (overrideFixedValue)
397  {
398  setPatchFields(pf);
399  }
400 }
401 
402 
403 template<>
404 void Foam::pointConstraints::constrainCorners<Foam::scalar>
405 (
406  GeometricField<scalar, pointPatchField, pointMesh>& pf
407 ) const
408 {}
409 
410 
411 template<>
412 void Foam::pointConstraints::constrainCorners<Foam::label>
413 (
414  GeometricField<label, pointPatchField, pointMesh>& pf
415 ) const
416 {}
417 
418 
419 // ************************************************************************* //
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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.
Definition: IOobject.H:295
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:226
HashTable< label, label, Hash< label > >::iterator iterator
Definition: Map.H:56
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
void constrainDisplacement(pointVectorField &displacement, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints),.
bool movePoints()
Correct weighting factors for moving mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
HashTable< label, label, Hash< label > >::const_iterator const_iterator
Definition: Map.H:59
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
pointConstraints(const pointMesh &)
Constructor from pointMesh.
static const twoDPointCorrector & New(const polyMesh &mesh)
Definition: MeshObject.C:44
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
Application of (multi-)patch point constraints.
Templated abstract base-class for optional mesh objects used to automate their allocation to the mesh...
Definition: MeshObject.H:86
void correctDisplacement(const pointField &p, vectorField &disp) const
Correct motion displacements.
dynamicFvMesh & mesh
const pointField & points
static void syncUntransformedData(const polyMesh &mesh, List< Type > &pointData, const CombineOp &cop)
Helper: sync data on collocated points only.
void constrainCorners(GeometricField< Type, pointPatchField, pointMesh > &pf) const
Apply patch-patch constraints only.
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.H:297
static void setPatchFields(GeometricField< Type, pointPatchField, pointMesh > &)
Helper: set patchField values from internal values (on.
List< label > labelList
A List of labels.
Definition: labelList.H:56
Internal::FieldType & primitiveFieldRef()
Return a reference to the internal field.
const Mesh & mesh() const
Return mesh.
defineTypeNameAndDebug(combustionModel, 0)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:240
void updateMesh(const mapPolyMesh &)
Update mesh topology using the morph engine.
label patchi
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
void correctBoundaryConditions()
Correct boundary field.
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:233
Namespace for OpenFOAM.
~pointConstraints()
Destructor.