pointEdgePointI.H
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) 2020 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 "polyMesh.H"
30#include "transform.H"
31
32// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33
34// Update this with w2 if w2 nearer to pt.
35template<class TrackingData>
37(
38 const point& pt,
39 const pointEdgePoint& w2,
40 const scalar tol,
41 TrackingData& td
42)
43{
44 const scalar dist2 = magSqr(pt - w2.origin());
45
46 if (!valid(td))
47 {
48 // current not yet set so use any value
49 distSqr_ = dist2;
50 origin_ = w2.origin();
51
52 return true;
53 }
54
55 const scalar diff = distSqr_ - dist2;
56
57 if (diff < 0)
58 {
59 // already nearer to pt
60 return false;
61 }
62
63 if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
64 {
65 // don't propagate small changes
66 return false;
67 }
68 else
69 {
70 // update with new values
71 distSqr_ = dist2;
72 origin_ = w2.origin();
73
74 return true;
75 }
76}
77
78
79// Update this with w2 (information on same point)
80template<class TrackingData>
82(
83 const pointEdgePoint& w2,
84 const scalar tol,
85 TrackingData& td
86)
87{
88 if (!valid(td))
89 {
90 // current not yet set so use any value
91 distSqr_ = w2.distSqr();
92 origin_ = w2.origin();
93
94 return true;
95 }
96
97 const scalar diff = distSqr_ - w2.distSqr();
98
99 if (diff < 0)
100 {
101 // already nearer to pt
102 return false;
103 }
104
105 if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
106 {
107 // don't propagate small changes
108 return false;
109 }
110 else
111 {
112 // update with new values
113 distSqr_ = w2.distSqr();
114 origin_ = w2.origin();
115
116 return true;
117 }
118}
119
120
121// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
122
124:
125 origin_(point::max),
126 distSqr_(GREAT)
127{}
128
129
131(
132 const point& origin,
133 const scalar distSqr
134)
135:
136 origin_(origin),
137 distSqr_(distSqr)
138{}
139
140
141// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
142
143template<class TrackingData>
144inline bool Foam::pointEdgePoint::valid(TrackingData& td) const
145{
146 return origin_ != point::max;
147}
148
149
150// Checks for cyclic points
151template<class TrackingData>
153(
154 const pointEdgePoint& w2,
155 const scalar tol,
156 TrackingData& td
157) const
158{
159 const scalar diff = Foam::mag(distSqr() - w2.distSqr());
160
161 if (diff < SMALL)
162 {
163 return true;
164 }
165 else
166 {
167 if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
168 {
169 return true;
170 }
171 else
172 {
173 return false;
174 }
175 }
176}
177
178
179template<class TrackingData>
181(
182 const polyPatch& patch,
183 const label patchPointi,
184 const point& coord,
185 TrackingData& td
186)
187{
188 origin_ -= coord;
189}
190
191
192template<class TrackingData>
194(
195 const tensor& rotTensor,
196 TrackingData& td
197)
198{
199 origin_ = Foam::transform(rotTensor, origin_);
200}
201
202
203// Update absolute geometric quantities. Note that distance (distSqr_)
204// is not affected by leaving/entering domain.
205template<class TrackingData>
207(
208 const polyPatch& patch,
209 const label patchPointi,
210 const point& coord,
211 TrackingData& td
212)
213{
214 // back to absolute form
215 origin_ += coord;
216}
217
218
219// Update this with information from connected edge
220template<class TrackingData>
222(
223 const polyMesh& mesh,
224 const label pointi,
225 const label edgeI,
226 const pointEdgePoint& edgeInfo,
227 const scalar tol,
228 TrackingData& td
229)
230{
231 return update(mesh.points()[pointi], edgeInfo, tol, td);
232}
233
234
235// Update this with new information on same point
236template<class TrackingData>
238(
239 const polyMesh& mesh,
240 const label pointi,
241 const pointEdgePoint& newPointInfo,
242 const scalar tol,
243 TrackingData& td
244)
245{
246 return update(mesh.points()[pointi], newPointInfo, tol, td);
247}
248
249
250// Update this with new information on same point. No extra information.
251template<class TrackingData>
253(
254 const pointEdgePoint& newPointInfo,
255 const scalar tol,
256 TrackingData& td
257)
258{
259 return update(newPointInfo, tol, td);
260}
261
262
263// Update this with information from connected point
264template<class TrackingData>
266(
267 const polyMesh& mesh,
268 const label edgeI,
269 const label pointi,
270 const pointEdgePoint& pointInfo,
271 const scalar tol,
272 TrackingData& td
273)
274{
275 const edge& e = mesh.edges()[edgeI];
276 return update(e.centre(mesh.points()), pointInfo, tol, td);
277}
278
279
280template<class TrackingData>
282(
283 const pointEdgePoint& rhs,
284 TrackingData& td
285) const
286{
287 return operator==(rhs);
288}
289
290
291// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
292
294(
295 const pointEdgePoint& rhs
296) const
297{
298 return origin_ == rhs.origin_ && distSqr_ == rhs.distSqr_;
299}
300
301
303(
304 const pointEdgePoint& rhs
305) const
306{
307 return !(*this == rhs);
308}
309
310
311// ************************************************************************* //
#define w2
Definition: blockCreate.C:35
bool valid() const
True if all internal ids are non-negative.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:66
friend Ostream & operator(Ostream &, const faMatrix< Type > &)
virtual bool update()
Update the mesh for both mesh motion and topology change.
Default transformation behaviour.
static const complex max
complex (VGREAT,VGREAT)
Definition: complex.H:293
Holds information regarding nearest wall point. Used in PointEdgeWave. (so not standard FaceCellWave)...
bool sameGeometry(const pointEdgePoint &, const scalar tol, TrackingData &td) const
Check for identical geometrical data (eg, cyclics checking)
bool equal(const pointEdgePoint &, TrackingData &td) const
Test for equality, with TrackingData.
pointEdgePoint()
Default construct.
void leaveDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert origin to relative vector to leaving point.
void enterDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert relative origin to absolute by adding entering point.
bool valid(TrackingData &td) const
Changed or contains original (invalid) value.
bool updatePoint(const polyMesh &mesh, const label pointi, const label edgeI, const pointEdgePoint &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on point.
bool updateEdge(const polyMesh &mesh, const label edgeI, const label pointi, const pointEdgePoint &pointInfo, const scalar tol, TrackingData &td)
Influence of point on edge.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:75
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
mesh update()
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
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:536
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:378
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
volScalarField & e
Definition: createFields.H:11
3D tensor transformation operations.