cyclicAMIPolyPatch.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 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 Class
25  Foam::cyclicAMIPolyPatch
26 
27 Description
28  Cyclic patch for Arbitrary Mesh Interface (AMI)
29 
30 SourceFiles
31  cyclicAMIPolyPatch.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef cyclicAMIPolyPatch_H
36 #define cyclicAMIPolyPatch_H
37 
38 #include "coupledPolyPatch.H"
40 #include "polyBoundaryMesh.H"
41 #include "coupleGroupIdentifier.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 /*---------------------------------------------------------------------------*\
49  Class cyclicAMIPolyPatch Declaration
50 \*---------------------------------------------------------------------------*/
51 
53 :
54  public coupledPolyPatch
55 {
56 
57 private:
58 
59  // Private Member Functions
60 
61  //- Return normal of face at max distance from rotation axis
62  vector findFaceNormalMaxRadius(const pointField& faceCentres) const;
63 
64  void calcTransforms
65  (
66  const primitivePatch& half0,
67  const pointField& half0Ctrs,
68  const vectorField& half0Areas,
69  const pointField& half1Ctrs,
70  const vectorField& half1Areas
71  );
72 
73 
74 protected:
75 
76  // Protected data
77 
78  //- Name of other half
79  mutable word nbrPatchName_;
80 
81  //- Optional patchGroup to find neighbPatch
83 
84  //- Index of other half
85  mutable label nbrPatchID_;
86 
87 
88  // Transformations
89 
90  // For rotation
91 
92  //- Axis of rotation for rotational cyclics
94 
95  //- Point on axis of rotation for rotational cyclics
97 
98  //- Flag to show whether the rotation angle is defined
100 
101  //- Rotation angle
102  scalar rotationAngle_;
103 
104 
105  // For translation
106 
107  //- Translation vector
109 
110 
111  //- AMI interpolation class
113 
114  //- Flag to indicate that slave patch should be reversed for AMI
115  const bool AMIReverse_;
116 
117  //- Flag to indicate that patches should match/overlap
118  bool AMIRequireMatch_;
119 
120  //- Low weight correction threshold for AMI
121  const scalar AMILowWeightCorrection_;
122 
123  //- Projection surface
125 
126  //- Dictionary used during projection surface construction
127  const dictionary surfDict_;
128 
129 
130  // Protected Member Functions
131 
132  //- Reset the AMI interpolator
133  virtual void resetAMI
134  (
137  ) const;
138 
139  //- Recalculate the transformation tensors
140  virtual void calcTransforms();
141 
142  //- Initialise the calculation of the patch geometry
143  virtual void initGeometry(PstreamBuffers&);
144 
145  //- Calculate the patch geometry
146  virtual void calcGeometry(PstreamBuffers&);
147 
148  //- Initialise the patches for moving points
149  virtual void initMovePoints(PstreamBuffers& pBufs, const pointField&);
150 
151  //- Correct patches after moving points
152  virtual void movePoints(PstreamBuffers& pBufs, const pointField&);
153 
154  //- Initialise the update of the patch topology
155  virtual void initUpdateMesh(PstreamBuffers&);
156 
157  //- Update of the patch topology
158  virtual void updateMesh(PstreamBuffers&);
159 
160  //- Clear geometry
161  virtual void clearGeom();
162 
163 
164 public:
165 
166  //- Runtime type information
167  TypeName("cyclicAMI");
168 
169 
170  // Constructors
171 
172  //- Construct from (base couped patch) components
174  (
175  const word& name,
176  const label size,
177  const label start,
178  const label index,
179  const polyBoundaryMesh& bm,
180  const word& patchType,
182  );
183 
184  //- Construct from dictionary
186  (
187  const word& name,
188  const dictionary& dict,
189  const label index,
190  const polyBoundaryMesh& bm,
191  const word& patchType
192  );
193 
194  //- Construct as copy, resetting the boundary mesh
196 
197  //- Construct given the original patch and resetting the
198  // face list and boundary mesh information
200  (
201  const cyclicAMIPolyPatch& pp,
202  const polyBoundaryMesh& bm,
203  const label index,
204  const label newSize,
205  const label newStart,
206  const word& nbrPatchName
207  );
208 
209  //- Construct given the original patch and a map
211  (
212  const cyclicAMIPolyPatch& pp,
213  const polyBoundaryMesh& bm,
214  const label index,
215  const labelUList& mapAddressing,
216  const label newStart
217  );
218 
219 
220  //- Construct and return a clone, resetting the boundary mesh
221  virtual autoPtr<polyPatch> clone(const polyBoundaryMesh& bm) const
222  {
223  return autoPtr<polyPatch>(new cyclicAMIPolyPatch(*this, bm));
224  }
225 
226  //- Construct and return a clone, resetting the face list
227  // and boundary mesh
229  (
230  const polyBoundaryMesh& bm,
231  const label index,
232  const label newSize,
233  const label newStart
234  ) const
235  {
236  return autoPtr<polyPatch>
237  (
239  (
240  *this,
241  bm,
242  index,
243  newSize,
244  newStart,
245  nbrPatchName_
246  )
247  );
248  }
249 
250  //- Construct and return a clone, resetting the face list
251  // and boundary mesh
253  (
254  const polyBoundaryMesh& bm,
255  const label index,
256  const labelUList& mapAddressing,
257  const label newStart
258  ) const
259  {
260  return autoPtr<polyPatch>
261  (
263  (
264  *this,
265  bm,
266  index,
267  mapAddressing,
268  newStart
269  )
270  );
271  }
272 
273 
274  //- Destructor
275  virtual ~cyclicAMIPolyPatch();
276 
277 
278  // Member Functions
279 
280  // Access
281 
282  //- Is patch 'coupled'. Note that on AMI the geometry is not
283  // coupled but the fields are!
284  virtual bool coupled() const
285  {
286  return false;
287  }
288 
289  //- Neighbour patch name
290  inline const word& neighbPatchName() const;
291 
292  //- Neighbour patch ID
293  virtual label neighbPatchID() const;
294 
295  //- Does this side own the patch?
296  virtual bool owner() const;
297 
298  //- Return a reference to the neighbour patch
299  virtual const cyclicAMIPolyPatch& neighbPatch() const;
300 
301  //- Return a reference to the projection surface
302  const autoPtr<searchableSurface>& surfPtr() const;
303 
304  //- Return a reference to the AMI interpolator
305  const AMIPatchToPatchInterpolation& AMI() const;
306 
307  //- Return true if applying the low weight correction
308  bool applyLowWeightCorrection() const;
309 
310 
311  // Transformations
312 
313  //- Axis of rotation for rotational cyclic AMI
314  inline const vector& rotationAxis() const;
315 
316  //- Point on axis of rotation for rotational cyclic AMI
317  inline const point& rotationCentre() const;
318 
319  //- Translation vector for translational cyclic AMI
320  inline const vector& separationVector() const;
321 
322  //- Transform patch-based positions from nbr side to this side
323  virtual void transformPosition(pointField&) const;
324 
325  //- Transform a patch-based position from nbr side to this side
326  virtual void transformPosition
327  (
328  point& l,
329  const label faceI
330  ) const;
331 
332  //- Transform a patch-based position from this side to nbr side
333  virtual void reverseTransformPosition
334  (
335  point& l,
336  const label faceI
337  ) const;
338 
339  //- Transform a patch-based direction from this side to nbr side
340  virtual void reverseTransformDirection
341  (
342  vector& d,
343  const label faceI
344  ) const;
345 
346 
347  // Interpolations
348 
349  //- Interpolate field
350  template<class Type>
352  (
353  const Field<Type>& fld,
354  const UList<Type>& defaultValues = UList<Type>()
355  ) const;
356 
357  //- Interpolate tmp field
358  template<class Type>
360  (
361  const tmp<Field<Type> >& tFld,
362  const UList<Type>& defaultValues = UList<Type>()
363  ) const;
364 
365  //- Low-level interpolate List
366  template<class Type, class CombineOp>
367  void interpolate
368  (
369  const UList<Type>& fld,
370  const CombineOp& cop,
371  List<Type>& result,
372  const UList<Type>& defaultValues = UList<Type>()
373  ) const;
374 
375 
376  //- Calculate the patch geometry
377  virtual void calcGeometry
378  (
379  const primitivePatch& referPatch,
380  const pointField& thisCtrs,
381  const vectorField& thisAreas,
382  const pointField& thisCc,
383  const pointField& nbrCtrs,
384  const vectorField& nbrAreas,
385  const pointField& nbrCc
386  );
387 
388  //- Initialize ordering for primitivePatch. Does not
389  // refer to *this (except for name() and type() etc.)
390  virtual void initOrder
391  (
393  const primitivePatch&
394  ) const;
395 
396  //- Return new ordering for primitivePatch.
397  // Ordering is -faceMap: for every face
398  // index of the new face -rotation:for every new face the clockwise
399  // shift of the original face. Return false if nothing changes
400  // (faceMap is identity, rotation is 0), true otherwise.
401  virtual bool order
402  (
404  const primitivePatch&,
406  labelList& rotation
407  ) const;
408 
409  //- Return face index on neighbour patch which shares point p
410  // following trajectory vector n
412  (
413  const label faceI,
414  const vector& n,
415  point& p
416  ) const;
417 
418  //- Write the polyPatch data as a dictionary
419  virtual void write(Ostream&) const;
420 };
421 
422 
423 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
424 
425 } // End namespace Foam
426 
427 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
428 
429 #include "cyclicAMIPolyPatchI.H"
430 
431 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
432 
433 #ifdef NoRepository
435 #endif
436 
437 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
438 
439 #endif
440 
441 // ************************************************************************* //
virtual bool owner() const
Does this side own the patch?
virtual void updateMesh(PstreamBuffers &)
Update of the patch topology.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
bool rotationAngleDefined_
Flag to show whether the rotation angle is defined.
const autoPtr< searchableSurface > & surfPtr() const
Return a reference to the projection surface.
Cyclic patch for Arbitrary Mesh Interface (AMI)
const point & rotationCentre() const
Point on axis of rotation for rotational cyclic AMI.
virtual void calcGeometry(PstreamBuffers &)
Calculate the patch geometry.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
virtual void reverseTransformPosition(point &l, const label faceI) const
Transform a patch-based position from this side to nbr side.
const word & name() const
Return name.
virtual void clearGeom()
Clear geometry.
virtual void initMovePoints(PstreamBuffers &pBufs, const pointField &)
Initialise the patches for moving points.
vector separationVector_
Translation vector.
virtual void initGeometry(PstreamBuffers &)
Initialise the calculation of the patch geometry.
autoPtr< AMIPatchToPatchInterpolation > AMIPtr_
AMI interpolation class.
Base class for Arbitrary Mesh Interface (AMI) methods.
Definition: AMIMethod.H:55
virtual void initUpdateMesh(PstreamBuffers &)
Initialise the update of the patch topology.
cyclicAMIPolyPatch(const word &name, const label size, const label start, const label index, const polyBoundaryMesh &bm, const word &patchType, const transformType transform=UNKNOWN)
Construct from (base couped patch) components.
virtual ~cyclicAMIPolyPatch()
Destructor.
TypeName("cyclicAMI")
Runtime type information.
label index() const
Return the index of this patch in the boundaryMesh.
A class for handling words, derived from string.
Definition: word.H:59
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
virtual void movePoints(PstreamBuffers &pBufs, const pointField &)
Correct patches after moving points.
autoPtr< searchableSurface > surfPtr_
Projection surface.
virtual label neighbPatchID() const
Neighbour patch ID.
bool AMIRequireMatch_
Flag to indicate that patches should match/overlap.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
word nbrPatchName_
Name of other half.
Namespace for OpenFOAM.
scalar rotationAngle_
Rotation angle.
label n
interpolationMethod
Enumeration specifying interpolation method.
point rotationCentre_
Point on axis of rotation for rotational cyclics.
dictionary dict
virtual const cyclicAMIPolyPatch & neighbPatch() const
Return a reference to the neighbour patch.
virtual void transformPosition(pointField &) const
Transform patch-based positions from nbr side to this side.
Foam::polyBoundaryMesh.
volScalarField & p
Definition: createFields.H:51
const coupleGroupIdentifier coupleGroup_
Optional patchGroup to find neighbPatch.
virtual void calcTransforms()
Recalculate the transformation tensors.
Interpolation class dealing with transfer of data between two primitive patches with an arbitrary mes...
const scalar AMILowWeightCorrection_
Low weight correction threshold for AMI.
Encapsulates using patchGroups to specify coupled patch.
label pointFace(const label faceI, const vector &n, point &p) const
Return face index on neighbour patch which shares point p.
virtual bool coupled() const
Is patch &#39;coupled&#39;. Note that on AMI the geometry is not.
tmp< Field< Type > > interpolate(const Field< Type > &fld, const UList< Type > &defaultValues=UList< Type >()) const
Interpolate field.
virtual autoPtr< polyPatch > clone(const polyBoundaryMesh &bm) const
Construct and return a clone, resetting the boundary mesh.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
A list of faces which address into the list of points.
bool applyLowWeightCorrection() const
Return true if applying the low weight correction.
label start() const
Return start label of this patch in the polyMesh face list.
Definition: polyPatch.H:300
const vector & separationVector() const
Translation vector for translational cyclic AMI.
const AMIPatchToPatchInterpolation & AMI() const
Return a reference to the AMI interpolator.
virtual void resetAMI(const AMIPatchToPatchInterpolation::interpolationMethod &AMIMethod=AMIPatchToPatchInterpolation::imFaceAreaWeight) const
Reset the AMI interpolator.
virtual transformType transform() const
Type of transform.
const bool AMIReverse_
Flag to indicate that slave patch should be reversed for AMI.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
virtual void initOrder(PstreamBuffers &, const primitivePatch &) const
Initialize ordering for primitivePatch. Does not.
label nbrPatchID_
Index of other half.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
virtual void reverseTransformDirection(vector &d, const label faceI) const
Transform a patch-based direction from this side to nbr side.
const vector & rotationAxis() const
Axis of rotation for rotational cyclic AMI.
virtual void write(Ostream &) const
Write the polyPatch data as a dictionary.
The coupledPolyPatch is an abstract base class for patches that couple regions of the computational d...
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
A class for managing temporary objects.
Definition: PtrList.H:118
vector rotationAxis_
Axis of rotation for rotational cyclics.
const dictionary surfDict_
Dictionary used during projection surface construction.
virtual bool order(PstreamBuffers &, const primitivePatch &, labelList &faceMap, labelList &rotation) const
Return new ordering for primitivePatch.
const word & neighbPatchName() const
Neighbour patch name.