patchToPatch.H
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-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 Class
25  Foam::patchToPatch
26 
27 Description
28  Class to generate coupling geometry between two primitive patches
29 
30 SourceFiles
31  patchToPatch.C
32  patchToPatchParallelOps.C
33  patchToPatchTemplates.C
34 
35 \*---------------------------------------------------------------------------*/
36 
37 #ifndef patchToPatch_H
38 #define patchToPatch_H
39 
40 #include "remote.H"
41 #include "distributionMap.H"
42 #include "primitivePatch.H"
43 #include "primitiveOldTimePatch.H"
44 #include "runTimeSelectionTables.H"
45 #include "treeBoundBox.H"
46 
47 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
48 
49 namespace Foam
50 {
51 
52 /*---------------------------------------------------------------------------*\
53  Class patchToPatch Declaration
54 \*---------------------------------------------------------------------------*/
55 
56 class patchToPatch
57 {
58 protected:
59 
60  // Private Data
61 
62  //- Flag to indicate that the two patches are co-directional and
63  // that the orientation of one should therefore be reversed
64  const bool reverse_;
65 
66  //- Index of the processor holding all faces of the patchToPatch, or -1
67  // if spread across multiple processors
69 
70  //- For each source face, the coupled local target faces
72 
73  //- For each target face, the coupled local source faces
75 
76  //- Map from source patch faces to target-local source patch faces
78 
79  //- Map from target patch faces to source-local target patch faces
81 
82  //- When running in parallel, a map from local source face index to
83  // source processor and face index
85 
86  //- When running in parallel, a map from local target face index to
87  // target processor and face index
89 
90 
91  // Private Member Functions
92 
93  // Searching
94 
95  //- Get the bound box for a source face
96  virtual treeBoundBox srcBox
97  (
98  const face& srcFace,
99  const pointField& srcPoints,
100  const vectorField& srcPointNormals
101  ) const = 0;
102 
103  //- Get the bound box for a source face
105  (
106  const primitiveOldTimePatch& srcPatch,
107  const vectorField& srcPointNormals,
108  const vectorField& srcPointNormals0,
109  const label srcFacei
110  ) const;
111 
112  //- Get the bound box for the source patch
114  (
115  const primitiveOldTimePatch& srcPatch,
116  const vectorField& srcPointNormals,
117  const vectorField& srcPointNormals0
118  ) const;
119 
120  //- Get the bound box for a target face
122  (
123  const primitiveOldTimePatch& tgtPatch,
124  const label tgtFacei
125  ) const;
126 
127  //- Get the bound box for the target patch
129  (
130  const primitiveOldTimePatch& tgtPatch
131  ) const;
132 
133 
134  // Intersection
135 
136  //- Intersect two faces
137  virtual bool intersectFaces
138  (
139  const primitiveOldTimePatch& srcPatch,
140  const vectorField& srcPointNormals,
141  const vectorField& srcPointNormals0,
142  const primitiveOldTimePatch& tgtPatch,
143  const label srcFacei,
144  const label tgtFacei
145  ) = 0;
146 
147  //- Intersect two faces
149  (
150  const primitiveOldTimePatch& srcPatch,
151  const vectorField& srcPointNormals,
152  const vectorField& srcPointNormals0,
153  const primitiveOldTimePatch& tgtPatch,
154  const label srcFacei,
155  const label tgtFacei
156  );
157 
158  //- Intersect a queue of source-target face pairs. Update completion
159  // lists and form a new queue of target-source face pairs.
161  (
162  const primitiveOldTimePatch& srcPatch,
163  const vectorField& srcPointNormals,
164  const vectorField& srcPointNormals0,
165  const primitiveOldTimePatch& tgtPatch,
166  const bool isSrc,
167  const DynamicList<labelPair>& queue,
168  labelList& faceComplete,
169  DynamicList<labelPair>& otherQueue,
170  const labelList& otherFaceComplete,
171  boolList& otherFaceQueued,
172  boolList& otherFaceVisited
173  );
174 
175  //- Intersect the patches
176  void intersectPatches
177  (
178  const primitiveOldTimePatch& srcPatch,
179  const vectorField& srcPointNormals,
180  const vectorField& srcPointNormals0,
181  const primitiveOldTimePatch& tgtPatch
182  );
183 
184 
185  // Parallel functionality
186 
187  //- Determine which target faces need to be sent to the source.
188  // This is done before intersection. Bound boxes are used to
189  // estimate what faces will intersect.
191  (
192  const primitiveOldTimePatch& srcPatch,
193  const vectorField& srcPointNormals,
194  const vectorField& srcPointNormals0,
195  const primitiveOldTimePatch& tgtPatch
196  ) const;
197 
198  //- Distribute a patch given its distribution map
200  (
201  const distributionMap& map,
202  const primitiveOldTimePatch& patch,
204  localPatchPtr
205  );
206 
207 
208  // Hooks
209 
210  //- Initialisation
211  virtual void initialise
212  (
213  const primitiveOldTimePatch& srcPatch,
214  const vectorField& srcPointNormals,
215  const vectorField& srcPointNormals0,
216  const primitiveOldTimePatch& tgtPatch
217  );
218 
219  //- Finalise the intersection locally. Trims the local target patch
220  // so that only parts that actually intersect the source remain.
221  // Returns new-to-old map for trimming target-associated data.
222  virtual labelList finaliseLocal
223  (
224  const primitiveOldTimePatch& srcPatch,
225  const vectorField& srcPointNormals,
226  const vectorField& srcPointNormals0,
227  const primitiveOldTimePatch& tgtPatch
228  );
229 
230  //- Distribute the source patch so that any information needed by
231  // the target is locally available
232  virtual void distributeSrc(const primitiveOldTimePatch& srcPatch);
233 
234  //- Send data that resulted from an intersection between the source
235  // patch and a distributed source-local-target patch back to the
236  // original target processes
237  virtual void rDistributeTgt(const primitiveOldTimePatch& tgtPatch);
238 
239  //- Finalising
240  virtual label finalise
241  (
242  const primitiveOldTimePatch& srcPatch,
243  const vectorField& srcPointNormals,
244  const vectorField& srcPointNormals0,
245  const primitiveOldTimePatch& tgtPatch,
246  const transformer& tgtToSrc
247  );
248 
249 
250  // Interpolation
251 
252  //- For each source face, the coupled target weights
253  virtual tmpNrc<List<DynamicList<scalar>>> srcWeights() const = 0;
254 
255  //- For each target face, the coupled source weights
256  virtual tmpNrc<List<DynamicList<scalar>>> tgtWeights() const = 0;
257 
258 
259 public:
260 
261  //- Runtime type information
262  TypeName("patchToPatch");
263 
264 
265  //- Declare runtime constructor selection table
267  (
268  autoPtr,
269  patchToPatch,
270  bool,
271  (const bool reverse),
272  (reverse)
273  );
274 
275 
276  // Constructors
277 
278  //- Construct from components
279  patchToPatch(const bool reverse);
280 
281  //- Disallow default bitwise copy construction
282  patchToPatch(const patchToPatch&) = delete;
283 
284 
285  //- Destructor
286  virtual ~patchToPatch();
287 
288 
289  // Selector
290 
291  //- Select from name
293  (
294  const word& patchToPatchType,
295  const bool reverse
296  );
297 
298 
299  // Member Functions
300 
301  // Access
302 
303  //- Flag to indicate that the two patches are co-directional and
304  // that the orientation of one should therefore be reversed
305  inline bool reverse() const;
306 
307  //- Index of the processor holding all faces of the patchToPatch,
308  // or -1 if spread across multiple processors
309  inline label singleProcess() const;
310 
311  //- Is this intersection on a single process?
312  inline bool isSingleProcess() const;
313 
314  //- Return a list indicating which source faces are coupled
315  PackedBoolList srcCoupled() const;
316 
317  //- Return a list indicating which target faces are coupled
318  PackedBoolList tgtCoupled() const;
319 
320  //- For each source face, the coupled target procs and faces
322 
323  //- For each target face, the coupled source procs and faces
325 
326 
327  // Interpolation
328 
329  //- Interpolate a source patch field to the target with no left
330  // over values specified. If the interpolation weight sum is less
331  // than one for a face then they will be normalised. If the
332  // interpolation weight sum is zero for a face then that face's
333  // value will be NaN.
334  template<class Type>
335  tmp<Field<Type>> srcToTgt(const Field<Type>& srcFld) const;
336 
337  //- Interpolate a source patch field to the target with left over
338  // values specified. If the interpolation weight sum is less than
339  // one for a face then the average will include the left over
340  // value multiplied by one minus the weight sum.
341  template<class Type>
343  (
344  const Field<Type>& srcFld,
345  const Field<Type>& leftOverTgtFld
346  ) const;
347 
348  //- Interpolate a target patch field to the source with no left
349  // over values specified. As the corresponding srcToTgt.
350  template<class Type>
351  tmp<Field<Type>> tgtToSrc(const Field<Type>& tgtFld) const;
352 
353  //- Interpolate a target patch field to the source with left
354  // over values specified. As the corresponding srcToTgt.
355  template<class Type>
357  (
358  const Field<Type>& tgtFld,
359  const Field<Type>& leftOverSrcFld
360  ) const;
361 
362 
363  // Manipulation
364 
365  //- Update addressing and weights for the given patches
366  void update
367  (
368  const primitiveOldTimePatch& srcPatch,
369  const vectorField& srcPointNormals,
370  const vectorField& srcPointNormals0,
371  const primitiveOldTimePatch& tgtPatch,
372  const transformer& tgtToSrc = NullObjectRef<transformer>()
373  );
374 
375  //- Update addressing and weights for the given patches
376  void update
377  (
378  const primitivePatch& srcPatch,
379  const vectorField& srcPointNormals,
380  const primitivePatch& tgtPatch,
381  const transformer& tgtToSrc = NullObjectRef<transformer>()
382  );
383 
384 
385  // Member Operators
386 
387  //- Disallow default bitwise assignment
388  void operator=(const patchToPatch&) = delete;
389 };
390 
391 
392 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
393 
394 } // End namespace Foam
395 
396 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
397 
398 #include "patchToPatchI.H"
399 
400 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
401 
402 #ifdef NoRepository
403  #include "patchToPatchTemplates.C"
404 #endif
405 
406 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
407 
408 #endif
409 
410 // ************************************************************************* //
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:78
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:91
A bit-packed bool list.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Class containing processor-to-processor mapping information.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
Class to generate coupling geometry between two primitive patches.
Definition: patchToPatch.H:56
virtual tmpNrc< List< DynamicList< scalar > > > srcWeights() const =0
For each source face, the coupled target weights.
static List< remote > distributePatch(const distributionMap &map, const primitiveOldTimePatch &patch, autoPtr< PrimitiveOldTimePatch< faceList, pointField >> &localPatchPtr)
Distribute a patch given its distribution map.
const bool reverse_
Flag to indicate that the two patches are co-directional and.
Definition: patchToPatch.H:63
tmp< Field< Type > > tgtToSrc(const Field< Type > &tgtFld) const
Interpolate a target patch field to the source with no left.
void operator=(const patchToPatch &)=delete
Disallow default bitwise assignment.
label singleProcess_
Index of the processor holding all faces of the patchToPatch, or -1.
Definition: patchToPatch.H:67
List< DynamicList< label > > tgtLocalSrcFaces_
For each target face, the coupled local source faces.
Definition: patchToPatch.H:73
autoPtr< distributionMap > srcMapPtr_
Map from source patch faces to target-local source patch faces.
Definition: patchToPatch.H:76
bool isSingleProcess() const
Is this intersection on a single process?
Definition: patchToPatchI.H:43
virtual treeBoundBox srcBox(const face &srcFace, const pointField &srcPoints, const vectorField &srcPointNormals) const =0
Get the bound box for a source face.
virtual ~patchToPatch()
Destructor.
Definition: patchToPatch.C:779
autoPtr< List< remote > > localSrcProcFacesPtr_
When running in parallel, a map from local source face index to.
Definition: patchToPatch.H:83
autoPtr< List< remote > > localTgtProcFacesPtr_
When running in parallel, a map from local target face index to.
Definition: patchToPatch.H:87
labelListList tgtPatchSendFaces(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch) const
Determine which target faces need to be sent to the source.
virtual void rDistributeTgt(const primitiveOldTimePatch &tgtPatch)
Send data that resulted from an intersection between the source.
Definition: patchToPatch.C:724
tmp< Field< Type > > srcToTgt(const Field< Type > &srcFld) const
Interpolate a source patch field to the target with no left.
virtual labelList finaliseLocal(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Finalise the intersection locally. Trims the local target patch.
Definition: patchToPatch.C:662
virtual tmpNrc< List< DynamicList< scalar > > > tgtWeights() const =0
For each target face, the coupled source weights.
virtual void distributeSrc(const primitiveOldTimePatch &srcPatch)
Distribute the source patch so that any information needed by.
Definition: patchToPatch.C:710
treeBoundBox tgtBox(const primitiveOldTimePatch &tgtPatch, const label tgtFacei) const
Get the bound box for a target face.
Definition: patchToPatch.C:144
static autoPtr< patchToPatch > New(const word &patchToPatchType, const bool reverse)
Select from name.
Definition: patchToPatch.C:786
List< List< remote > > tgtSrcProcFaces() const
For each target face, the coupled source procs and faces.
Definition: patchToPatch.C:847
List< DynamicList< label > > srcLocalTgtFaces_
For each source face, the coupled local target faces.
Definition: patchToPatch.H:70
PackedBoolList srcCoupled() const
Return a list indicating which source faces are coupled.
Definition: patchToPatch.C:810
virtual bool intersectFaces(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const label srcFacei, const label tgtFacei)=0
Intersect two faces.
bool reverse() const
Flag to indicate that the two patches are co-directional and.
Definition: patchToPatchI.H:31
bool findOrIntersectFaces(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const label srcFacei, const label tgtFacei)
Intersect two faces.
Definition: patchToPatch.C:182
List< List< remote > > srcTgtProcFaces() const
For each source face, the coupled target procs and faces.
Definition: patchToPatch.C:833
virtual void initialise(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Initialisation.
Definition: patchToPatch.C:640
TypeName("patchToPatch")
Runtime type information.
void intersectPatches(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch)
Intersect the patches.
Definition: patchToPatch.C:374
label singleProcess() const
Index of the processor holding all faces of the patchToPatch,.
Definition: patchToPatchI.H:37
label intersectPatchQueue(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const bool isSrc, const DynamicList< labelPair > &queue, labelList &faceComplete, DynamicList< labelPair > &otherQueue, const labelList &otherFaceComplete, boolList &otherFaceQueued, boolList &otherFaceVisited)
Intersect a queue of source-target face pairs. Update completion.
Definition: patchToPatch.C:232
void update(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const transformer &tgtToSrc=NullObjectRef< transformer >())
Update addressing and weights for the given patches.
Definition: patchToPatch.C:861
PackedBoolList tgtCoupled() const
Return a list indicating which target faces are coupled.
Definition: patchToPatch.C:821
patchToPatch(const bool reverse)
Construct from components.
Definition: patchToPatch.C:764
autoPtr< distributionMap > tgtMapPtr_
Map from target patch faces to source-local target patch faces.
Definition: patchToPatch.H:79
virtual label finalise(const primitiveOldTimePatch &srcPatch, const vectorField &srcPointNormals, const vectorField &srcPointNormals0, const primitiveOldTimePatch &tgtPatch, const transformer &tgtToSrc)
Finalising.
Definition: patchToPatch.C:737
declareRunTimeSelectionTable(autoPtr, patchToPatch, bool,(const bool reverse),(reverse))
Declare runtime constructor selection table.
A class for managing temporary objects without reference counting.
Definition: tmpNrc.H:53
A class for managing temporary objects.
Definition: tmp.H:55
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space.
Definition: transformer.H:84
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:90
A class for handling words, derived from string.
Definition: word.H:62
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
Macros to ease declaration of run-time selection tables.