patchPatchDist.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) 2011-2018 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 
26 #include "patchPatchDist.H"
27 #include "PatchEdgeFaceWave.H"
28 #include "syncTools.H"
29 #include "polyMesh.H"
30 #include "patchEdgeFaceInfo.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
35 (
36  const polyPatch& patch,
37  const labelHashSet& nbrPatchIDs
38 )
39 :
40  patch_(patch),
41  nbrPatchIDs_(nbrPatchIDs),
42  nUnset_(0)
43 {
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
49 
51 {}
52 
53 
54 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
55 
57 {
58  // Mark all edge connected to a nbrPatch.
59  label nBnd = 0;
60  forAllConstIter(labelHashSet, nbrPatchIDs_, iter)
61  {
62  label nbrPatchi = iter.key();
63  const polyPatch& nbrPatch = patch_.boundaryMesh()[nbrPatchi];
64  nBnd += nbrPatch.nEdges()-nbrPatch.nInternalEdges();
65  }
66 
67  // Mark all edges. Note: should use HashSet but have no syncTools
68  // functionality for these.
69  EdgeMap<label> nbrEdges(2*nBnd);
70 
71  forAllConstIter(labelHashSet, nbrPatchIDs_, iter)
72  {
73  label nbrPatchi = iter.key();
74  const polyPatch& nbrPatch = patch_.boundaryMesh()[nbrPatchi];
75  const labelList& nbrMp = nbrPatch.meshPoints();
76 
77  for
78  (
79  label edgeI = nbrPatch.nInternalEdges();
80  edgeI < nbrPatch.nEdges();
81  edgeI++
82  )
83  {
84  const edge& e = nbrPatch.edges()[edgeI];
85  const edge meshE = edge(nbrMp[e[0]], nbrMp[e[1]]);
86  nbrEdges.insert(meshE, nbrPatchi);
87  }
88  }
89 
90 
91  // Make sure these boundary edges are marked everywhere.
93  (
94  patch_.boundaryMesh().mesh(),
95  nbrEdges,
97  );
98 
99 
100  // Data on all edges and faces
101  List<patchEdgeFaceInfo> allEdgeInfo(patch_.nEdges());
102  List<patchEdgeFaceInfo> allFaceInfo(patch_.size());
103 
104  // Initial seed
105  label nBndEdges = patch_.nEdges() - patch_.nInternalEdges();
106  DynamicList<label> initialEdges(2*nBndEdges);
107  DynamicList<patchEdgeFaceInfo> initialEdgesInfo(2*nBndEdges);
108 
109 
110  // Seed all my edges that are also nbrEdges
111 
112  const labelList& mp = patch_.meshPoints();
113 
114  for
115  (
116  label edgeI = patch_.nInternalEdges();
117  edgeI < patch_.nEdges();
118  edgeI++
119  )
120  {
121  const edge& e = patch_.edges()[edgeI];
122  const edge meshE = edge(mp[e[0]], mp[e[1]]);
123  EdgeMap<label>::const_iterator edgeFnd = nbrEdges.find(meshE);
124  if (edgeFnd != nbrEdges.end())
125  {
126  initialEdges.append(edgeI);
127  initialEdgesInfo.append
128  (
130  (
131  e.centre(patch_.localPoints()),
132  0.0
133  )
134  );
135  }
136  }
137 
138 
139  // Walk
141  <
144  > calc
145  (
146  patch_.boundaryMesh().mesh(),
147  patch_,
148  initialEdges,
149  initialEdgesInfo,
150  allEdgeInfo,
151  allFaceInfo,
152  returnReduce(patch_.nEdges(), sumOp<label>())
153  );
154 
155 
156  // Extract into *this
157  setSize(patch_.size());
158  nUnset_ = 0;
159  forAll(allFaceInfo, facei)
160  {
161  if (allFaceInfo[facei].valid(calc.data()))
162  {
163  operator[](facei) = Foam::sqrt(allFaceInfo[facei].distSqr());
164  }
165  else
166  {
167  nUnset_++;
168  }
169  }
170 }
171 
172 
173 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:112
Type & operator[](const label)
Return element of UList.
Definition: UListI.H:167
Wave propagation of information along patch. Every iteration information goes through one layer of fa...
const polyBoundaryMesh & boundaryMesh() const
Return boundaryMesh reference.
Definition: polyPatch.C:291
dimensionedScalar sqrt(const dimensionedScalar &ds)
label nInternalEdges() const
Number of internal edges.
const Field< PointType > & localPoints() const
Return pointField of points in patch.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
bool insert(const Key &, const T &newElmt)
Insert a new hashedEntry.
Definition: HashTableI.H:80
iterator find(const Key &)
Find and return an iterator set at the hashedEntry.
Definition: HashTable.C:142
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-9/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
virtual ~patchPatchDist()
Destructor.
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:296
const polyMesh & mesh() const
Return the mesh reference.
const dimensionedScalar mp
Proton mass.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
label nEdges() const
Return number of edges in patch.
static void syncEdgeMap(const polyMesh &, EdgeMap< T > &edgeValues, const CombineOp &cop, const TransformOp &top)
Synchronise values on selected edges.
PrimitivePatch< SubList< face >, const pointField & > primitivePatch
Addressing for a faceList slice.
patchPatchDist(const polyPatch &pp, const labelHashSet &nbrPatchIDs)
Construct from patch and neighbour patches.
void setSize(const label)
Reset size of List.
Definition: List.C:281
virtual void correct()
Correct for mesh geom/topo changes.
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66