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-2023 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 {
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, pointConstraint());
94 
95  // From constraint index to mesh point
96  labelList patchPatchPoints(nPatchPatchPoints);
97 
98  label pppi = 0;
99 
100  forAll(pbm, patchi)
101  {
102  if (!isA<emptyPointPatch>(pbm[patchi]) && !pbm[patchi].coupled())
103  {
104  const labelList& bp = bm[patchi].boundaryPoints();
105  const labelList& meshPoints = pbm[patchi].meshPoints();
106 
107  forAll(bp, pointi)
108  {
109  label ppp = meshPoints[bp[pointi]];
110 
111  Map<label>::iterator iter = patchPatchPointSet.find(ppp);
112 
113  label constraintI = -1;
114 
115  if (iter == patchPatchPointSet.end())
116  {
117  patchPatchPointSet.insert(ppp, pppi);
118  patchPatchPoints[pppi] = ppp;
119  constraintI = pppi++;
120  }
121  else
122  {
123  constraintI = iter();
124  }
125 
126  // Apply to patch constraints
127  pbm[patchi].applyConstraint
128  (
129  bp[pointi],
130  patchPatchPointConstraints_[constraintI]
131  );
132  }
133  }
134  }
135 
136  if (debug)
137  {
138  Pout<< indent << "Have (local) constrained points:"
139  << nPatchPatchPoints << endl;
140  }
141 
142 
143  // Extend set with constraints across coupled points
144  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
145 
146  {
147  const globalMeshData& gd = mesh.globalData();
148  const labelListList& globalPointSlaves = gd.globalPointSlaves();
149  const distributionMap& globalPointSlavesMap = gd.globalPointSlavesMap();
150  const Map<label>& cpPointMap = gd.coupledPatch().meshPointMap();
151  const labelList& cpMeshPoints = gd.coupledPatch().meshPoints();
152 
153  // Constraints on coupled points
154  List<pointConstraint> constraints
155  (
156  globalPointSlavesMap.constructSize()
157  );
158 
159  // Copy from patchPatch constraints into coupledConstraints.
160  forAll(pbm, patchi)
161  {
162  if (!isA<emptyPointPatch>(pbm[patchi]) && !pbm[patchi].coupled())
163  {
164  const labelList& bp = bm[patchi].boundaryPoints();
165  const labelList& meshPoints = pbm[patchi].meshPoints();
166 
167  forAll(bp, pointi)
168  {
169  label ppp = meshPoints[bp[pointi]];
170 
171  Map<label>::const_iterator fnd = cpPointMap.find(ppp);
172  if (fnd != cpPointMap.end())
173  {
174  // Can just copy (instead of apply) constraint
175  // will already be consistent across multiple patches.
176  constraints[fnd()] = patchPatchPointConstraints_
177  [
178  patchPatchPointSet[ppp]
179  ];
180  }
181  }
182  }
183  }
184 
185  // Exchange data
186  globalPointSlavesMap.distribute(constraints);
187 
188  // Combine master with slave constraints
189  forAll(globalPointSlaves, pointi)
190  {
191  const labelList& slaves = globalPointSlaves[pointi];
192 
193  // Combine master constraint with slave constraints
194  forAll(slaves, i)
195  {
196  constraints[pointi].combine(constraints[slaves[i]]);
197  }
198  // Duplicate master constraint into slave slots
199  forAll(slaves, i)
200  {
201  constraints[slaves[i]] = constraints[pointi];
202  }
203  }
204 
205  // Send back
206  globalPointSlavesMap.reverseDistribute
207  (
208  cpMeshPoints.size(),
209  constraints
210  );
211 
212  // Add back into patchPatch constraints
213  forAll(constraints, coupledPointi)
214  {
215  if (constraints[coupledPointi].first() != 0)
216  {
217  label meshPointi = cpMeshPoints[coupledPointi];
218 
219  Map<label>::iterator iter = patchPatchPointSet.find(meshPointi);
220 
221  label constraintI = -1;
222 
223  if (iter == patchPatchPointSet.end())
224  {
225  // Pout<< indent << "on meshpoint:" << meshPointi
226  // << " coupled:" << coupledPointi
227  // << " at:" << mesh.points()[meshPointi]
228  // << " have new constraint:"
229  // << constraints[coupledPointi]
230  // << endl;
231 
232  // Allocate new constraint
233  if (patchPatchPoints.size() <= pppi)
234  {
235  patchPatchPoints.setSize(pppi+100);
236  }
237  patchPatchPointSet.insert(meshPointi, pppi);
238  patchPatchPoints[pppi] = meshPointi;
239  constraintI = pppi++;
240  }
241  else
242  {
243  // Pout<< indent << "on meshpoint:" << meshPointi
244  // << " coupled:" << coupledPointi
245  // << " at:" << mesh.points()[meshPointi]
246  // << " have possibly extended constraint:"
247  // << constraints[coupledPointi]
248  // << endl;
249 
250  constraintI = iter();
251  }
252 
253  // Extend the patchPatchPointConstraints_ array if necessary
254  if (patchPatchPointConstraints_.size() <= constraintI)
255  {
256  patchPatchPointConstraints_.setSize
257  (
258  constraintI + 1,
259  pointConstraint()
260  );
261  }
262 
263  // Combine (new or existing) constraint with one
264  // on coupled.
265  patchPatchPointConstraints_[constraintI].combine
266  (
267  constraints[coupledPointi]
268  );
269  }
270  }
271  }
272 
273 
274  nPatchPatchPoints = pppi;
275  patchPatchPoints.setSize(nPatchPatchPoints);
276  patchPatchPointConstraints_.setSize(nPatchPatchPoints);
277 
278 
279  if (debug)
280  {
281  Pout<< indent << "Have (global) constrained points:"
282  << nPatchPatchPoints << endl;
283  }
284 
285 
286  // Copy out all non-trivial constraints
287  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
288 
289  patchPatchPointConstraintPoints_.setSize(nPatchPatchPoints);
290  patchPatchPointConstraintTensors_.setSize(nPatchPatchPoints);
291 
292  label nConstraints = 0;
293 
294  forAll(patchPatchPointConstraints_, i)
295  {
296  // Note: check for more than zero constraints. (could check for
297  // more than one constraint but what about coupled points that
298  // inherit the constraintness)
299  if (patchPatchPointConstraints_[i].first() != 0)
300  {
301  patchPatchPointConstraintPoints_[nConstraints] =
302  patchPatchPoints[i];
303 
304  patchPatchPointConstraintTensors_[nConstraints] =
305  patchPatchPointConstraints_[i].constraintTransformation();
306 
307  nConstraints++;
308  }
309  }
310 
311  if (debug)
312  {
313  Pout<< indent << "Have non-trivial constrained points:"
314  << nConstraints << endl;
315  }
316 
317  patchPatchPointConstraints_.setSize(nConstraints);
318  patchPatchPointConstraintPoints_.setSize(nConstraints);
319  patchPatchPointConstraintTensors_.setSize(nConstraints);
320 
321 
322  if (debug)
323  {
324  Pout<< decrIndent
325  << "pointConstraints::makePatchPatchAddressing() : "
326  << "finished constructing boundary addressing"
327  << endl;
328  }
329 }
330 
331 
332 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
333 
335 :
337  <
338  pointMesh,
341  >(pm)
342 {
343  if (debug)
344  {
345  Pout<< "pointConstraints::pointConstraints(const pointMesh&): "
346  << "Constructing from pointMesh " << pm.name()
347  << endl;
348  }
349 
350  makePatchPatchAddressing();
351 }
352 
353 
354 // * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * * //
355 
357 {
358  if (debug)
359  {
360  Pout<< "pointConstraints::~pointConstraints()" << endl;
361  }
362 }
363 
364 
365 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
366 
368 {
369  return true;
370 }
371 
372 
374 {
375  makePatchPatchAddressing();
376 }
377 
378 
380 {
381  makePatchPatchAddressing();
382 }
383 
384 
386 {
388  makePatchPatchAddressing();
389 }
390 
391 
393 (
394  pointVectorField& pf,
395  const bool overrideFixedValue
396 ) const
397 {
398  // Override constrained pointPatchField types with the constraint value.
399  // This relies on only constrained pointPatchField implementing the evaluate
400  // function
402 
403  // Sync any dangling points
404  syncUntransformedData
405  (
406  pf.mesh()(),
407  pf.primitiveFieldRef(),
409  );
410 
411  // Apply multiple constraints on edge/corner points
412  constrainCorners(pf);
413 
414  // Apply any 2D motion constraints (or should they go before
415  // corner constraints?)
417  (
418  mesh()().points(),
419  pf.primitiveFieldRef()
420  );
421 
422  if (overrideFixedValue)
423  {
424  setPatchFields(pf);
425  }
426 }
427 
428 
429 template<>
430 void Foam::pointConstraints::constrainCorners<Foam::scalar>
431 (
432  PointField<scalar>& pf
433 ) const
434 {}
435 
436 
437 template<>
438 void Foam::pointConstraints::constrainCorners<Foam::label>
439 (
440  PointField<label>& pf
441 ) const
442 {}
443 
444 
445 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Templated abstract base-class for demand-driven mesh objects used to automate their allocation to the...
static twoDPointCorrector & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
const Mesh & mesh() const
Return mesh.
Internal::FieldType & primitiveFieldRef()
Return a reference to the primitive field.
void correctBoundaryConditions()
Correct boundary field.
friend class iterator
Declare friendship with the iterator.
Definition: HashTable.H:194
friend class const_iterator
Declare friendship with the const_iterator.
Definition: HashTable.H:197
const word & name() const
Return name.
Definition: IOobject.H:310
virtual const fileName & name() const
Return the name of the stream.
Definition: IOstream.H:294
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
void setSize(const label)
Reset size of List.
Definition: List.C:281
Application of (multi-)patch point constraints.
virtual bool movePoints()
Correct weighting factors for moving mesh.
pointConstraints(const pointMesh &)
Constructor from pointMesh.
virtual void topoChange(const polyTopoChangeMap &)
Update mesh topology using the morph engine.
virtual void distribute(const polyDistributionMap &)
Update corresponding to the given distribution map.
virtual void mapMesh(const polyMeshMap &)
Update from another mesh using the given map.
~pointConstraints()
Destructor.
void constrainDisplacement(pointVectorField &displacement, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints),.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:53
const globalMeshData & globalData() const
Return parallel info.
Definition: pointMesh.H:131
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information.
Definition: polyMeshMap.H:51
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
void correctDisplacement(const pointField &p, vectorField &disp) const
Correct motion displacements.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
const pointField & points
Namespace for OpenFOAM.
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
Definition: Ostream.H:241
List< label > labelList
A List of labels.
Definition: labelList.H:56
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
errorManip< error > abort(error &err)
Definition: errorManip.H:131
Ostream & incrIndent(Ostream &os)
Increment the indent level.
Definition: Ostream.H:234
labelList first(const UList< labelPair > &p)
Definition: patchToPatch.C:39
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:57
defineTypeNameAndDebug(combustionModel, 0)
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
error FatalError
Ostream & indent(Ostream &os)
Indent stream.
Definition: Ostream.H:227