wallPointI.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 "wallPoint.H"
27 #include "transformer.H"
28 #include "SubField.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class TrackingData>
34 (
35  const point& pt,
36  const wallPoint& w2,
37  const scalar tol,
38  TrackingData& td
39 )
40 {
41  const scalar dist2 = magSqr(pt - w2.origin());
42 
43  if (valid(td))
44  {
45  const scalar diff = distSqr() - dist2;
46 
47  if (diff < 0)
48  {
49  // already nearer to pt
50  return false;
51  }
52 
53  if ((diff < small) || ((distSqr() > small) && (diff/distSqr() < tol)))
54  {
55  // don't propagate small changes
56  return false;
57  }
58  }
59 
60  // Either *this is not yet valid or w2 is closer
61  distSqr() = dist2;
62  origin() = w2.origin();
63 
64  return true;
65 }
66 
67 
68 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
69 
71 :
72  origin_(point::max),
73  distSqr_(-great)
74 {}
75 
76 
78 (
79  const point& origin,
80  const scalar distSqr
81 )
82 :
83  origin_(origin),
84  distSqr_(distSqr)
85 {}
86 
87 
89 (
90  const face& f,
91  const pointField& points,
92  const point& centre,
93  const scalar distSqr
94 )
95 :
96  origin_(centre),
97  distSqr_(distSqr)
98 {}
99 
100 
101 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
104 {
105  return origin_;
106 }
107 
108 
110 {
111  return origin_;
112 }
113 
114 
115 inline Foam::scalar Foam::wallPoint::distSqr() const
116 {
117  return distSqr_;
118 }
119 
120 
121 inline Foam::scalar& Foam::wallPoint::distSqr()
122 {
123  return distSqr_;
124 }
125 
126 
127 template<class TrackingData>
128 inline Foam::scalar Foam::wallPoint::dist(TrackingData& td) const
129 {
130  return valid(td) ? sqrt(distSqr_) : great;
131 }
132 
133 
134 template<class TrackingData>
135 inline bool Foam::wallPoint::valid(TrackingData& td) const
136 {
137  return distSqr_ > -small;
138 }
139 
140 
141 template<class TrackingData>
143 (
144  const wallPoint& w2,
145  const scalar tol,
146  TrackingData& td
147 ) const
148 {
149  const scalar diff = mag(distSqr() - w2.distSqr());
150 
151  return
152  diff < small
153  || ((distSqr() > small) && (diff/distSqr() < tol));
154 }
155 
156 
157 template<class TrackingData>
159 (
160  const transformer& transform,
161  TrackingData& td
162 )
163 {
164  // Note that distSqr_ is not affected by crossing an interface
165  origin_ = transform.transformPosition(origin_);
166 }
167 
168 
169 template<class TrackingData>
171 (
172  const wallPoint& rhs,
173  TrackingData& td
174 ) const
175 {
176  return operator==(rhs);
177 }
178 
179 
180 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
181 
182 inline bool Foam::wallPoint::operator==
183 (
184  const Foam::wallPoint& rhs
185 ) const
186 {
187  return origin() == rhs.origin();
188 }
189 
190 
191 inline bool Foam::wallPoint::operator!=
192 (
193  const Foam::wallPoint& rhs
194 ) const
195 {
196  return !(*this == rhs);
197 }
198 
199 
200 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
201 
203 {
204  return os << w.origin() << token::SPACE << w.distSqr();
205 }
206 
207 
209 {
210  return is >> w.origin() >> w.distSqr();
211 }
212 
213 
214 // ************************************************************************* //
#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
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:76
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: wallPoint.H:59
scalar dist(TrackingData &td) const
void transform(const transformer &transform, TrackingData &td)
Transform across an interface.
Definition: wallPointI.H:159
scalar distSqr() const
Definition: wallPointI.H:115
bool equal(const wallPoint &, TrackingData &td) const
Same (like operator==)
Definition: wallPointI.H:171
wallPoint()
Construct null.
Definition: wallPointI.H:70
bool valid(TrackingData &td) const
Check whether origin has been changed at all or.
Definition: wallPointI.H:135
bool sameGeometry(const wallPoint &, const scalar, TrackingData &td) const
Check for identical geometrical data. Used for cyclics checking.
Definition: wallPointI.H:143
const point & origin() const
Definition: wallPointI.H:103
bool update(const point &, const wallPoint &w2, const scalar tol, TrackingData &td)
Evaluate distance to point. Update distSqr, origin from whomever.
Definition: wallPointI.H:34
const pointField & points
bool valid(const PtrList< ModelType > &l)
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 > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
Ostream & operator<<(Ostream &os, const fvConstraints &constraints)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
labelList f(nPoints)