polyMeshTools.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) 2012-2016 OpenFOAM Foundation
9 Copyright (C) 2021-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 "polyMeshTools.H"
30#include "syncTools.H"
31#include "pyramidPointFaceRef.H"
32#include "primitiveMeshTools.H"
33
34// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
35
37(
38 const polyMesh& mesh,
39 const vectorField& areas,
40 const vectorField& cc
41)
42{
43 const labelList& own = mesh.faceOwner();
44 const labelList& nei = mesh.faceNeighbour();
45 const polyBoundaryMesh& pbm = mesh.boundaryMesh();
46
47 auto tortho = tmp<scalarField>::New(mesh.nFaces(), scalar(1));
48 auto& ortho = tortho.ref();
49
50 // Internal faces
51 forAll(nei, facei)
52 {
54 (
55 cc[own[facei]],
56 cc[nei[facei]],
57 areas[facei]
58 );
59 }
60
61
62 // Coupled faces
63
64 pointField neighbourCc;
66
67 for (const polyPatch& pp : pbm)
68 {
69 if (pp.coupled())
70 {
71 forAll(pp, i)
72 {
73 label facei = pp.start() + i;
74 label bFacei = facei - mesh.nInternalFaces();
75
77 (
78 cc[own[facei]],
79 neighbourCc[bFacei],
80 areas[facei]
81 );
82 }
83 }
84 }
85
86 return tortho;
87}
88
89
91(
92 const polyMesh& mesh,
93 const pointField& p,
94 const vectorField& fCtrs,
95 const vectorField& fAreas,
96 const vectorField& cellCtrs
97)
98{
99 const labelList& own = mesh.faceOwner();
100 const labelList& nei = mesh.faceNeighbour();
101 const faceList& faces = mesh.faces();
102 const polyBoundaryMesh& pbm = mesh.boundaryMesh();
103
104 auto tskew = tmp<scalarField>::New(mesh.nFaces());
105 auto& skew = tskew.ref();
106
107 forAll(nei, facei)
108 {
110 (
111 faces,
112 p,
113 fCtrs,
114 fAreas,
115
116 facei,
117 cellCtrs[own[facei]],
118 cellCtrs[nei[facei]]
119 );
120 }
121
122
123 // Boundary faces: consider them to have only skewness error.
124 // (i.e. treat as if mirror cell on other side)
125
126 pointField neighbourCc;
127 syncTools::swapBoundaryCellPositions(mesh, cellCtrs, neighbourCc);
128
129 for (const polyPatch& pp : pbm)
130 {
131 if (pp.coupled())
132 {
133 forAll(pp, i)
134 {
135 label facei = pp.start() + i;
136 label bFacei = facei - mesh.nInternalFaces();
137
139 (
140 faces,
141 p,
142 fCtrs,
143 fAreas,
144
145 facei,
146 cellCtrs[own[facei]],
147 neighbourCc[bFacei]
148 );
149 }
150 }
151 else
152 {
153 forAll(pp, i)
154 {
155 label facei = pp.start() + i;
156
158 (
159 faces,
160 p,
161 fCtrs,
162 fAreas,
163
164 facei,
165 cellCtrs[own[facei]]
166 );
167 }
168 }
169 }
170
171 return tskew;
172}
173
174
176(
177 const polyMesh& mesh,
178 const vectorField& fCtrs,
179 const vectorField& fAreas,
180 const vectorField& cellCtrs
181)
182{
183 const labelList& own = mesh.faceOwner();
184 const labelList& nei = mesh.faceNeighbour();
185 const polyBoundaryMesh& pbm = mesh.boundaryMesh();
186
187 auto tweight = tmp<scalarField>::New(mesh.nFaces(), scalar(1));
188 auto& weight = tweight.ref();
189
190 // Internal faces
191 forAll(nei, facei)
192 {
193 const point& fc = fCtrs[facei];
194 const vector& fa = fAreas[facei];
195
196 scalar dOwn = mag(fa & (fc-cellCtrs[own[facei]]));
197 scalar dNei = mag(fa & (cellCtrs[nei[facei]]-fc));
198
199 weight[facei] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
200 }
201
202
203 // Coupled faces
204
205 pointField neiCc;
207
208 for (const polyPatch& pp : pbm)
209 {
210 if (pp.coupled())
211 {
212 forAll(pp, i)
213 {
214 label facei = pp.start() + i;
215 label bFacei = facei - mesh.nInternalFaces();
216
217 const point& fc = fCtrs[facei];
218 const vector& fa = fAreas[facei];
219
220 scalar dOwn = mag(fa & (fc-cellCtrs[own[facei]]));
221 scalar dNei = mag(fa & (neiCc[bFacei]-fc));
222
223 weight[facei] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
224 }
225 }
226 }
227
228 return tweight;
229}
230
231
233(
234 const polyMesh& mesh,
235 const scalarField& vol
236)
237{
238 const labelList& own = mesh.faceOwner();
239 const labelList& nei = mesh.faceNeighbour();
240 const polyBoundaryMesh& pbm = mesh.boundaryMesh();
241
242 auto tratio = tmp<scalarField>::New(mesh.nFaces(), scalar(1));
243 auto& ratio = tratio.ref();
244
245 // Internal faces
246 forAll(nei, facei)
247 {
248 scalar volOwn = vol[own[facei]];
249 scalar volNei = vol[nei[facei]];
250
251 ratio[facei] = min(volOwn,volNei)/(max(volOwn, volNei)+VSMALL);
252 }
253
254
255 // Coupled faces
256
257 scalarField neiVol;
259
260 for (const polyPatch& pp : pbm)
261 {
262 if (pp.coupled())
263 {
264 forAll(pp, i)
265 {
266 label facei = pp.start() + i;
267 label bFacei = facei - mesh.nInternalFaces();
268
269 scalar volOwn = vol[own[facei]];
270 scalar volNei = neiVol[bFacei];
271
272 ratio[facei] = min(volOwn,volNei)/(max(volOwn, volNei)+VSMALL);
273 }
274 }
275 }
276
277 return tratio;
278}
279
280
282(
283 const polyMesh::readUpdateState& state0,
284 const polyMesh::readUpdateState& state1
285)
286{
287 if
288 (
289 (
290 state0 == polyMesh::UNCHANGED
291 && state1 != polyMesh::UNCHANGED
292 )
293 || (
294 state0 == polyMesh::POINTS_MOVED
295 && (state1 != polyMesh::UNCHANGED && state1 != polyMesh::POINTS_MOVED)
296 )
297 || (
298 state0 == polyMesh::TOPO_CHANGE
299 && state1 == polyMesh::TOPO_PATCH_CHANGE
300 )
301 )
302 {
303 return state1;
304 }
305
306 return state0;
307}
308
309
310// ************************************************************************* //
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
tmp< scalarField > faceSkewness() const
Return face skewness.
Definition: cellQuality.C:243
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
static tmp< scalarField > faceOrthogonality(const polyMesh &mesh, const vectorField &fAreas, const vectorField &cellCtrs)
Definition: polyMeshTools.C:37
static tmp< scalarField > faceWeights(const polyMesh &mesh, const vectorField &fCtrs, const vectorField &fAreas, const vectorField &cellCtrs)
Generate interpolation factors field.
static tmp< scalarField > volRatio(const polyMesh &mesh, const scalarField &vol)
Generate volume ratio field.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1108
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:91
@ TOPO_PATCH_CHANGE
Definition: polyMesh.H:95
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1121
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:456
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1127
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
static scalar boundaryFaceSkewness(const UList< face > &faces, const pointField &p, const vectorField &fCtrs, const vectorField &fAreas, const label facei, const point &ownCc)
Skewness of single boundary face.
static tmp< scalarField > faceOrthogonality(const primitiveMesh &mesh, const vectorField &fAreas, const vectorField &cellCtrs)
Generate non-orthogonality field (internal faces only)
label nInternalFaces() const noexcept
Number of internal faces.
label nFaces() const noexcept
Number of mesh faces.
static void swapBoundaryCellList(const polyMesh &mesh, const UList< T > &cellData, List< T > &neighbourCellData)
Swap to obtain neighbour cell values for all boundary faces.
static void swapBoundaryCellPositions(const polyMesh &mesh, const UList< point > &cellData, List< point > &neighbourCellData)
Swap to obtain neighbour cell positions for all boundary faces.
Definition: syncTools.C:34
A class for managing temporary objects.
Definition: tmp.H:65
volScalarField & p
dynamicFvMesh & mesh
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
dimensionedTensor skew(const dimensionedTensor &dt)
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333