FreeFOAM The Cross-Platform CFD Toolkit
wallPointI.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 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 // Update this with w2 if w2 nearer to pt.
32 inline bool Foam::wallPoint::update
33 (
34  const point& pt,
35  const wallPoint& w2,
36  const scalar tol
37 )
38 {
39  //Already done in calling algorithm
40  //if (w2.origin() == origin_)
41  //{
42  // // Shortcut. Same input so same distance.
43  // return false;
44  //}
45 
46  scalar dist2 = magSqr(pt - w2.origin());
47 
48  if (!valid())
49  {
50  // current not yet set so use any value
51  distSqr_ = dist2;
52  origin_ = w2.origin();
53 
54  return true;
55  }
56 
57  scalar diff = distSqr_ - dist2;
58 
59  if (diff < 0)
60  {
61  // already nearer to pt
62  return false;
63  }
64 
65  if ((diff < SMALL) || ((distSqr_ > SMALL) && (diff/distSqr_ < tol)))
66  {
67  // don't propagate small changes
68  return false;
69  }
70  else
71  {
72  // update with new values
73  distSqr_ = dist2;
74  origin_ = w2.origin();
75 
76  return true;
77  }
78 }
79 
80 
81 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
82 
83 // Null constructor
85 :
86  origin_(greatPoint),
87  distSqr_(GREAT)
88 {}
89 
90 
91 // Construct from origin, distance
92 inline Foam::wallPoint::wallPoint(const point& origin, const scalar distSqr)
93 :
94  origin_(origin), distSqr_(distSqr)
95 {}
96 
97 
98 // Construct as copy
100 :
101  origin_(wpt.origin()), distSqr_(wpt.distSqr())
102 {}
103 
104 
105 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
106 
108 {
109  return origin_;
110 }
111 
112 
114 {
115  return origin_;
116 }
117 
118 
119 inline Foam::scalar Foam::wallPoint::distSqr() const
120 {
121  return distSqr_;
122 }
123 
124 
125 inline Foam::scalar& Foam::wallPoint::distSqr()
126 {
127  return distSqr_;
128 }
129 
130 
131 inline bool Foam::wallPoint::valid() const
132 {
133  return origin_ != greatPoint;
134 }
135 
136 
137 // Checks for cyclic faces
139 (
140  const polyMesh&,
141  const wallPoint& w2,
142  const scalar tol
143 )
144  const
145 {
146  scalar diff = mag(distSqr() - w2.distSqr());
147 
148  if (diff < SMALL)
149  {
150  return true;
151  }
152  else
153  {
154  if ((distSqr() > SMALL) && ((diff/distSqr()) < tol))
155  {
156  return true;
157  }
158  else
159  {
160  return false;
161  }
162  }
163 }
164 
165 
167 (
168  const polyMesh&,
169  const polyPatch&,
170  const label,
171  const point& faceCentre
172 )
173 {
174  origin_ -= faceCentre;
175 }
176 
177 
178 inline void Foam::wallPoint::transform
179 (
180  const polyMesh&,
181  const tensor& rotTensor
182 )
183 {
184  origin_ = Foam::transform(rotTensor, origin_);
185 }
186 
187 
188 // Update absolute geometric quantities. Note that distance (distSqr_)
189 // is not affected by leaving/entering domain.
191 (
192  const polyMesh&,
193  const polyPatch&,
194  const label,
195  const point& faceCentre
196 )
197 {
198  // back to absolute form
199  origin_ += faceCentre;
200 }
201 
202 
203 // Update this with w2 if w2 nearer to pt.
204 inline bool Foam::wallPoint::updateCell
205 (
206  const polyMesh& mesh,
207  const label thisCellI,
208  const label neighbourFaceI,
209  const wallPoint& neighbourWallInfo,
210  const scalar tol
211 )
212 {
213  return
214  update
215  (
216  mesh.cellCentres()[thisCellI],
217  neighbourWallInfo,
218  tol
219  );
220  }
221 
222 
223 // Update this with w2 if w2 nearer to pt.
224 inline bool Foam::wallPoint::updateFace
225 (
226  const polyMesh& mesh,
227  const label thisFaceI,
228  const label neighbourCellI,
229  const wallPoint& neighbourWallInfo,
230  const scalar tol
231 )
232 {
233  return
234  update
235  (
236  mesh.faceCentres()[thisFaceI],
237  neighbourWallInfo,
238  tol
239  );
240 }
241 
242 // Update this with w2 if w2 nearer to pt.
243 inline bool Foam::wallPoint::updateFace
244 (
245  const polyMesh& mesh,
246  const label thisFaceI,
247  const wallPoint& neighbourWallInfo,
248  const scalar tol
249 )
250 {
251  return
252  update
253  (
254  mesh.faceCentres()[thisFaceI],
255  neighbourWallInfo,
256  tol
257  );
258 }
259 
260 
261 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
262 
263 inline bool Foam::wallPoint::operator==(const Foam::wallPoint& rhs) const
264 {
265  return origin() == rhs.origin();
266 }
267 
268 
269 inline bool Foam::wallPoint::operator!=(const Foam::wallPoint& rhs) const
270 {
271  return !(*this == rhs);
272 }
273 
274 
275 // ************************ vim: set sw=4 sts=4 et: ************************ //