wallFaceI.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 "wallFace.H"
27 #include "polyMesh.H"
28 #include "transformer.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class TrackingData>
34 (
35  const point& pt,
36  const wallFace& w2,
37  const scalar tol,
38  TrackingData& td
39 )
40 {
41  const scalar dist2 =
42  sqr
43  (
44  face(identityMap(w2.points().size()))
45  .nearestPoint(pt, w2.points())
46  .distance()
47  );
48 
49  if (valid(td))
50  {
51  scalar diff = distSqr() - dist2;
52 
53  if (diff < 0)
54  {
55  // already nearer to pt
56  return false;
57  }
58 
59  if ((diff < small) || ((distSqr() > small) && (diff/distSqr() < tol)))
60  {
61  // don't propagate small changes
62  return false;
63  }
64  }
65 
66  // Either *this is not yet valid or w2 is closer
67  distSqr() = dist2;
68  points() = w2.points();
69 
70  return true;
71 }
72 
73 
74 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
75 
77 :
78  points_(0),
79  distSqr_(-great)
80 {}
81 
82 
84 (
85  const face& f,
86  const pointField& points,
87  const scalar distSqr
88 )
89 :
90  points_(f.points(points)),
91  distSqr_(distSqr)
92 {}
93 
94 
96 (
97  const face& f,
98  const pointField& points,
99  const point& centre,
100  const scalar distSqr
101 )
102 :
103  points_(f.points(points)),
104  distSqr_(distSqr)
105 {}
106 
107 
108 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
109 
111 {
112  return points_;
113 }
114 
115 
117 {
118  return points_;
119 }
120 
121 
122 inline Foam::scalar Foam::wallFace::distSqr() const
123 {
124  return distSqr_;
125 }
126 
127 
128 inline Foam::scalar& Foam::wallFace::distSqr()
129 {
130  return distSqr_;
131 }
132 
133 
134 template<class TrackingData>
135 inline Foam::scalar Foam::wallFace::dist(TrackingData& td) const
136 {
137  return valid(td) ? sqrt(distSqr_) : great;
138 }
139 
140 
141 template<class TrackingData>
142 inline bool Foam::wallFace::valid(TrackingData& td) const
143 {
144  return distSqr_ > -small;
145 }
146 
147 
148 template<class TrackingData>
150 (
151  const wallFace& w2,
152  const scalar tol,
153  TrackingData& td
154 ) const
155 {
156  const scalar diff = mag(distSqr() - w2.distSqr());
157 
158  return
159  diff < small
160  || ((distSqr() > small) && (diff/distSqr() < tol));
161 }
162 
163 
164 template<class TrackingData>
166 (
167  const transformer& transform,
168  TrackingData& td
169 )
170 {
171  // Note that distSqr_ is not affected by crossing an interface
172  transform.transformPosition(points_, points_);
173 }
174 
175 
176 template<class TrackingData>
178 (
179  const wallFace& rhs,
180  TrackingData& td
181 ) const
182 {
183  return operator==(rhs);
184 }
185 
186 
187 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
188 
189 inline bool Foam::wallFace::operator==
190 (
191  const Foam::wallFace& rhs
192 ) const
193 {
194  return points() == rhs.points();
195 }
196 
197 
198 inline bool Foam::wallFace::operator!=
199 (
200  const Foam::wallFace& rhs
201 ) const
202 {
203  return !(*this == rhs);
204 }
205 
206 
207 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
208 
210 {
211  return os << w.points() << token::SPACE << w.distSqr();
212 }
213 
214 
216 {
217  return is >> w.points() >> w.distSqr();
218 }
219 
220 
221 // ************************************************************************* //
#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
scalar distance() const
Return distance to hit.
Definition: PointHit.H:139
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
pointHit nearestPoint(const point &p, const pointField &) const
Return nearest point to face.
Vector-tensor class used to perform translations, rotations and scaling operations in 3D space.
Definition: transformer.H:84
Holds information regarding nearest wall point. Used in wall distance calculation.
Definition: wallFace.H:59
wallFace()
Construct null.
Definition: wallFaceI.H:76
scalar dist(TrackingData &td) const
bool equal(const wallFace &, TrackingData &td) const
Test equality.
Definition: wallFaceI.H:178
void transform(const transformer &transform, TrackingData &td)
Transform across an interface.
Definition: wallFaceI.H:166
scalar distSqr() const
Definition: wallFaceI.H:122
bool update(const point &pt, const wallFace &w2, const scalar tol, TrackingData &td)
Evaluate distance to point. Update distSqr, origin from whomever.
Definition: wallFaceI.H:34
bool valid(TrackingData &td) const
Check whether the wallFace has been changed at all or still.
Definition: wallFaceI.H:142
const pointField & points() const
Definition: wallFaceI.H:110
bool sameGeometry(const wallFace &, const scalar, TrackingData &td) const
Check for identical geometrical data. Used for checking.
Definition: wallFaceI.H:150
const pointField & points
bool valid(const PtrList< ModelType > &l)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Istream & operator>>(Istream &, pistonPointEdgeData &)
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
dimensioned< scalar > mag(const dimensioned< Type > &)
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
Definition: ListOps.C:104
labelList f(nPoints)