externalPointEdgePointI.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 | Copyright (C) 2013-2016 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 "polyMesh.H"
27 #include "transform.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class TrackingData>
32 inline bool Foam::externalPointEdgePoint::update
33 (
34  const point& pt,
35  const externalPointEdgePoint& w2,
36  const scalar tol,
37  TrackingData& td
38 )
39 {
40  scalar dist2 = magSqr(pt - w2.origin());
41 
42  if (!valid(td))
43  {
44  // current not yet set so use any value
45  distSqr_ = dist2;
46  origin_ = w2.origin();
47 
48  return true;
49  }
50 
51  scalar diff = distSqr_ - dist2;
52 
53  if (diff < 0)
54  {
55  // already nearer to pt
56  return false;
57  }
58 
59  if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
60  {
61  // don't propagate small changes
62  return false;
63  }
64  else
65  {
66  // update with new values
67  distSqr_ = dist2;
68  origin_ = w2.origin();
69 
70  return true;
71  }
72 }
73 
74 
75 template<class TrackingData>
76 inline bool Foam::externalPointEdgePoint::update
77 (
78  const externalPointEdgePoint& w2,
79  const scalar tol,
80  TrackingData& td
81 )
82 {
83  if (!valid(td))
84  {
85  // current not yet set so use any value
86  distSqr_ = w2.distSqr();
87  origin_ = w2.origin();
88 
89  return true;
90  }
91 
92  scalar diff = distSqr_ - w2.distSqr();
93 
94  if (diff < 0)
95  {
96  // already nearer to pt
97  return false;
98  }
99 
100  if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
101  {
102  // don't propagate small changes
103  return false;
104  }
105  else
106  {
107  // update with new values
108  distSqr_ = w2.distSqr();
109  origin_ = w2.origin();
110 
111  return true;
112  }
113 }
114 
115 
116 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
117 
119 :
120  origin_(point::max),
121  distSqr_(GREAT)
122 {}
123 
124 
126 (
127  const point& origin,
128  const scalar distSqr
129 )
130 :
131  origin_(origin),
132  distSqr_(distSqr)
133 {}
134 
135 
137 (
138  const externalPointEdgePoint& wpt
139 )
140 :
141  origin_(wpt.origin()),
142  distSqr_(wpt.distSqr())
143 {}
144 
145 
146 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
147 
149 {
150  return origin_;
151 }
152 
153 
154 inline Foam::scalar Foam::externalPointEdgePoint::distSqr() const
155 {
156  return distSqr_;
157 }
158 
159 
160 template<class TrackingData>
161 inline bool Foam::externalPointEdgePoint::valid(TrackingData& td) const
162 {
163  return origin_ != point::max;
164 }
165 
166 
167 // Checks for cyclic points
168 template<class TrackingData>
170 (
171  const externalPointEdgePoint& w2,
172  const scalar tol,
173  TrackingData& td
174 ) const
175 {
176  scalar diff = Foam::mag(distSqr() - w2.distSqr());
177 
178  if (diff < SMALL)
179  {
180  return true;
181  }
182  else
183  {
184  if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
185  {
186  return true;
187  }
188  else
189  {
190  return false;
191  }
192  }
193 }
194 
195 
196 template<class TrackingData>
198 (
199  const polyPatch& patch,
200  const label patchPointi,
201  const point& coord,
202  TrackingData& td
203 )
204 {
205  origin_ -= coord;
206 }
207 
208 
209 template<class TrackingData>
211 (
212  const tensor& rotTensor,
213  TrackingData& td
214 )
215 {
216  origin_ = Foam::transform(rotTensor, origin_);
217 }
218 
219 
220 template<class TrackingData>
222 (
223  const polyPatch& patch,
224  const label patchPointi,
225  const point& coord,
226  TrackingData& td
227 )
228 {
229  // back to absolute form
230  origin_ += coord;
231 }
232 
233 
234 template<class TrackingData>
236 (
237  const polyMesh& mesh,
238  const label pointi,
239  const label edgeI,
240  const externalPointEdgePoint& edgeInfo,
241  const scalar tol,
242  TrackingData& td
243 )
244 {
245  return update(td.points_[pointi], edgeInfo, tol, td);
246 }
247 
248 
249 template<class TrackingData>
251 (
252  const polyMesh& mesh,
253  const label pointi,
254  const externalPointEdgePoint& newPointInfo,
255  const scalar tol,
256  TrackingData& td
257 )
258 {
259  return update(td.points_[pointi], newPointInfo, tol, td);
260 }
261 
262 
263 template<class TrackingData>
265 (
266  const externalPointEdgePoint& newPointInfo,
267  const scalar tol,
268  TrackingData& td
269 )
270 {
271  return update(newPointInfo, tol, td);
272 }
273 
274 
275 template<class TrackingData>
277 (
278  const polyMesh& mesh,
279  const label edgeI,
280  const label pointi,
281  const externalPointEdgePoint& pointInfo,
282  const scalar tol,
283  TrackingData& td
284 )
285 {
286  const edge& e = mesh.edges()[edgeI];
287  return update(e.centre(td.points_), pointInfo, tol, td);
288 }
289 
290 
291 template<class TrackingData>
293 (
294  const externalPointEdgePoint& rhs,
295  TrackingData& td
296 ) const
297 {
298  return operator==(rhs);
299 }
300 
301 
302 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
303 
304 inline bool Foam::externalPointEdgePoint::operator==
305 (
307 )
308 const
309 {
310  return (origin() == rhs.origin()) && (distSqr() == rhs.distSqr());
311 }
312 
313 
314 inline bool Foam::externalPointEdgePoint::operator!=
315 (
317 )
318 const
319 {
320  return !(*this == rhs);
321 }
322 
323 
324 // ************************************************************************* //
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
const double e
Elementary charge.
Definition: doubleFloat.H:78
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
point centre(const pointField &) const
Return centre (centroid)
Definition: edgeI.H:151
bool operator==(const externalPointEdgePoint &) const
bool updatePoint(const polyMesh &mesh, const label pointi, const label edgeI, const externalPointEdgePoint &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on point.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
Holds information regarding nearest wall point. Used in PointEdgeWave. (so not standard FaceCellWave)...
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
3D tensor transformation operations.
void transform(const tensor &rotTensor, TrackingData &td)
Apply rotation matrix to origin.
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
bool sameGeometry(const externalPointEdgePoint &, const scalar tol, TrackingData &td) const
Check for identical geometrical data. Used for cyclics checking.
vector point
Point is a vector.
Definition: point.H:41
dimensioned< scalar > mag(const dimensioned< Type > &)
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
void leaveDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert origin to relative vector to leaving point.
bool equal(const externalPointEdgePoint &, TrackingData &td) const
Equivalent to operator== with TrackingData.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465
bool updateEdge(const polyMesh &mesh, const label edgeI, const label pointi, const externalPointEdgePoint &pointInfo, const scalar tol, TrackingData &td)
Influence of point on edge.
void enterDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert relative origin to absolute by adding entering point.