slidingInterface.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-2022 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 "slidingInterface.H"
27 #include "polyTopoChanger.H"
28 #include "polyMesh.H"
29 #include "polyTopoChange.H"
31 #include "plane.H"
32 
33 // Index of debug signs:
34 // p - adjusting a projection point
35 // * - adjusting edge intersection
36 
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
38 
39 namespace Foam
40 {
41  defineTypeNameAndDebug(slidingInterface, 0);
43  (
44  polyMeshModifier,
45  slidingInterface,
46  dictionary
47  );
48 
49  template<>
50  const char* Foam::NamedEnum
51  <
53  2
54  >::names[] =
55  {
56  "integral",
57  "partial"
58  };
59 }
60 
61 
64 
65 
66 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
67 
68 void Foam::slidingInterface::checkDefinition()
69 {
70  const polyMesh& mesh = topoChanger().mesh();
71 
72  if
73  (
74  !masterFaceZoneID_.active()
75  || !slaveFaceZoneID_.active()
76  || !cutPointZoneID_.active()
77  || !cutFaceZoneID_.active()
78  || !masterPatchID_.active()
79  || !slavePatchID_.active()
80  )
81  {
83  << "Not all zones and patches needed in the definition "
84  << "have been found. Please check your mesh definition."
85  << abort(FatalError);
86  }
87 
88  // Check the sizes and set up state
89  if
90  (
91  mesh.faceZones()[masterFaceZoneID_.index()].empty()
92  || mesh.faceZones()[slaveFaceZoneID_.index()].empty()
93  )
94  {
96  << "Please check your mesh definition."
97  << abort(FatalError);
98  }
99 
100  if (debug)
101  {
102  Pout<< "Sliding interface object " << name() << " :" << nl
103  << " master face zone: " << masterFaceZoneID_.index() << nl
104  << " slave face zone: " << slaveFaceZoneID_.index() << endl;
105  }
106 }
107 
108 
109 void Foam::slidingInterface::clearOut() const
110 {
111  clearPointProjection();
112  clearAttachedAddressing();
113  clearAddressing();
114 }
115 
116 
117 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
118 
119 
120 // Construct from components
122 (
123  const word& name,
124  const label index,
125  const polyTopoChanger& mme,
126  const word& masterFaceZoneName,
127  const word& slaveFaceZoneName,
128  const word& cutPointZoneName,
129  const word& cutFaceZoneName,
130  const word& masterPatchName,
131  const word& slavePatchName,
132  const typeOfMatch tom,
133  const bool coupleDecouple,
134  const intersection::algorithm algo
135 )
136 :
137  polyMeshModifier(name, index, mme, true),
138  masterFaceZoneID_
139  (
140  masterFaceZoneName,
141  mme.mesh().faceZones()
142  ),
143  slaveFaceZoneID_
144  (
145  slaveFaceZoneName,
146  mme.mesh().faceZones()
147  ),
148  cutPointZoneID_
149  (
150  cutPointZoneName,
151  mme.mesh().pointZones()
152  ),
153  cutFaceZoneID_
154  (
155  cutFaceZoneName,
156  mme.mesh().faceZones()
157  ),
158  masterPatchID_
159  (
160  masterPatchName,
161  mme.mesh().boundaryMesh()
162  ),
163  slavePatchID_
164  (
165  slavePatchName,
166  mme.mesh().boundaryMesh()
167  ),
168  matchType_(tom),
169  coupleDecouple_(coupleDecouple),
170  attached_(false),
171  projectionAlgo_(algo),
172  trigger_(false),
173  pointMergeTol_(pointMergeTolDefault_),
174  edgeMergeTol_(edgeMergeTolDefault_),
175  nFacesPerSlaveEdge_(nFacesPerSlaveEdgeDefault_),
176  edgeFaceEscapeLimit_(edgeFaceEscapeLimitDefault_),
177  integralAdjTol_(integralAdjTolDefault_),
178  edgeMasterCatchFraction_(edgeMasterCatchFractionDefault_),
179  edgeCoPlanarTol_(edgeCoPlanarTolDefault_),
180  edgeEndCutoffTol_(edgeEndCutoffTolDefault_),
181  cutFaceMasterPtr_(nullptr),
182  cutFaceSlavePtr_(nullptr),
183  masterFaceCellsPtr_(nullptr),
184  slaveFaceCellsPtr_(nullptr),
185  masterStickOutFacesPtr_(nullptr),
186  slaveStickOutFacesPtr_(nullptr),
187  retiredPointMapPtr_(nullptr),
188  cutPointEdgePairMapPtr_(nullptr),
189  slavePointPointHitsPtr_(nullptr),
190  slavePointEdgeHitsPtr_(nullptr),
191  slavePointFaceHitsPtr_(nullptr),
192  masterPointEdgeHitsPtr_(nullptr),
193  projectedSlavePointsPtr_(nullptr)
194 {
195  checkDefinition();
196 
197  if (attached_)
198  {
200  << "Creation of a sliding interface from components "
201  << "in attached state not supported."
202  << abort(FatalError);
203  }
204  else
205  {
206  calcAttachedAddressing();
207  }
208 }
209 
210 
211 // Construct from components
213 (
214  const word& name,
215  const dictionary& dict,
216  const label index,
217  const polyTopoChanger& mme
218 )
219 :
220  polyMeshModifier(name, index, mme, Switch(dict.lookup("active"))),
221  masterFaceZoneID_
222  (
223  dict.lookup("masterFaceZoneName"),
224  mme.mesh().faceZones()
225  ),
226  slaveFaceZoneID_
227  (
228  dict.lookup("slaveFaceZoneName"),
229  mme.mesh().faceZones()
230  ),
231  cutPointZoneID_
232  (
233  dict.lookup("cutPointZoneName"),
234  mme.mesh().pointZones()
235  ),
236  cutFaceZoneID_
237  (
238  dict.lookup("cutFaceZoneName"),
239  mme.mesh().faceZones()
240  ),
241  masterPatchID_
242  (
243  dict.lookup("masterPatchName"),
244  mme.mesh().boundaryMesh()
245  ),
246  slavePatchID_
247  (
248  dict.lookup("slavePatchName"),
249  mme.mesh().boundaryMesh()
250  ),
251  matchType_(typeOfMatchNames_.read((dict.lookup("typeOfMatch")))),
252  coupleDecouple_(dict.lookup("coupleDecouple")),
253  attached_(dict.lookup("attached")),
254  projectionAlgo_
255  (
256  intersection::algorithmNames_.read(dict.lookup("projection"))
257  ),
258  trigger_(false),
259  cutFaceMasterPtr_(nullptr),
260  cutFaceSlavePtr_(nullptr),
261  masterFaceCellsPtr_(nullptr),
262  slaveFaceCellsPtr_(nullptr),
263  masterStickOutFacesPtr_(nullptr),
264  slaveStickOutFacesPtr_(nullptr),
265  retiredPointMapPtr_(nullptr),
266  cutPointEdgePairMapPtr_(nullptr),
267  slavePointPointHitsPtr_(nullptr),
268  slavePointEdgeHitsPtr_(nullptr),
269  slavePointFaceHitsPtr_(nullptr),
270  masterPointEdgeHitsPtr_(nullptr),
271  projectedSlavePointsPtr_(nullptr)
272 {
273  // Optionally default tolerances from dictionary
274  setTolerances(dict);
275 
276  checkDefinition();
277 
278  // If the interface is attached, the master and slave face zone addressing
279  // needs to be read in; otherwise it will be created
280  if (attached_)
281  {
282  if (debug)
283  {
284  Pout<< "slidingInterface::slidingInterface(...) "
285  << " for object " << name << " : "
286  << "Interface attached. Reading master and slave face zones "
287  << "and retired point lookup." << endl;
288  }
289 
290  // The face zone addressing is written out in the definition dictionary
291  masterFaceCellsPtr_ = new labelList(dict.lookup("masterFaceCells"));
292  slaveFaceCellsPtr_ = new labelList(dict.lookup("slaveFaceCells"));
293 
294  masterStickOutFacesPtr_ =
295  new labelList(dict.lookup("masterStickOutFaces"));
296  slaveStickOutFacesPtr_ =
297  new labelList(dict.lookup("slaveStickOutFaces"));
298 
299  retiredPointMapPtr_ = new Map<label>(dict.lookup("retiredPointMap"));
300  cutPointEdgePairMapPtr_ =
301  new Map<Pair<edge>>(dict.lookup("cutPointEdgePairMap"));
302  }
303  else
304  {
305  calcAttachedAddressing();
306  }
307 }
308 
309 
310 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
311 
313 {
314  clearOut();
315 }
316 
317 
318 void Foam::slidingInterface::clearAddressing() const
319 {
320  deleteDemandDrivenData(cutFaceMasterPtr_);
321  deleteDemandDrivenData(cutFaceSlavePtr_);
322 }
323 
324 
325 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
326 
328 {
329  return masterFaceZoneID_;
330 }
331 
332 
334 {
335  return slaveFaceZoneID_;
336 }
337 
338 
340 {
341  if (coupleDecouple_)
342  {
343  // Always changes. If not attached, project points
344  if (debug)
345  {
346  Pout<< "bool slidingInterface::changeTopology() const "
347  << "for object " << name() << " : "
348  << "Couple-decouple mode." << endl;
349  }
350 
351  if (!attached_)
352  {
353  projectPoints();
354  }
355  else
356  {
357  }
358 
359  return true;
360  }
361 
362  if
363  (
364  attached_
365  && !topoChanger().mesh().changing()
366  )
367  {
368  // If the mesh is not moving or morphing and the interface is
369  // already attached, the topology will not change
370  return false;
371  }
372  else
373  {
374  // Check if the motion changes point projection
375  return projectPoints();
376  }
377 }
378 
379 
381 {
382  if (coupleDecouple_)
383  {
384  if (attached_)
385  {
386  // Attached, coupling
387  decoupleInterface(ref);
388  }
389  else
390  {
391  // Detached, coupling
392  coupleInterface(ref);
393  }
394 
395  return;
396  }
397 
398  if (trigger_)
399  {
400  if (attached_)
401  {
402  // Clear old coupling data
403  clearCouple(ref);
404  }
405 
406  coupleInterface(ref);
407 
408  trigger_ = false;
409  }
410 }
411 
412 
414 {
415  if (debug)
416  {
417  Pout<< "void slidingInterface::modifyMotionPoints("
418  << "pointField& motionPoints) const for object " << name() << " : "
419  << "Adjusting motion points." << endl;
420  }
421 
422  const polyMesh& mesh = topoChanger().mesh();
423 
424  // Get point from the cut zone
425  const labelList& cutPoints = mesh.pointZones()[cutPointZoneID_.index()];
426 
427  if (cutPoints.size() && !projectedSlavePointsPtr_)
428  {
429  return;
430  }
431  else
432  {
433  const pointField& projectedSlavePoints = *projectedSlavePointsPtr_;
434 
435  const Map<label>& rpm = retiredPointMap();
436 
437  const Map<Pair<edge>>& cpepm = cutPointEdgePairMap();
438 
439  const Map<label>& slaveZonePointMap =
440  mesh.faceZones()[slaveFaceZoneID_.index()]().meshPointMap();
441 
442  const primitiveFacePatch& masterPatch =
443  mesh.faceZones()[masterFaceZoneID_.index()]();
444  const edgeList& masterEdges = masterPatch.edges();
445  const pointField& masterLocalPoints = masterPatch.localPoints();
446 
447  const primitiveFacePatch& slavePatch =
448  mesh.faceZones()[slaveFaceZoneID_.index()]();
449  const edgeList& slaveEdges = slavePatch.edges();
450  const pointField& slaveLocalPoints = slavePatch.localPoints();
451  const vectorField& slavePointNormals = slavePatch.pointNormals();
452 
453  forAll(cutPoints, pointi)
454  {
455  // Try to find the cut point in retired points
456  Map<label>::const_iterator rpmIter = rpm.find(cutPoints[pointi]);
457 
458  if (rpmIter != rpm.end())
459  {
460  if (debug)
461  {
462  Pout<< "p";
463  }
464 
465  // Cut point is a retired point
466  motionPoints[cutPoints[pointi]] =
467  projectedSlavePoints[slaveZonePointMap.find(rpmIter())()];
468  }
469  else
470  {
471  // A cut point is not a projected slave point. Therefore, it
472  // must be an edge-to-edge intersection.
473 
474  Map<Pair<edge>>::const_iterator cpepmIter =
475  cpepm.find(cutPoints[pointi]);
476 
477  if (cpepmIter != cpepm.end())
478  {
479  // Pout<< "Need to re-create hit for point "
480  // << cutPoints[pointi]
481  // << " lookup: " << cpepmIter()
482  // << endl;
483 
484  // Note.
485  // The edge cutting code is repeated in
486  // slidingInterface::coupleInterface. This is done for
487  // efficiency reasons and avoids multiple creation of
488  // cutting planes. Please update both simultaneously.
489  //
490  const edge& globalMasterEdge = cpepmIter().first();
491 
492  const label curMasterEdgeIndex =
493  masterPatch.whichEdge
494  (
495  edge
496  (
497  masterPatch.whichPoint
498  (
499  globalMasterEdge.start()
500  ),
501  masterPatch.whichPoint
502  (
503  globalMasterEdge.end()
504  )
505  )
506  );
507 
508  const edge& cme = masterEdges[curMasterEdgeIndex];
509 
510  // Pout<< "curMasterEdgeIndex: " << curMasterEdgeIndex
511  // << " cme: " << cme
512  // << endl;
513 
514  const edge& globalSlaveEdge = cpepmIter().second();
515 
516  const label curSlaveEdgeIndex =
517  slavePatch.whichEdge
518  (
519  edge
520  (
521  slavePatch.whichPoint
522  (
523  globalSlaveEdge.start()
524  ),
525  slavePatch.whichPoint
526  (
527  globalSlaveEdge.end()
528  )
529  )
530  );
531 
532  const edge& curSlaveEdge = slaveEdges[curSlaveEdgeIndex];
533  // Pout<< "curSlaveEdgeIndex: " << curSlaveEdgeIndex
534  // << " curSlaveEdge: " << curSlaveEdge
535  // << endl;
536  const point& a = projectedSlavePoints[curSlaveEdge.start()];
537  const point& b = projectedSlavePoints[curSlaveEdge.end()];
538 
539  point c =
540  0.5*
541  (
542  slaveLocalPoints[curSlaveEdge.start()]
543  + slavePointNormals[curSlaveEdge.start()]
544  + slaveLocalPoints[curSlaveEdge.end()]
545  + slavePointNormals[curSlaveEdge.end()]
546  );
547 
548  // Create the plane
549  plane cutPlane(a, b, c);
550 
551  linePointRef curSlaveLine =
552  curSlaveEdge.line(slaveLocalPoints);
553  const scalar curSlaveLineMag = curSlaveLine.mag();
554 
555  scalar cutOnMaster =
556  cutPlane.lineIntersect
557  (
558  cme.line(masterLocalPoints)
559  );
560 
561  if
562  (
563  cutOnMaster > edgeEndCutoffTol_
564  && cutOnMaster < 1.0 - edgeEndCutoffTol_
565  )
566  {
567  // Master is cut, check the slave
568  point masterCutPoint =
569  masterLocalPoints[cme.start()]
570  + cutOnMaster*cme.vec(masterLocalPoints);
571 
572  pointHit slaveCut =
573  curSlaveLine.nearestDist(masterCutPoint);
574 
575  if (slaveCut.hit())
576  {
577  // Strict checking of slave cut to avoid capturing
578  // end points.
579  scalar cutOnSlave =
580  (
581  (
582  slaveCut.hitPoint()
583  - curSlaveLine.start()
584  ) & curSlaveLine.vec()
585  )/sqr(curSlaveLineMag);
586 
587  // Calculate merge tolerance from the
588  // target edge length
589  scalar mergeTol =
590  edgeCoPlanarTol_*mag(b - a);
591 
592  if
593  (
594  cutOnSlave > edgeEndCutoffTol_
595  && cutOnSlave < 1.0 - edgeEndCutoffTol_
596  && slaveCut.distance() < mergeTol
597  )
598  {
599  // Cut both master and slave.
600  motionPoints[cutPoints[pointi]] =
601  masterCutPoint;
602  }
603  }
604  else
605  {
606  Pout<< "Missed slave edge!!! This is an error. "
607  << "Master edge: "
608  << cme.line(masterLocalPoints)
609  << " slave edge: " << curSlaveLine
610  << " point: " << masterCutPoint
611  << " weight: " <<
612  (
613  (
614  slaveCut.missPoint()
615  - curSlaveLine.start()
616  ) & curSlaveLine.vec()
617  )/sqr(curSlaveLineMag)
618  << endl;
619  }
620  }
621  else
622  {
623  Pout<< "Missed master edge!!! This is an error"
624  << endl;
625  }
626  }
627  else
628  {
630  << "Cut point " << cutPoints[pointi]
631  << " not recognised as either the projected "
632  << "or as intersection point. Error in point "
633  << "projection or data mapping"
634  << abort(FatalError);
635  }
636  }
637  }
638  if (debug)
639  {
640  Pout<< endl;
641  }
642  }
643 }
644 
645 
647 {
648  if (debug)
649  {
650  Pout<< "void slidingInterface::topoChange(const polyTopoChangeMap& m)"
651  << " const for object " << name() << " : "
652  << "Updating topology." << endl;
653  }
654 
655  // Mesh has changed topologically. Update local topological data
656  const polyMesh& mesh = topoChanger().mesh();
657 
658  masterFaceZoneID_.update(mesh.faceZones());
659  slaveFaceZoneID_.update(mesh.faceZones());
660  cutPointZoneID_.update(mesh.pointZones());
661  cutFaceZoneID_.update(mesh.faceZones());
662 
663  masterPatchID_.update(mesh.boundaryMesh());
664  slavePatchID_.update(mesh.boundaryMesh());
665 
666 //MJ:Disabled updating
667 // if (!attached())
668 // {
669 // calcAttachedAddressing();
670 // }
671 // else
672 // {
673 // renumberAttachedAddressing(m);
674 // }
675 }
676 
677 
679 {
680  if (!projectedSlavePointsPtr_)
681  {
682  projectPoints();
683  }
684 
685  return *projectedSlavePointsPtr_;
686 }
687 
689 {
690  pointMergeTol_ = dict.lookupOrDefault<scalar>
691  (
692  "pointMergeTol",
693  pointMergeTol_
694  );
695  edgeMergeTol_ = dict.lookupOrDefault<scalar>
696  (
697  "edgeMergeTol",
698  edgeMergeTol_
699  );
700  nFacesPerSlaveEdge_ = dict.lookupOrDefault<label>
701  (
702  "nFacesPerSlaveEdge",
703  nFacesPerSlaveEdge_
704  );
705  edgeFaceEscapeLimit_ = dict.lookupOrDefault<label>
706  (
707  "edgeFaceEscapeLimit",
708  edgeFaceEscapeLimit_
709  );
710  integralAdjTol_ = dict.lookupOrDefault<scalar>
711  (
712  "integralAdjTol",
713  integralAdjTol_
714  );
715  edgeMasterCatchFraction_ = dict.lookupOrDefault<scalar>
716  (
717  "edgeMasterCatchFraction",
718  edgeMasterCatchFraction_
719  );
720  edgeCoPlanarTol_ = dict.lookupOrDefault<scalar>
721  (
722  "edgeCoPlanarTol",
723  edgeCoPlanarTol_
724  );
725  edgeEndCutoffTol_ = dict.lookupOrDefault<scalar>
726  (
727  "edgeEndCutoffTol",
728  edgeEndCutoffTol_
729  );
730 
731  if (report)
732  {
733  Info<< "Sliding interface parameters:" << nl
734  << "pointMergeTol : " << pointMergeTol_ << nl
735  << "edgeMergeTol : " << edgeMergeTol_ << nl
736  << "nFacesPerSlaveEdge : " << nFacesPerSlaveEdge_ << nl
737  << "edgeFaceEscapeLimit : " << edgeFaceEscapeLimit_ << nl
738  << "integralAdjTol : " << integralAdjTol_ << nl
739  << "edgeMasterCatchFraction : " << edgeMasterCatchFraction_ << nl
740  << "edgeCoPlanarTol : " << edgeCoPlanarTol_ << nl
741  << "edgeEndCutoffTol : " << edgeEndCutoffTol_ << endl;
742  }
743 }
744 
745 
747 {
748  os << nl << type() << nl
749  << name()<< nl
750  << masterFaceZoneID_.name() << nl
751  << slaveFaceZoneID_.name() << nl
752  << cutPointZoneID_.name() << nl
753  << cutFaceZoneID_.name() << nl
754  << masterPatchID_.name() << nl
755  << slavePatchID_.name() << nl
756  << typeOfMatchNames_[matchType_] << nl
757  << coupleDecouple_ << nl
758  << attached_ << endl;
759 }
760 
761 
762 // To write out all those tolerances
763 #define WRITE_NON_DEFAULT(name) \
764  if ( name ## _ != name ## Default_ )\
765  { \
766  os << " " #name " " << name ## _ << token::END_STATEMENT << nl; \
767  }
768 
769 
771 {
772  os << nl << name() << nl << token::BEGIN_BLOCK << nl
773  << " type " << type() << token::END_STATEMENT << nl
774  << " masterFaceZoneName " << masterFaceZoneID_.name()
776  << " slaveFaceZoneName " << slaveFaceZoneID_.name()
778  << " cutPointZoneName " << cutPointZoneID_.name()
780  << " cutFaceZoneName " << cutFaceZoneID_.name()
782  << " masterPatchName " << masterPatchID_.name()
784  << " slavePatchName " << slavePatchID_.name()
786  << " typeOfMatch " << typeOfMatchNames_[matchType_]
788  << " coupleDecouple " << coupleDecouple_
790  << " projection " << intersection::algorithmNames_[projectionAlgo_]
792  << " attached " << attached_
794  << " active " << active()
795  << token::END_STATEMENT << nl;
796 
797  if (attached_)
798  {
799  writeEntry(os, "masterFaceCells", *masterFaceCellsPtr_);
800  writeEntry(os, "slaveFaceCells", *slaveFaceCellsPtr_);
801  writeEntry(os, "masterStickOutFaces", *masterStickOutFacesPtr_);
802  writeEntry(os, "slaveStickOutFaces", *slaveStickOutFacesPtr_);
803 
804  os << " retiredPointMap " << retiredPointMap()
805  << token::END_STATEMENT << nl
806  << " cutPointEdgePairMap " << cutPointEdgePairMap()
807  << token::END_STATEMENT << nl;
808  }
809 
810  WRITE_NON_DEFAULT(pointMergeTol)
811  WRITE_NON_DEFAULT(edgeMergeTol)
812  WRITE_NON_DEFAULT(nFacesPerSlaveEdge)
813  WRITE_NON_DEFAULT(edgeFaceEscapeLimit)
814  WRITE_NON_DEFAULT(integralAdjTol)
815  WRITE_NON_DEFAULT(edgeMasterCatchFraction)
816  WRITE_NON_DEFAULT(edgeCoPlanarTol)
817  WRITE_NON_DEFAULT(edgeEndCutoffTol)
818 
819  os << token::END_BLOCK << endl;
820 }
821 
822 
823 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:453
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
A line primitive.
Definition: line.H:56
bool update()
Update the mesh for topology change, mesh to mesh mapping.
Definition: fvMesh.C:584
const faceZoneID & slaveFaceZoneID() const
Return slave face zone ID.
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
scalar mag() const
Return scalar magnitude.
Definition: lineI.H:80
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
static const NamedEnum< algorithm, 3 > algorithmNames_
Projection algorithm names.
Definition: intersection.H:80
virtual void modifyMotionPoints(pointField &motionPoints) const
Modify motion points to comply with the topological change.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
bool empty() const
Return true if the UPtrList is empty (ie, size() is zero)
Definition: UPtrListI.H:36
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
Geometric class that creates a 2D plane and can return the intersection point between a line and the ...
Definition: plane.H:60
virtual void write(Ostream &) const
Write.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
fvMesh & mesh
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:51
const dimensionedScalar c
Speed of light in a vacuum.
virtual void setRefinement(polyTopoChange &) const
Insert the layer addition/removal instructions.
const Point & missPoint() const
Return miss point.
Definition: PointHit.H:145
Macros for easy insertion into run-time selection tables.
const pointField & pointProjection() const
Return projected points for a slave patch.
const meshPointZones & pointZones() const
Return point zones.
Definition: polyMesh.H:489
label whichEdge(const edge &) const
Given an edge in local point labels, return its.
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
A list of faces which address into the list of points.
bool read(const char *, int32_t &)
Definition: int32IO.C:85
linePointRef line(const pointField &) const
Return edge line.
Definition: edgeI.H:187
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
List of mesh modifiers defining the mesh dynamics.
A class for handling words, derived from string.
Definition: word.H:59
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
#define WRITE_NON_DEFAULT(name)
const Point & hitPoint() const
Return hit point.
Definition: PointHit.H:126
PointRef start() const
Return first vertex.
Definition: lineI.H:60
bool hit() const
Is there a hit.
Definition: PointHit.H:120
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
List< label > labelList
A List of labels.
Definition: labelList.H:56
slidingInterface(const word &name, const label index, const polyTopoChanger &mme, const word &masterFaceZoneName, const word &slaveFaceZoneName, const word &cutPointZoneName, const word &cutFaceZoneName, const word &masterPatchName, const word &slavePatchName, const typeOfMatch tom, const bool coupleDecouple=false, const intersection::algorithm algo=intersection::algorithm::visible)
Construct from components.
Virtual base class for mesh modifiers.
errorManip< error > abort(error &err)
Definition: errorManip.H:131
label whichPoint(const label gp) const
Given a global point index, return the local point index.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
vector vec(const pointField &) const
Return the vector (end - start)
Definition: edgeI.H:175
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
static const char nl
Definition: Ostream.H:260
defineTypeNameAndDebug(combustionModel, 0)
static const NamedEnum< typeOfMatch, 2 > typeOfMatchNames_
Direction names.
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
typeOfMatch
Type of match.
scalar lineIntersect(const line< Point, PointRef > &l) const
Return the cutting point between the plane and.
Definition: plane.H:186
const Field< PointType > & pointNormals() const
Return point normals for patch.
Point vec() const
Return start-end vector.
Definition: lineI.H:87
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
label end() const
Return end vertex label.
Definition: edgeI.H:92
prefixOSstream Pout(cout, "Pout")
Definition: IOstreams.H:53
const meshFaceZones & faceZones() const
Return face zones.
Definition: polyMesh.H:495
Direct mesh changes based on v1.3 polyTopoChange syntax.
virtual void topoChange(const polyTopoChangeMap &)
Force recalculation of locally stored data on topological change.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
bool update(const point &, const FvWallInfoDataBase< WallInfo, Type, Derived > &w2, const scalar tol, TrackingData &td)
Evaluate distance to point. Update distSqr, origin from whomever.
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
const faceZoneID & masterFaceZoneID() const
Return master face zone ID.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
void setTolerances(const dictionary &, bool report=false)
Set the tolerances from the values in a dictionary.
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
rDeltaT ref()
const polyMesh & mesh() const
Return the mesh reference.
void deleteDemandDrivenData(DataPtr &dataPtr)
virtual ~slidingInterface()
Destructor.
virtual void writeDict(Ostream &) const
Write dictionary.
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95
virtual bool changeTopology() const
Check for topology change.
label start() const
Return start vertex label.
Definition: edgeI.H:81
Namespace for OpenFOAM.
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:864