inversePointDistanceDiffusivity.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) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2018 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
31#include "HashSet.H"
32#include "pointEdgePoint.H"
33#include "PointEdgeWave.H"
34
35// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36
37namespace Foam
38{
40
42 (
46 );
47}
48
49
50// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51
53(
54 const fvMesh& mesh,
55 Istream& mdData
56)
57:
58 uniformDiffusivity(mesh, mdData),
59 patchNames_(mdData)
60{
61 correct();
62}
63
64
65// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
66
68{
69 const polyBoundaryMesh& bdry = mesh().boundaryMesh();
70
71 labelHashSet patchSet(bdry.patchSet(patchNames_));
72
73 label nPatchEdges = 0;
74
75 for (const label patchi : patchSet)
76 {
77 nPatchEdges += bdry[patchi].nEdges();
78 }
79
80 // Distance to wall on points and edges.
81 List<pointEdgePoint> pointWallDist(mesh().nPoints());
82 List<pointEdgePoint> edgeWallDist(mesh().nEdges());
83
84 int dummyTrackData = 0;
85
86
87 {
88 // Seeds
89 List<pointEdgePoint> seedInfo(nPatchEdges);
90 labelList seedPoints(nPatchEdges);
91
92 nPatchEdges = 0;
93
94 for (const label patchi : patchSet)
95 {
96 const polyPatch& patch = bdry[patchi];
97
98 const labelList& meshPoints = patch.meshPoints();
99
100 for (const label pointi : meshPoints)
101 {
102 if (!pointWallDist[pointi].valid(dummyTrackData))
103 {
104 // Not yet seeded
105 seedInfo[nPatchEdges] = pointEdgePoint
106 (
107 mesh().points()[pointi],
108 0.0
109 );
110 seedPoints[nPatchEdges] = pointi;
111 pointWallDist[pointi] = seedInfo[nPatchEdges];
112
113 nPatchEdges++;
114 }
115 }
116 }
117 seedInfo.setSize(nPatchEdges);
118 seedPoints.setSize(nPatchEdges);
119
120 // Do calculations
122 (
123 mesh(),
124 seedPoints,
125 seedInfo,
126
127 pointWallDist,
128 edgeWallDist,
129 mesh().globalData().nTotalPoints(),// max iterations
130 dummyTrackData
131 );
132 }
133
134
135 for (label facei=0; facei<mesh().nInternalFaces(); ++facei)
136 {
137 const face& f = mesh().faces()[facei];
138
139 scalar dist = 0;
140
141 forAll(f, fp)
142 {
143 dist += sqrt(pointWallDist[f[fp]].distSqr());
144 }
145 dist /= f.size();
146
147 faceDiffusivity_[facei] = 1.0/dist;
148 }
149
150
151 surfaceScalarField::Boundary& faceDiffusivityBf =
152 faceDiffusivity_.boundaryFieldRef();
153
154 forAll(faceDiffusivityBf, patchi)
155 {
156 fvsPatchScalarField& bfld = faceDiffusivityBf[patchi];
157
158 if (patchSet.found(patchi))
159 {
160 const labelUList& faceCells = bfld.patch().faceCells();
161
162 forAll(bfld, i)
163 {
164 const cell& ownFaces = mesh().cells()[faceCells[i]];
165
166 labelHashSet cPoints(4*ownFaces.size());
167
168 scalar dist = 0;
169
170 forAll(ownFaces, ownFacei)
171 {
172 const face& f = mesh().faces()[ownFaces[ownFacei]];
173
174 forAll(f, fp)
175 {
176 if (cPoints.insert(f[fp]))
177 {
178 dist += sqrt(pointWallDist[f[fp]].distSqr());
179 }
180 }
181 }
182 dist /= cPoints.size();
183
184 bfld[i] = 1.0/dist;
185 }
186 }
187 else
188 {
189 const label start = bfld.patch().start();
190
191 forAll(bfld, i)
192 {
193 const face& f = mesh().faces()[start+i];
194
195 scalar dist = 0;
196
197 forAll(f, fp)
198 {
199 dist += sqrt(pointWallDist[f[fp]].distSqr());
200 }
201 dist /= f.size();
202
203 bfld[i] = 1.0/dist;
204 }
205 }
206 }
207}
208
209
210// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
bool insert(const Key &key)
Insert a new entry, not overwriting existing entries.
Definition: HashSet.H:191
bool found(const Key &key) const
Return true if hashed entry is found in table.
Definition: HashTableI.H:100
label size() const noexcept
The number of elements in table.
Definition: HashTableI.H:52
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:64
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
void setSize(const label n)
Alias for resize()
Definition: List.H:218
Wave propagation of information through grid. Every iteration information goes through one layer of e...
Definition: PointEdgeWave.H:91
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:57
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:59
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
virtual label start() const
Return start label of this patch in the polyMesh face list.
Definition: fvPatch.H:179
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:113
const fvPatch & patch() const
Return patch.
Inverse distance to the given patches motion diffusivity.
virtual void correct()
Correct the motion diffusivity.
Abstract base class for cell-centre mesh motion diffusivity.
Holds information regarding nearest wall point. Used in PointEdgeWave. (so not standard FaceCellWave)...
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
labelHashSet patchSet(const UList< wordRe > &patchNames, const bool warnNotFound=true, const bool useGroups=true) const
Return the set of patch IDs corresponding to the given names.
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1108
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
label nInternalFaces() const noexcept
Number of internal faces.
const cellList & cells() const
Uniform uniform finite volume mesh motion diffusivity.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
const pointField & points
label nPoints
Namespace for OpenFOAM.
dimensionedScalar sqrt(const dimensionedScalar &ds)
labelList f(nPoints)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333