nonConformalCyclicPolyPatch.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) 2021-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 
29 #include "polyBoundaryMesh.H"
30 #include "polyMesh.H"
31 #include "SubField.H"
32 #include "nonConformalBoundary.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
39 
42  (
43  polyPatch,
46  );
47 }
48 
49 
50 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
51 
53 {
55  intersectionIsValid_ = false;
56  raysIsValid_ = false;
57 }
58 
59 
61 {
62  static_cast<cyclicTransform&>(*this) =
64  (
65  name(),
66  origPatch().faceAreas(),
67  *this,
68  nbrPatchName(),
69  nbrPatch(),
70  matchTolerance()
71  );
72 }
73 
74 
76 (
77  PstreamBuffers& pBufs,
78  const pointField& p
79 )
80 {
82  intersectionIsValid_ = false;
83  raysIsValid_ = false;
84 }
85 
86 
88 {
90  intersectionIsValid_ = false;
91  raysIsValid_ = false;
92 }
93 
94 
96 {
98  intersectionIsValid_ = false;
99  raysIsValid_ = false;
100 }
101 
102 
104 {
105  cyclicPolyPatch::rename(newNames);
107 }
108 
109 
111 {
112  cyclicPolyPatch::reorder(newToOldIndex);
114 }
115 
116 
117 // * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * * //
118 
120 (
121  const word& name,
122  const label size,
123  const label start,
124  const label index,
125  const polyBoundaryMesh& bm,
126  const word& patchType
127 )
128 :
129  cyclicPolyPatch(name, size, start, index, bm, patchType),
130  nonConformalCoupledPolyPatch(static_cast<const polyPatch&>(*this)),
131  intersectionIsValid_(false),
132  intersection_(false),
133  raysIsValid_(false),
134  rays_(false)
135 {}
136 
137 
139 (
140  const word& name,
141  const label size,
142  const label start,
143  const label index,
144  const polyBoundaryMesh& bm,
145  const word& patchType,
146  const word& nbrPatchName,
147  const word& origPatchName,
149 )
150 :
152  (
153  name,
154  size,
155  start,
156  index,
157  bm,
158  patchType,
159  nbrPatchName,
160  transform
161  ),
162  nonConformalCoupledPolyPatch(*this, origPatchName),
163  intersectionIsValid_(false),
164  intersection_(false),
165  raysIsValid_(false),
166  rays_(false)
167 {}
168 
169 
171 (
172  const word& name,
173  const dictionary& dict,
174  const label index,
175  const polyBoundaryMesh& bm,
176  const word& patchType
177 )
178 :
179  cyclicPolyPatch(name, dict, index, bm, patchType, true),
181  intersectionIsValid_(false),
182  intersection_(false),
183  raysIsValid_(false),
184  rays_(false)
185 {}
186 
187 
189 (
190  const nonConformalCyclicPolyPatch& pp,
191  const polyBoundaryMesh& bm
192 )
193 :
194  cyclicPolyPatch(pp, bm),
195  nonConformalCoupledPolyPatch(*this, pp),
196  intersectionIsValid_(false),
197  intersection_(false),
198  raysIsValid_(false),
199  rays_(false)
200 {}
201 
202 
204 (
205  const nonConformalCyclicPolyPatch& pp,
206  const polyBoundaryMesh& bm,
207  const label index,
208  const label newSize,
209  const label newStart,
210  const word& nbrPatchName,
211  const word& origPatchName
212 )
213 :
214  cyclicPolyPatch(pp, bm, index, newSize, newStart, nbrPatchName),
215  nonConformalCoupledPolyPatch(*this, origPatchName),
216  intersectionIsValid_(false),
217  intersection_(false),
218  raysIsValid_(false),
219  rays_(false)
220 {}
221 
222 
223 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
224 
226 {}
227 
228 
229 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
230 
233 {
234  return
235  refCast<const nonConformalCyclicPolyPatch>
236  (
238  );
239 }
240 
241 
243 {
244  return false;
245 }
246 
247 
250 {
251  if (!owner())
252  {
254  << "The non-conformal cyclic intersection is only available to "
255  << "the owner patch" << abort(FatalError);
256  }
257 
258  if (!intersectionIsValid_)
259  {
260  const polyMesh& mesh = boundaryMesh().mesh();
261 
263 
264  intersection_.update
265  (
266  origPatch(),
267  ncb.patchPointNormals(origPatchIndex()),
268  nbrPatch().origPatch(),
269  transform()
270  );
271 
272  intersectionIsValid_ = true;
273  }
274 
275  return intersection_;
276 }
277 
278 
281 {
282  if (!owner())
283  {
285  << "The non-conformal cyclic rays is only available to "
286  << "the owner patch" << abort(FatalError);
287  }
288 
289  if (!raysIsValid_)
290  {
291  const polyMesh& mesh = boundaryMesh().mesh();
292 
294 
295  rays_.update
296  (
298  (
299  origPatch(),
300  mesh.points(),
301  mesh.oldPoints()
302  ),
303  ncb.patchPointNormals(origPatchIndex()),
304  ncb.patchPointNormals0(origPatchIndex()),
306  (
307  nbrPatch().origPatch(),
308  mesh.points(),
309  mesh.oldPoints()
310  ),
311  transform()
312  );
313 
314  raysIsValid_ = true;
315  }
316 
317  return rays_;
318 }
319 
320 
322 (
323  const scalar fraction,
324  const label origFacei,
325  const vector& p,
326  const vector& n,
327  point& nbrP
328 ) const
329 {
330  const polyMesh& mesh = boundaryMesh().mesh();
331 
332  const nonConformalCyclicPolyPatch& ownerPatch =
333  owner() ? *this : nbrPatch();
334 
335  auto ownerRaysMethod =
336  owner()
339 
340  return
341  (ownerPatch.rays().*ownerRaysMethod)
342  (
344  (
345  nbrPatch().origPatch(),
346  mesh.points(),
347  mesh.oldPoints()
348  ),
349  fraction,
350  origFacei,
351  transform().invTransformPosition(p),
352  transform().invTransform(n),
353  nbrP
354  );
355 }
356 
357 
359 {
362 }
363 
364 
365 // ************************************************************************* //
label n
Macros for easy insertion into run-time selection tables.
static nonConformalBoundary & New(const word &name, const polyMesh &mesh)
Construct and return the named DemandDrivenMeshObject.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Cyclic plane patch.
const cyclicPolyPatch & nbrPatch() const
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void initMovePoints(PstreamBuffers &, const pointField &)
Initialise the patches for moving points.
virtual void rename(const wordList &newNames)
Reset the patch name.
virtual void reorder(const labelUList &newToOldIndex)
Reset the patch index.
virtual void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.
Cyclic plane transformation.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Mesh object that stores an all boundary patch and mapping to and from it and the mesh and the individ...
tmp< vectorField > patchPointNormals(const label patchi) const
Get parallel consistent point normals for the patch.
tmp< vectorField > patchPointNormals0(const label patchi) const
Get parallel consistent old-time point normals for the patch.
Non-conformal coupled poly patch. As nonConformalPolyPatch, but this patch is coupled to another non-...
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
virtual void rename(const wordList &newNames)
Reset the patch name.
virtual void reorder(const labelUList &newToOldIndex)
Reset the patch index.
Non-conformal cyclic poly patch. As nonConformalCoupledPolyPatch, but the neighbouring patch is local...
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
remote ray(const scalar fraction, const label origFacei, const vector &p, const vector &n, point &nbrP) const
Compute a ray intersection across the coupling.
virtual void rename(const wordList &newNames)
Reset the patch name.
virtual void clearGeom()
Clear geometry.
bool intersectionIsValid_
Is the intersection engine up to date?
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
const patchToPatches::rays & rays() const
Access the rays engine.
virtual void reorder(const labelUList &newToOldIndex)
Reset the patch index.
nonConformalCyclicPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType)
Construct from components.
virtual void initCalcGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
virtual bool coupled() const
Is this patch coupled? Returns false. For NCC patches the poly.
const patchToPatches::intersection & intersection() const
Access the intersection engine.
bool raysIsValid_
Is the rays engine up to date?
const nonConformalCyclicPolyPatch & nbrPatch() const
Neighbour patch.
virtual void initTopoChange(PstreamBuffers &)
Initialise the update of the patch topology.
Class to generate patchToPatch coupling geometry. A full geometric intersection is done between a fac...
remote srcToTgtRay(const primitiveOldTimePatch &tgtPatch, const scalar fraction, const label srcFacei, const vector &srcP, const vector &srcN, point &tgtP) const
Compute a ray intersection from the source side to the target.
remote tgtToSrcRay(const primitiveOldTimePatch &srcPatch, const scalar fraction, const label tgtFacei, const vector &tgtP, const vector &tgtN, point &srcP) const
Compute a ray intersection from the target side to the source.
Foam::polyBoundaryMesh.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
virtual const pointField & oldPoints() const
Return old points for mesh motion.
Definition: polyMesh.C:1351
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1313
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:70
virtual void clearGeom()
Clear geometry.
Definition: polyPatch.C:76
Struct for keeping processor, element (cell, face, point) index.
Definition: remote.H:57
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
Namespace for OpenFOAM.
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
addToRunTimeSelectionTable(polyPatch, mergedCyclicPolyPatch, word)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
errorManip< error > abort(error &err)
Definition: errorManip.H:131
PrimitiveOldTimePatch< SubList< face >, const pointField & > primitiveOldTimePatch
Addressing for a faceList slice.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
defineTypeNameAndDebug(combustionModel, 0)
error FatalError
dictionary dict
volScalarField & p