faPatch.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) 2016-2017 Wikki Ltd
9 Copyright (C) 2019-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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 "faPatch.H"
31#include "faBoundaryMesh.H"
32#include "faMesh.H"
33#include "areaFields.H"
34#include "edgeFields.H"
35#include "edgeHashes.H"
36#include "polyMesh.H"
37#include "polyPatch.H"
38//#include "pointPatchField.H"
39#include "demandDrivenData.H"
40
41// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
42
43namespace Foam
44{
48}
49
50
51// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
52
53bool Foam::faPatch::constraintType(const word& patchType)
54{
55 // Reasonable to expect any faPatch constraint has an identically
56 // named polyPatch/pointPatch equivalent
57
58 return polyPatch::constraintType(patchType);
59}
60
61
63{
64 const auto& cnstrTable = *dictionaryConstructorTablePtr_;
65
66 wordList cTypes(cnstrTable.size());
67
68 label i = 0;
69
70 forAllConstIters(cnstrTable, iter)
71 {
72 if (constraintType(iter.key()))
73 {
74 cTypes[i++] = iter.key();
75 }
76 }
77
78 cTypes.resize(i);
79
80 return cTypes;
81}
82
83
84// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
85
86void Foam::faPatch::clearOut()
87{
88 deleteDemandDrivenData(edgeFacesPtr_);
89 deleteDemandDrivenData(pointLabelsPtr_);
90 deleteDemandDrivenData(pointEdgesPtr_);
91}
92
93
94// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
95
97(
98 const word& name,
99 const labelUList& edgeLabels,
100 const label index,
101 const faBoundaryMesh& bm,
102 const label nbrPolyPatchi,
103 const word& patchType
104)
105:
106 patchIdentifier(name, index),
107 labelList(edgeLabels),
108 nbrPolyPatchId_(nbrPolyPatchi),
109 boundaryMesh_(bm),
110 edgeFacesPtr_(nullptr),
111 pointLabelsPtr_(nullptr),
112 pointEdgesPtr_(nullptr)
113{
114 if (!patchType.empty() && constraintType(patchType))
115 {
116 inGroups().appendUniq(patchType);
117 }
118}
119
120
122(
123 const word& name,
124 const dictionary& dict,
125 const label index,
126 const faBoundaryMesh& bm,
127 const word& patchType
128)
129:
130 patchIdentifier(name, dict, index),
131 labelList(dict.get<labelList>("edgeLabels")),
132 nbrPolyPatchId_(dict.get<label>("ngbPolyPatchIndex")),
133 boundaryMesh_(bm),
134 edgeFacesPtr_(nullptr),
135 pointLabelsPtr_(nullptr),
136 pointEdgesPtr_(nullptr)
137{
138 if (!patchType.empty() && constraintType(patchType))
139 {
140 inGroups().appendUniq(patchType);
141 }
142}
143
144
146(
147 const faPatch& p,
148 const faBoundaryMesh& bm,
149 const label index,
150 const labelUList& edgeLabels,
151 const label nbrPolyPatchi
152)
153:
154 patchIdentifier(p, index),
155 labelList(edgeLabels),
156 nbrPolyPatchId_(p.nbrPolyPatchId_),
157 boundaryMesh_(bm),
158 edgeFacesPtr_(nullptr),
159 pointLabelsPtr_(nullptr),
160 pointEdgesPtr_(nullptr)
161{}
162
163
165(
166 const faPatch& p,
167 const faBoundaryMesh& bm
168)
169:
170 faPatch
171 (
172 p,
173 bm,
174 p.index(),
175 p.edgeLabels(),
176 p.nbrPolyPatchId_
177 )
178{}
179
180
181// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
182
184{
185 clearOut();
186}
187
188
189// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
190
192{
193 return boundaryMesh_;
194}
195
196
197Foam::label Foam::faPatch::start() const
198{
199 return boundaryMesh().mesh().patchStarts()[index()];
200}
201
202
204{
205 const auto& connections = boundaryMesh().mesh().boundaryConnections();
206 const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
207
209
210 // Like an IndirectList but removing the nInternalEdges offset
211 label count = 0;
212 for (const label patchEdgei : this->edgeLabels())
213 {
214 const label bndEdgei = (patchEdgei - nInternalEdges);
215 output[count] = connections[bndEdgei];
216 ++count;
217 }
218
219 return output;
220}
221
222
224{
225 const auto& connections = boundaryMesh().mesh().boundaryConnections();
226 const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
227
228 labelHashSet procsUsed(2*Pstream::nProcs());
229
230 for (const label patchEdgei : this->edgeLabels())
231 {
232 const label bndEdgei = (patchEdgei - nInternalEdges);
233 const label proci = connections[bndEdgei].first();
234 procsUsed.insert(proci);
235 }
236 procsUsed.erase(-1); // placeholder value
237 procsUsed.erase(Pstream::myProcNo());
238
239 return procsUsed.sortedToc();
240}
241
242
244{
245 const auto& connections = boundaryMesh().mesh().boundaryConnections();
246 const label nInternalEdges = boundaryMesh().mesh().nInternalEdges();
247
248 Map<label> procCount(2*Pstream::nProcs());
249
250 for (const label patchEdgei : this->edgeLabels())
251 {
252 const label bndEdgei = (patchEdgei - nInternalEdges);
253 const label proci = connections[bndEdgei].first();
254 ++procCount(proci);
255 }
256 procCount.erase(-1); // placeholder value
257 procCount.erase(Pstream::myProcNo());
258
259 // Flatten as list
260 List<labelPair> output(procCount.size());
261 label count = 0;
262 for (const label proci : procCount.sortedToc())
263 {
264 output[count].first() = proci;
265 output[count].second() = procCount[proci]; // size
266 ++count;
267 }
268
269 return output;
270}
271
272
274{
275 if (!pointLabelsPtr_)
276 {
277 calcPointLabels();
278 }
279
280 return *pointLabelsPtr_;
281}
282
283
285{
286 if (!pointEdgesPtr_)
287 {
288 calcPointEdges();
289 }
290
291 return *pointEdgesPtr_;
292}
293
294
296{
297 const edgeList::subList edges = patchSlice(boundaryMesh().mesh().edges());
298
299 // Walk boundary edges.
300 // The edge orientation corresponds to the face orientation
301 // (outwards normal).
302
303 // Note: could combine this with calcPointEdges for more efficiency
304
305 // Map<label> markedPoints(4*edges.size());
306 labelHashSet markedPoints(4*edges.size());
307 DynamicList<label> dynEdgePoints(2*edges.size());
308
309 for (const edge& e : edges)
310 {
311 // if (markedPoints.insert(e.first(), markedPoints.size()))
312 if (markedPoints.insert(e.first()))
313 {
314 dynEdgePoints.append(e.first());
315 }
316 // if (markedPoints.insert(e.second(), markedPoints.size()))
317 if (markedPoints.insert(e.second()))
318 {
319 dynEdgePoints.append(e.second());
320 }
321 }
322
323 // Transfer to plain list (reuse storage)
324 pointLabelsPtr_ = new labelList(std::move(dynEdgePoints));
349}
350
351
353{
354 const edgeList::subList edges = patchSlice(boundaryMesh().mesh().edges());
355
356 const labelList& edgePoints = pointLabels();
357
358 // Cannot use invertManyToMany - we have non-local edge numbering
359
360 // Intermediate storage for pointEdges.
361 // Points on the boundary will normally connect 1 or 2 edges only.
362 List<DynamicList<label,2>> dynPointEdges(edgePoints.size());
363
364 forAll(edges, edgei)
365 {
366 const edge& e = edges[edgei];
367
368 dynPointEdges[edgePoints.find(e.first())].append(edgei);
369 dynPointEdges[edgePoints.find(e.second())].append(edgei);
370 }
371
372 // Flatten to regular list
373 pointEdgesPtr_ = new labelListList(edgePoints.size());
374 auto& pEdges = *pointEdgesPtr_;
375
376 forAll(pEdges, pointi)
377 {
378 pEdges[pointi] = std::move(dynPointEdges[pointi]);
379 }
380}
381
382
384{
385 if (nbrPolyPatchId_ < 0)
386 {
387 return tmp<vectorField>::New();
388 }
389
390 return boundaryMesh().mesh().haloFaceNormals(this->index());
391}
392
393
395{
396 if (nbrPolyPatchId_ < 0)
397 {
398 return tmp<vectorField>::New();
399 }
400
401 // Unit normals for the neighbour patch faces
403 (
404 boundaryMesh().mesh().haloFaceNormals(this->index())
405 );
406
407 const labelListList& pntEdges = pointEdges();
408
409 auto tpointNorm = tmp<vectorField>::New(pntEdges.size());
410 auto& pointNorm = tpointNorm.ref();
411
412 forAll(pointNorm, pointi)
413 {
414 vector& n = pointNorm[pointi];
415 n = Zero;
416
417 for (const label bndEdgei : pntEdges[pointi])
418 {
419 n += faceNormals[bndEdgei];
420 }
421
422 n.normalise();
423 }
424
425 return tpointNorm;
426}
427
428
430{
431 if (!edgeFacesPtr_)
432 {
433 edgeFacesPtr_ = new labelList::subList
434 (
435 patchSlice(boundaryMesh().mesh().edgeOwner())
436 );
437 }
438
439 return *edgeFacesPtr_;
440}
441
442
444{
445 return boundaryMesh().mesh().edgeCentres().boundaryField()[index()];
446}
447
448
450{
451 return boundaryMesh().mesh().Le().boundaryField()[index()];
452}
453
454
456{
457 return boundaryMesh().mesh().magLe().boundaryField()[index()];
458}
459
460
462{
463 auto tedgeNorm = tmp<vectorField>::New(edgeLengths());
464
465 for (vector& n : tedgeNorm.ref())
466 {
467 n.normalise();
468 }
469
470 return tedgeNorm;
471}
472
473
475{
476 auto tfc = tmp<vectorField>::New(size());
477 auto& fc = tfc.ref();
478
479 // get reference to global face centres
480 const vectorField& gfc =
481 boundaryMesh().mesh().areaCentres().internalField();
482
483 const labelUList& faceLabels = edgeFaces();
484
485 forAll(faceLabels, edgeI)
486 {
487 fc[edgeI] = gfc[faceLabels[edgeI]];
488 }
489
490 return tfc;
491}
492
493
495{
496 return edgeNormals()*(edgeNormals() & (edgeCentres() - edgeFaceCentres()));
497}
498
499
501{
502 dc = scalar(1)/(edgeNormals() & delta());
503}
504
505
507{
508 vectorField unitDelta(delta()/mag(delta()));
509 vectorField edgeNormMag(edgeNormals()/mag(edgeNormals()));
510 scalarField dn(edgeNormals() & delta());
511
512 k = edgeNormMag - (scalar(1)/(unitDelta & edgeNormMag))*unitDelta;
513}
514
515
517{
518 return boundaryMesh().mesh().deltaCoeffs().boundaryField()[index()];
519}
520
521
523{
524 w = scalar(1);
525}
526
527
529{
530 return boundaryMesh().mesh().weights().boundaryField()[index()];
531}
532
533
535{}
536
537
539{
540 clearOut();
541 static_cast<labelList&>(*this) = newEdges;
542}
543
544
546{
547 clearOut();
548 static_cast<labelList&>(*this) = std::move(newEdges);
549}
550
551
553{
554 os.writeEntry("type", type());
555
557
558 os.writeEntry("ngbPolyPatchIndex", nbrPolyPatchId_);
559 static_cast<const labelList&>(*this).writeEntry("edgeLabels", os);
560}
561
562
563// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
564
566{
567 p.write(os);
569 return os;
570}
571
572
573// ************************************************************************* //
label k
scalar delta
label n
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
void append(const T &val)
Copy append an element to the end of this list.
Definition: DynamicListI.H:503
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
List< Key > sortedToc() const
The table of contents (the keys) in sorted order.
Definition: HashTable.C:137
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
bool erase(const iterator &iter)
Erase an entry specified by given iterator.
Definition: HashTable.C:440
virtual bool check(const char *operation) const
Check IOstream status for given operation.
Definition: IOstream.C:58
SubList< label > subList
Declare type of subList.
Definition: List.H:111
void append(const T &val)
Append an element at the end of the list.
Definition: ListI.H:175
label appendUniq(const T &val)
Append an element if not already in the list.
Definition: ListI.H:232
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:239
label nInternalEdges() const
Number of internal edges.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
label nProcs() const noexcept
Number of ranks associated with PstreamBuffers.
A List obtained as a section of another List.
Definition: SubList.H:70
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
T & first()
Return the first element of the list.
Definition: UListI.H:202
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
Definition: UList.C:212
void writeEntry(Ostream &os) const
Write the UList with its compound type.
Definition: UListIO.C:38
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface....
Definition: boundaryMesh.H:63
const bMesh & mesh() const
Definition: boundaryMesh.H:206
virtual const word & constraintType() const
Return the constraint type this pointPatch implements.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
Finite area boundary mesh.
Finite area patch class. Used for 2-D non-Euclidian finite area method.
Definition: faPatch.H:78
static wordList constraintTypes()
Return a list of all the constraint patch types.
Definition: faPatch.C:62
const labelListList & pointEdges() const
Return patch point-edge addressing.
Definition: faPatch.C:284
tmp< vectorField > ngbPolyPatchPointNormals() const
Return normals of neighbour polyPatch joined points.
Definition: faPatch.C:394
const labelList & pointLabels() const
Return patch point labels.
Definition: faPatch.C:273
virtual ~faPatch()
Destructor.
Definition: faPatch.C:183
void resetEdges(const labelUList &newEdges)
Reset the list of edges (use with caution)
Definition: faPatch.C:538
virtual void makeWeights(scalarField &) const
Make patch weighting factors.
Definition: faPatch.C:522
virtual void makeDeltaCoeffs(scalarField &) const
Make patch edge - neighbour face distances.
Definition: faPatch.C:500
List< labelPair > boundaryConnections() const
Definition: faPatch.C:203
const scalarField & magEdgeLengths() const
Return edge length magnitudes.
Definition: faPatch.C:455
tmp< vectorField > ngbPolyPatchFaceNormals() const
Return normals of neighbour polyPatch faces.
Definition: faPatch.C:383
const vectorField & edgeLengths() const
Return edge length vectors.
Definition: faPatch.C:449
virtual tmp< vectorField > delta() const
Return cell-centre to face-centre vector.
Definition: faPatch.C:494
tmp< vectorField > edgeFaceCentres() const
Return neighbour face centres.
Definition: faPatch.C:474
void calcPointEdges() const
Calculate patch point-edge addressing.
Definition: faPatch.C:352
const faBoundaryMesh & boundaryMesh() const noexcept
Return boundaryMesh reference.
Definition: faPatch.C:191
static bool constraintType(const word &pt)
Return true if the given type is a constraint type.
Definition: faPatch.C:53
const scalarField & weights() const
Return patch weighting factors.
Definition: faPatch.C:528
const labelUList & edgeFaces() const
Return edge-face addressing.
Definition: faPatch.C:429
const vectorField & edgeCentres() const
Return edge centres.
Definition: faPatch.C:443
labelList boundaryProcs() const
Definition: faPatch.C:223
tmp< vectorField > edgeNormals() const
Return edge normals.
Definition: faPatch.C:461
void calcPointLabels() const
Calculate patch point labels.
Definition: faPatch.C:295
List< labelPair > boundaryProcSizes() const
Definition: faPatch.C:243
const scalarField & deltaCoeffs() const
Return patch edge - neighbour face distances.
Definition: faPatch.C:516
void makeCorrectionVectors(vectorField &) const
Definition: faPatch.C:506
label start() const
Patch start in edge list.
Definition: faPatch.C:197
virtual bool write()
Write the output fields.
void movePoints()
Update for new mesh geometry.
Identifies a patch by name and index, with optional physical type and group information.
const wordList & inGroups() const noexcept
The (optional) groups that the patch belongs to.
int myProcNo() const noexcept
Return processor number.
virtual bool write(const bool valid=true) const
Write using setting from DB.
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
volScalarField & p
dynamicFvMesh & mesh
Template functions to aid in the implementation of demand driven data.
OBJstream os(runTime.globalPath()/outputName)
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
surfaceVectorField faceNormals(mesh.Sf()/mesh.magSf())
#define FUNCTION_NAME
Namespace for OpenFOAM.
List< label > labelList
A List of labels.
Definition: List.H:66
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & operator<<(Ostream &, const boundaryPatch &p)
Write boundaryPatch as dictionary entries (without surrounding braces)
Definition: boundaryPatch.C:83
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
List< labelList > labelListList
A List of labelList.
Definition: labelList.H:56
static Ostream & output(Ostream &os, const IntRange< T > &range)
Definition: IntRanges.C:66
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
const direction noexcept
Definition: Scalar.H:223
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
void deleteDemandDrivenData(DataPtr &dataPtr)
labelList pointLabels(nPoints, -1)
#define defineRunTimeSelectionTable(baseType, argNames)
Define run-time selection table.
dictionary dict
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Definition: stdFoam.H:278