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