patchEdgeFaceInfoI.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 if w2 nearer to pt.
32 template<class TrackingData>
33 inline bool Foam::patchEdgeFaceInfo::update
34 (
35  const point& pt,
36  const patchEdgeFaceInfo& w2,
37  const scalar tol,
38  TrackingData& td
39 )
40 {
41  scalar dist2 = magSqr(pt - w2.origin());
42 
43  if (!valid(td))
44  {
45  // current not yet set so use any value
46  distSqr_ = dist2;
47  origin_ = w2.origin();
48 
49  return true;
50  }
51 
52  scalar diff = distSqr_ - dist2;
53 
54  if (diff < 0)
55  {
56  // already nearer to pt
57  return false;
58  }
59 
60  if ((diff < small) || ((distSqr_ > small) && (diff/distSqr_ < tol)))
61  {
62  // don't propagate small changes
63  return false;
64  }
65  else
66  {
67  // update with new values
68  distSqr_ = dist2;
69  origin_ = w2.origin();
70 
71  return true;
72  }
73 }
74 
75 
76 // Update this with w2 (information on same edge)
77 template<class TrackingData>
78 inline bool Foam::patchEdgeFaceInfo::update
79 (
80  const patchEdgeFaceInfo& w2,
81  const scalar tol,
82  TrackingData& td
83 )
84 {
85  if (!valid(td))
86  {
87  // current not yet set so use any value
88  distSqr_ = w2.distSqr();
89  origin_ = w2.origin();
90 
91  return true;
92  }
93 
94  scalar diff = distSqr_ - w2.distSqr();
95 
96  if (diff < 0)
97  {
98  // already nearer to pt
99  return false;
100  }
101 
102  if ((diff < small) || ((distSqr_ > small) && (diff/distSqr_ < tol)))
103  {
104  // don't propagate small changes
105  return false;
106  }
107  else
108  {
109  // update with new values
110  distSqr_ = w2.distSqr();
111  origin_ = w2.origin();
112 
113  return true;
114  }
115 }
116 
117 
118 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
119 
120 // Null constructor
122 :
123  origin_(point::max),
124  distSqr_(sqr(great))
125 {}
126 
127 
128 // Construct from origin, distance
130 (
131  const point& origin,
132  const scalar distSqr
133 )
134 :
135  origin_(origin),
136  distSqr_(distSqr)
137 {}
138 
139 
140 // Construct as copy
142 :
143  origin_(wpt.origin()),
144  distSqr_(wpt.distSqr())
145 {}
146 
147 
148 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
149 
151 {
152  return origin_;
153 }
154 
155 
156 inline Foam::scalar Foam::patchEdgeFaceInfo::distSqr() const
157 {
158  return distSqr_;
159 }
160 
161 
162 template<class TrackingData>
163 inline bool Foam::patchEdgeFaceInfo::valid(TrackingData& td) const
164 {
165  return origin_ != point::max;
166 }
167 
168 
169 template<class TrackingData>
171 (
172  const polyMesh& mesh,
173  const primitivePatch& patch,
174  const tensor& rotTensor,
175  const scalar tol,
176  TrackingData& td
177 )
178 {
179  origin_ = Foam::transform(rotTensor, origin_);
180 }
181 
182 
183 template<class TrackingData>
185 (
186  const polyMesh& mesh,
187  const primitivePatch& patch,
188  const label edgeI,
189  const label facei,
190  const patchEdgeFaceInfo& faceInfo,
191  const scalar tol,
192  TrackingData& td
193 )
194 {
195  const edge& e = patch.edges()[edgeI];
196  point eMid =
197  0.5
198  * (
199  patch.points()[patch.meshPoints()[e[0]]]
200  + patch.points()[patch.meshPoints()[e[1]]]
201  );
202  return update(eMid, faceInfo, tol, td);
203 }
204 
205 
206 template<class TrackingData>
208 (
209  const polyMesh& mesh,
210  const primitivePatch& patch,
211  const patchEdgeFaceInfo& edgeInfo,
212  const bool sameOrientation,
213  const scalar tol,
214  TrackingData& td
215 )
216 {
217  return update(edgeInfo, tol, td);
218 }
219 
220 
221 template<class TrackingData>
223 (
224  const polyMesh& mesh,
225  const primitivePatch& patch,
226  const label facei,
227  const label edgeI,
228  const patchEdgeFaceInfo& edgeInfo,
229  const scalar tol,
230  TrackingData& td
231 )
232 {
233  return update(patch.faceCentres()[facei], edgeInfo, tol, td);
234 }
235 
236 
237 template<class TrackingData>
239 (
240  const patchEdgeFaceInfo& rhs,
241  TrackingData& td
242 ) const
243 {
244  return operator==(rhs);
245 }
246 
247 
248 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
249 
250 inline bool Foam::patchEdgeFaceInfo::operator==
251 (
252  const Foam::patchEdgeFaceInfo& rhs
253 ) const
254 {
255  return origin() == rhs.origin();
256 }
257 
258 
259 inline bool Foam::patchEdgeFaceInfo::operator!=
260 (
261  const Foam::patchEdgeFaceInfo& rhs
262 ) const
263 {
264  return !(*this == rhs);
265 }
266 
267 
268 // ************************************************************************* //
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
bool operator==(const patchEdgeFaceInfo &) const
patchEdgeFaceInfo()
Construct null.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const Field< PointType > & faceCentres() const
Return face centres for patch.
bool equal(const patchEdgeFaceInfo &, TrackingData &td) const
Same (like operator==)
bool updateFace(const polyMesh &mesh, const primitivePatch &patch, const label facei, const label edgeI, const patchEdgeFaceInfo &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on face.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
A list of faces which address into the list of points.
const point & origin() const
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.
bool updateEdge(const polyMesh &mesh, const primitivePatch &patch, const label edgeI, const label facei, const patchEdgeFaceInfo &faceInfo, const scalar tol, TrackingData &td)
Influence of face on edge.
const Field< PointType > & points() const
Return reference to global points.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
vector point
Point is a vector.
Definition: point.H:41
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:98
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477
void transform(const polyMesh &mesh, const primitivePatch &patch, const tensor &rotTensor, const scalar tol, TrackingData &td)
Apply rotation matrix.