pointPatchDist.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) 2013-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 "pointPatchDist.H"
27 #include "externalPointEdgePoint.H"
28 #include "pointMesh.H"
29 #include "PointEdgeWave.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 Foam::pointPatchDist::pointPatchDist
34 (
35  const pointMesh& pMesh,
36  const labelHashSet& patchIDs,
37  const pointField& points
38 )
39 :
41  (
42  IOobject
43  (
44  "pointDistance",
45  pMesh.db().time().timeName(),
46  pMesh.db()
47  ),
48  pMesh,
49  dimensionedScalar("y", dimLength, great)
50  ),
51  points_(points),
52  patchIDs_(patchIDs),
53  nUnset_(0)
54 {
55  correct();
56 }
57 
58 
59 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
60 
62 {}
63 
64 
65 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66 
68 {
69  const pointBoundaryMesh& pbm = mesh().boundary();
70 
71  label nPoints = 0;
72 
73  forAllConstIter(labelHashSet, patchIDs_, iter)
74  {
75  label patchi = iter.key();
76  nPoints += pbm[patchi].meshPoints().size();
77  }
78 
80 
81  // Set initial changed points to all the patch points(if patch present)
82  List<externalPointEdgePoint> wallInfo(nPoints);
83  labelList wallPoints(nPoints);
84  nPoints = 0;
85 
86  forAllConstIter(labelHashSet, patchIDs_, iter)
87  {
88  label patchi = iter.key();
89  // Retrieve the patch now we have its index in patches.
90 
91  const labelList& mp = pbm[patchi].meshPoints();
92 
93  forAll(mp, ppI)
94  {
95  label meshPointi = mp[ppI];
96  wallPoints[nPoints] = meshPointi;
98  (
99  td.points_[meshPointi],
100  0.0
101  );
102  nPoints++;
103  }
104  }
105 
106  // Current info on points
107  List<externalPointEdgePoint> allPointInfo(mesh()().nPoints());
108 
109  // Current info on edges
110  List<externalPointEdgePoint> allEdgeInfo(mesh()().nEdges());
111 
113  <
116  > wallCalc
117  (
118  mesh()(),
119  wallPoints,
120  wallInfo,
121 
122  allPointInfo,
123  allEdgeInfo,
124  mesh().globalData().nTotalPoints(), // max iterations
125  td
126  );
127 
128  pointScalarField& psf = *this;
129 
130 
131  forAll(allPointInfo, pointi)
132  {
133  if (allPointInfo[pointi].valid(td))
134  {
135  psf[pointi] = Foam::sqrt(allPointInfo[pointi].distSqr());
136  }
137  else
138  {
139  nUnset_++;
140  }
141  }
142 }
143 
144 
145 // ************************************************************************* //
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
Definition: pointFields.H:49
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
Foam::pointBoundaryMesh.
virtual ~pointPatchDist()
Destructor.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
dimensionedScalar sqrt(const dimensionedScalar &ds)
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:626
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
Holds information regarding nearest wall point. Used in PointEdgeWave. (so not standard FaceCellWave)...
label nPoints
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
const Mesh & mesh() const
Return mesh.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-6/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:72
const Time & time() const
Return time.
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
label patchi
Class used to pass data into container.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:361
Wave propagation of information through grid. Every iteration information goes through one layer of e...
Definition: PointEdgeWave.H:85
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
const dimensionedScalar mp
Proton mass.
void correct()
Correct for mesh geom/topo changes.