pointEdgeStructuredWalkI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2018 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 // Update this with w2.
32 template<class TrackingData>
33 inline bool Foam::pointEdgeStructuredWalk::update
34 (
35  const pointEdgeStructuredWalk& w2,
36  const scalar tol,
37  TrackingData& td
38 )
39 {
40  if (!valid(td))
41  {
42  // current not yet set. Walked from w2 to here (=point0)
43  dist_ = w2.dist_ + mag(point0_-w2.previousPoint_);
44  previousPoint_ = point0_;
45  data_ = w2.data_;
46 
47  return true;
48  }
49  else
50  {
51  return false;
52  }
53 }
54 
55 
56 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
57 
58 // Null constructor
60 :
61  point0_(vector::max),
62  previousPoint_(vector::max),
63  dist_(0),
64  data_(Zero)
65 {}
66 
67 
68 // Construct from origin, distance
70 (
71  const point& point0,
72  const point& previousPoint,
73  const scalar dist,
74  const vector& data
75 )
76 :
77  point0_(point0),
78  previousPoint_(previousPoint),
79  dist_(dist),
80  data_(data)
81 {}
82 
83 
84 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 
87 {
88  return point0_ != vector::max;
89 }
90 
91 
92 //inline const Foam::point& Foam::pointEdgeStructuredWalk::previousPoint() const
93 //{
94 // return previousPoint_;
95 //}
96 
97 
98 inline Foam::scalar Foam::pointEdgeStructuredWalk::dist() const
99 {
100  return dist_;
101 }
102 
103 
105 {
106  return data_;
107 }
108 
109 
110 template<class TrackingData>
111 inline bool Foam::pointEdgeStructuredWalk::valid(TrackingData& td) const
112 {
113  return previousPoint_ != vector::max;
114 }
115 
116 
117 // Checks for cyclic points
118 template<class TrackingData>
120 (
121  const pointEdgeStructuredWalk& w2,
122  const scalar tol,
123  TrackingData& td
124 ) const
125 {
126  scalar diff = Foam::mag(dist() - w2.dist());
127 
128  if (diff < small)
129  {
130  return true;
131  }
132  else
133  {
134  if ((dist() > small) && ((diff/dist()) < tol))
135  {
136  return true;
137  }
138  else
139  {
140  return false;
141  }
142  }
143 }
144 
145 
146 template<class TrackingData>
148 (
149  const polyPatch& patch,
150  const label patchPointi,
151  const point& coord,
152  TrackingData& td
153 )
154 {
155  previousPoint_ -= coord;
156 }
157 
158 
159 template<class TrackingData>
161 (
162  const tensor& rotTensor,
163  TrackingData& td
164 )
165 {
166  previousPoint_ = Foam::transform(rotTensor, previousPoint_);
167 }
168 
169 
170 // Update absolute geometric quantities. Note that distance (dist_)
171 // is not affected by leaving/entering domain.
172 template<class TrackingData>
174 (
175  const polyPatch& patch,
176  const label patchPointi,
177  const point& coord,
178  TrackingData& td
179 )
180 {
181  // back to absolute form
182  previousPoint_ += coord;
183 }
184 
185 
186 // Update this with information from connected edge
187 template<class TrackingData>
189 (
190  const polyMesh& mesh,
191  const label pointi,
192  const label edgeI,
193  const pointEdgeStructuredWalk& edgeInfo,
194  const scalar tol,
195  TrackingData& td
196 )
197 {
198  if (inZone())
199  {
200  return update(edgeInfo, tol, td);
201  }
202  else
203  {
204  return false;
205  }
206 }
207 
208 
209 // Update this with new information on same point
210 template<class TrackingData>
212 (
213  const polyMesh& mesh,
214  const label pointi,
215  const pointEdgeStructuredWalk& newPointInfo,
216  const scalar tol,
217  TrackingData& td
218 )
219 {
220  if (inZone())
221  {
222  return update(newPointInfo, tol, td);
223  }
224  else
225  {
226  return false;
227  }
228 }
229 
230 
231 // Update this with new information on same point. No extra information.
232 template<class TrackingData>
234 (
235  const pointEdgeStructuredWalk& newPointInfo,
236  const scalar tol,
237  TrackingData& td
238 )
239 {
240  return update(newPointInfo, tol, td);
241 }
242 
243 
244 // Update this with information from connected point
245 template<class TrackingData>
247 (
248  const polyMesh& mesh,
249  const label edgeI,
250  const label pointi,
251  const pointEdgeStructuredWalk& pointInfo,
252  const scalar tol,
253  TrackingData& td
254 )
255 {
256  if (inZone())
257  {
258  return update(pointInfo, tol, td);
259  }
260  else
261  {
262  return false;
263  }
264 }
265 
266 
267 template<class TrackingData>
269 (
270  const pointEdgeStructuredWalk& rhs,
271  TrackingData& td
272 ) const
273 {
274  return operator==(rhs);
275 }
276 
277 
278 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
279 
280 inline bool Foam::pointEdgeStructuredWalk::operator==
281 (
283 ) const
284 {
285  return previousPoint_ == rhs.previousPoint_;
286 }
287 
288 
289 inline bool Foam::pointEdgeStructuredWalk::operator!=
290 (
292 ) const
293 {
294  return !(*this == rhs);
295 }
296 
297 
298 // ************************************************************************* //
bool sameGeometry(const pointEdgeStructuredWalk &, const scalar tol, TrackingData &td) const
Check for identical geometrical data. Used for cyclics checking.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:403
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
bool updatePoint(const polyMesh &mesh, const label pointi, const label edgeI, const pointEdgeStructuredWalk &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on point.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
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
Check whether origin has been changed at all or.
void transform(const tensor &rotTensor, TrackingData &td)
Apply rotation matrix to origin.
Determines length of string of edges walked to point.
void leaveDomain(const polyPatch &patch, const label patchPointi, const point &pos, TrackingData &td)
Convert origin to relative vector to leaving point.
3D tensor transformation operations.
static const zero Zero
Definition: zero.H:97
bool equal(const pointEdgeStructuredWalk &, TrackingData &) const
Same (like operator==)
bool operator==(const pointEdgeStructuredWalk &) const
bool updateEdge(const polyMesh &mesh, const label edgeI, const label pointi, const pointEdgeStructuredWalk &pointInfo, const scalar tol, TrackingData &td)
Influence of point on edge.
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
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477