FreeFOAM The Cross-Platform CFD Toolkit
pointDataI.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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/polyMesh.H>
27 #include <OpenFOAM/transform.H>
28 #include <meshTools/wallPoint.H>
29 
30 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
31 
32 // Update this with w2 if w2 nearer to pt.
33 inline bool Foam::pointData::update
34 (
35  const point& pt,
36  const pointData& w2,
37  const scalar tol
38 )
39 {
40  scalar dist2 = magSqr(pt - w2.origin());
41 
42  if (!valid())
43  {
44  distSqr_ = dist2;
45  origin_ = w2.origin();
46  s_ = w2.s();
47  v_ = w2.v();
48 
49  return true;
50  }
51 
52 
53 // if (v_ != w2.v())
54 // {
55 // return false;
56 // }
57 
58 
59  scalar diff = distSqr_ - dist2;
60 
61  if (diff < 0)
62  {
63  // already nearer to pt
64  return false;
65  }
66 
67  if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
68  {
69  // don't propagate small changes
70  return false;
71  }
72  else
73  {
74  // update with new values
75  distSqr_ = dist2;
76  origin_ = w2.origin();
77  s_ = w2.s();
78  v_ = w2.v();
79 
80  return true;
81  }
82 }
83 
84 
85 // Update this with w2 (information on same point)
86 inline bool Foam::pointData::update
87 (
88  const pointData& w2,
89  const scalar tol
90 )
91 {
92  if (!valid())
93  {
94  // current not yet set so use any value
95  distSqr_ = w2.distSqr();
96  origin_ = w2.origin();
97  s_ = w2.s();
98  v_ = w2.v();
99  return true;
100  }
101 
102 
103 // if (v_ != w2.v())
104 // {
105 // return false;
106 // }
107 
108 
109  scalar diff = distSqr_ - w2.distSqr();
110 
111  if (diff < 0)
112  {
113  // already nearer to pt
114  return false;
115  }
116 
117  if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
118  {
119  // don't propagate small changes
120  return false;
121  }
122  else
123  {
124  // update with new values
125  distSqr_ = w2.distSqr();
126  origin_ = w2.origin();
127  s_ = w2.s();
128  v_ = w2.v();
129 
130  return true;
131  }
132 }
133 
134 
135 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
136 
137 // Null constructor
139 :
140  origin_(wallPoint::greatPoint),
141  distSqr_(GREAT),
142  s_(GREAT),
143  v_(wallPoint::greatPoint)
144 {}
145 
146 
147 // Construct from origin, distance
149 (
150  const point& origin,
151  const scalar distSqr,
152  const scalar s,
153  const vector& v
154 )
155 :
156  origin_(origin),
157  distSqr_(distSqr),
158  s_(s),
159  v_(v)
160 {}
161 
162 
163 // Construct as copy
165 :
166  origin_(wpt.origin()),
167  distSqr_(wpt.distSqr()),
168  s_(wpt.s()),
169  v_(wpt.v())
170 {}
171 
172 
173 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
174 
176 {
177  return origin_;
178 }
179 
180 
181 inline Foam::scalar Foam::pointData::distSqr() const
182 {
183  return distSqr_;
184 }
185 
186 
187 inline Foam::scalar Foam::pointData::s() const
188 {
189  return s_;
190 }
191 
192 
193 inline const Foam::vector& Foam::pointData::v() const
194 {
195  return v_;
196 }
197 
198 
199 inline bool Foam::pointData::valid() const
200 {
201  return origin_ != wallPoint::greatPoint;
202 }
203 
204 
205 // Checks for cyclic points
207 (
208  const pointData& w2,
209  const scalar tol
210 ) const
211 {
212  scalar diff = Foam::mag(distSqr() - w2.distSqr());
213 
214  if (diff < SMALL)
215  {
216  return true;
217  }
218  else
219  {
220  if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
221  {
222  return true;
223  }
224  else
225  {
226  return false;
227  }
228  }
229 }
230 
231 
233 (
234  const polyPatch& patch,
235  const label patchPointI,
236  const point& coord
237 )
238 {
239  origin_ -= coord;
240 }
241 
242 
243 inline void Foam::pointData::transform(const tensor& rotTensor)
244 {
245  origin_ = Foam::transform(rotTensor, origin_);
246 }
247 
248 
249 // Update absolute geometric quantities. Note that distance (distSqr_)
250 // is not affected by leaving/entering domain.
252 (
253  const polyPatch& patch,
254  const label patchPointI,
255  const point& coord
256 )
257 {
258  // back to absolute form
259  origin_ += coord;
260 }
261 
262 
263 // Update this with information from connected edge
265 (
266  const polyMesh& mesh,
267  const label pointI,
268  const label edgeI,
269  const pointData& edgeInfo,
270  const scalar tol
271 )
272 {
273  return
274  update
275  (
276  mesh.points()[pointI],
277  edgeInfo,
278  tol
279  );
280  }
281 
282 
283 // Update this with new information on same point
285 (
286  const polyMesh& mesh,
287  const label pointI,
288  const pointData& newPointInfo,
289  const scalar tol
290 )
291 {
292  return
293  update
294  (
295  mesh.points()[pointI],
296  newPointInfo,
297  tol
298  );
299 }
300 
301 
302 // Update this with new information on same point. No extra information.
304 (
305  const pointData& newPointInfo,
306  const scalar tol
307 )
308 {
309  return update(newPointInfo, tol);
310 }
311 
312 
313 // Update this with information from connected point
314 inline bool Foam::pointData::updateEdge
315 (
316  const polyMesh& mesh,
317  const label edgeI,
318  const label pointI,
319  const pointData& pointInfo,
320  const scalar tol
321 )
322 {
323  const pointField& points = mesh.points();
324 
325  const edge& e = mesh.edges()[edgeI];
326 
327  const point edgeMid(0.5*(points[e[0]] + points[e[1]]));
328 
329  return
330  update
331  (
332  edgeMid,
333  pointInfo,
334  tol
335  );
336 }
337 
338 
339 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
340 
341 inline bool Foam::pointData::operator==(const pointData& rhs) const
342 {
343  return origin() == rhs.origin();
344 }
345 
346 
347 inline bool Foam::pointData::operator!=(const pointData& rhs) const
348 {
349  return !(*this == rhs);
350 }
351 
352 
353 // ************************ vim: set sw=4 sts=4 et: ************************ //