patchEdgeFacePointI.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-2023 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::patchEdgeFacePoint::update
34 (
35  const point& pt,
36  const patchEdgeFacePoint& 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::patchEdgeFacePoint::update
79 (
80  const patchEdgeFacePoint& 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 
121 :
122  origin_(point::max),
123  distSqr_(sqr(great))
124 {}
125 
126 
128 (
129  const point& origin,
130  const scalar distSqr
131 )
132 :
133  origin_(origin),
134  distSqr_(distSqr)
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
141 {
142  return origin_;
143 }
144 
145 
146 inline Foam::scalar Foam::patchEdgeFacePoint::distSqr() const
147 {
148  return distSqr_;
149 }
150 
151 
152 template<class TrackingData>
153 inline bool Foam::patchEdgeFacePoint::valid(TrackingData& td) const
154 {
155  return origin_ != point::max;
156 }
157 
158 
159 template<class TrackingData>
161 (
162  const polyMesh& mesh,
163  const primitivePatch& patch,
164  const tensor& rotTensor,
165  const scalar tol,
166  TrackingData& td
167 )
168 {
169  origin_ = Foam::transform(rotTensor, origin_);
170 }
171 
172 
173 template<class TrackingData>
175 (
176  const polyMesh& mesh,
177  const primitivePatch& patch,
178  const label edgei,
179  const label facei,
180  const patchEdgeFacePoint& faceInfo,
181  const scalar tol,
182  TrackingData& td
183 )
184 {
185  const edge& e = patch.edges()[edgei];
186  point eMid =
187  0.5
188  * (
189  patch.points()[patch.meshPoints()[e[0]]]
190  + patch.points()[patch.meshPoints()[e[1]]]
191  );
192  return update(eMid, faceInfo, tol, td);
193 }
194 
195 
196 template<class TrackingData>
198 (
199  const polyMesh& mesh,
200  const primitivePatch& patch,
201  const patchEdgeFacePoint& edgeInfo,
202  const bool sameOrientation,
203  const scalar tol,
204  TrackingData& td
205 )
206 {
207  return update(edgeInfo, tol, td);
208 }
209 
210 
211 template<class TrackingData>
213 (
214  const polyMesh& mesh,
215  const primitivePatch& patch,
216  const label facei,
217  const label edgei,
218  const patchEdgeFacePoint& edgeInfo,
219  const scalar tol,
220  TrackingData& td
221 )
222 {
223  return update(patch.faceCentres()[facei], edgeInfo, tol, td);
224 }
225 
226 
227 template<class TrackingData>
229 (
230  const patchEdgeFacePoint& rhs,
231  TrackingData& td
232 ) const
233 {
234  return operator==(rhs);
235 }
236 
237 
238 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
239 
240 inline bool Foam::patchEdgeFacePoint::operator==
241 (
242  const Foam::patchEdgeFacePoint& rhs
243 ) const
244 {
245  return origin() == rhs.origin();
246 }
247 
248 
249 inline bool Foam::patchEdgeFacePoint::operator!=
250 (
251  const Foam::patchEdgeFacePoint& rhs
252 ) const
253 {
254  return !(*this == rhs);
255 }
256 
257 
258 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
259 
260 inline Foam::Ostream& Foam::operator<<
261 (
262  Foam::Ostream& os,
263  const Foam::patchEdgeFacePoint& wDist
264 )
265 {
266  return os << wDist.origin() << token::SPACE << wDist.distSqr();
267 }
268 
269 
270 inline Foam::Istream& Foam::operator>>
271 (
272  Foam::Istream& is,
274 )
275 {
276  return is >> wDist.origin_ >> wDist.distSqr_;
277 }
278 
279 
280 // ************************************************************************* //
#define w2
Definition: blockCreate.C:32
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
const Field< PointType > & points() const
Return reference to global points.
const Field< PointType > & faceCentres() const
Return face centres for patch.
static const Form max
Definition: VectorSpace.H:120
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:61
Transport of nearest point location for use in PatchEdgeFaceWave.
patchEdgeFacePoint()
Construct null.
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
bool equal(const patchEdgeFacePoint &, TrackingData &td) const
Same (like operator==)
bool updateFace(const polyMesh &mesh, const primitivePatch &patch, const label facei, const label edgei, const patchEdgeFacePoint &edgeInfo, const scalar tol, TrackingData &td)
Influence of edge on face.
const point & origin() const
bool updateEdge(const polyMesh &mesh, const primitivePatch &patch, const label edgei, const label facei, const patchEdgeFacePoint &faceInfo, const scalar tol, TrackingData &td)
Influence of face on edge.
void transform(const polyMesh &mesh, const primitivePatch &patch, const tensor &rotTensor, const scalar tol, TrackingData &td)
Apply rotation matrix.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:80
bool valid(const PtrList< ModelType > &l)
const doubleScalar e
Definition: doubleScalar.H:106
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
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
vector point
Point is a vector.
Definition: point.H:41
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:409
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:504
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
3D tensor transformation operations.