lineI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011 OpenFOAM Foundation
9  Copyright (C) 2018-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "zero.H"
30 #include "IOstreams.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class Point, class PointRef>
35 inline Foam::line<Point, PointRef>::line(const Point& from, const Point& to)
36 :
37  a_(from),
38  b_(to)
39 {}
40 
41 
42 template<class Point, class PointRef>
44 (
45  const UList<Point>& points,
46  const FixedList<label, 2>& indices
47 )
48 :
49  a_(points[indices[0]]),
50  b_(points[indices[1]])
51 {}
52 
53 
54 template<class Point, class PointRef>
56 {
57  is >> *this;
58 }
59 
60 
61 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
62 
63 template<class Point, class PointRef>
64 inline PointRef Foam::line<Point, PointRef>::first() const noexcept
65 {
66  return a_;
67 }
68 
69 
70 template<class Point, class PointRef>
71 inline PointRef Foam::line<Point, PointRef>::second() const noexcept
72 {
73  return b_;
74 }
75 
76 
77 template<class Point, class PointRef>
78 inline PointRef Foam::line<Point, PointRef>::last() const noexcept
79 {
80  return b_;
81 }
82 
83 
84 template<class Point, class PointRef>
85 inline PointRef Foam::line<Point, PointRef>::start() const noexcept
86 {
87  return first();
88 }
89 
90 template<class Point, class PointRef>
91 inline PointRef Foam::line<Point, PointRef>::end() const noexcept
92 {
93  return second();
94 }
95 
96 
97 template<class Point, class PointRef>
99 {
100  return 0.5*(a_ + b_);
101 }
102 
103 
104 template<class Point, class PointRef>
105 inline Foam::scalar Foam::line<Point, PointRef>::mag() const
106 {
107  return ::Foam::mag(vec());
108 }
109 
110 
111 template<class Point, class PointRef>
113 {
114  return b_ - a_;
115 }
116 
117 
118 template<class Point, class PointRef>
120 {
121  const Point v = (b_ - a_);
122  const scalar s(::Foam::mag(v));
123 
124  return s < ROOTVSMALL ? Zero : v/s;
125 }
126 
127 
128 template<class Point, class PointRef>
130 (
131  const Point& p
132 ) const
133 {
134  Point v = vec();
135 
136  Point w(p - a_);
137 
138  const scalar c1 = v & w;
139 
140  if (c1 <= 0)
141  {
142  return PointHit<Point>(false, a_, Foam::mag(p - a_), true);
143  }
144 
145  const scalar c2 = v & v;
146 
147  if (c2 <= c1)
148  {
149  return PointHit<Point>(false, b_, Foam::mag(p - b_), true);
150  }
151 
152  const scalar b = c1/c2;
153 
154  Point pb(a_ + b*v);
155 
156  return PointHit<Point>(true, pb, Foam::mag(p - pb), false);
157 }
158 
159 
160 template<class Point, class PointRef>
162 (
164  Point& thisPt,
165  Point& edgePt
166 ) const
167 {
168  // From Mathworld Line-Line distance/(Gellert et al. 1989, p. 538).
169  Point a(end() - start());
170  Point b(edge.end() - edge.start());
171  Point c(edge.start() - start());
172 
173  Point crossab = a ^ b;
174  const scalar magCrossSqr = Foam::magSqr(crossab);
175 
176  if (magCrossSqr > VSMALL)
177  {
178  scalar s = ((c ^ b) & crossab)/magCrossSqr;
179  scalar t = ((c ^ a) & crossab)/magCrossSqr;
180 
181  // Check for end points outside of range 0..1
182  if (s >= 0 && s <= 1 && t >= 0 && t <= 1)
183  {
184  // Both inside range 0..1
185  thisPt = start() + a*s;
186  edgePt = edge.start() + b*t;
187  }
188  else
189  {
190  // Do brute force. Distance of everything to everything.
191  // Can quite possibly be improved!
192 
193  // From edge endpoints to *this
194  PointHit<Point> this0(nearestDist(edge.start()));
195  PointHit<Point> this1(nearestDist(edge.end()));
196  scalar thisDist = min(this0.distance(), this1.distance());
197 
198  // From *this to edge
199  PointHit<Point> edge0(edge.nearestDist(start()));
200  PointHit<Point> edge1(edge.nearestDist(end()));
201  scalar edgeDist = min(edge0.distance(), edge1.distance());
202 
203  if (thisDist < edgeDist)
204  {
205  if (this0.distance() < this1.distance())
206  {
207  thisPt = this0.rawPoint();
208  edgePt = edge.start();
209  }
210  else
211  {
212  thisPt = this1.rawPoint();
213  edgePt = edge.end();
214  }
215  }
216  else
217  {
218  if (edge0.distance() < edge1.distance())
219  {
220  thisPt = start();
221  edgePt = edge0.rawPoint();
222  }
223  else
224  {
225  thisPt = end();
226  edgePt = edge1.rawPoint();
227  }
228  }
229  }
230  }
231  else
232  {
233  // Parallel lines. Find overlap of both lines by projecting onto
234  // direction vector (now equal for both lines).
235 
236  scalar edge0 = edge.start() & a;
237  scalar edge1 = edge.end() & a;
238  bool edgeOrder = edge0 < edge1;
239 
240  scalar minEdge = (edgeOrder ? edge0 : edge1);
241  scalar maxEdge = (edgeOrder ? edge1 : edge0);
242  const Point& minEdgePt = (edgeOrder ? edge.start() : edge.end());
243  const Point& maxEdgePt = (edgeOrder ? edge.end() : edge.start());
244 
245  scalar this0 = start() & a;
246  scalar this1 = end() & a;
247  bool thisOrder = this0 < this1;
248 
249  scalar minThis = min(this0, this1);
250  scalar maxThis = max(this1, this0);
251  const Point& minThisPt = (thisOrder ? start() : end());
252  const Point& maxThisPt = (thisOrder ? end() : start());
253 
254  if (maxEdge < minThis)
255  {
256  // edge completely below *this
257  edgePt = maxEdgePt;
258  thisPt = minThisPt;
259  }
260  else if (maxEdge < maxThis)
261  {
262  // maxEdge inside interval of *this
263  edgePt = maxEdgePt;
264  thisPt = nearestDist(edgePt).rawPoint();
265  }
266  else
267  {
268  // maxEdge outside. Check if minEdge inside.
269  if (minEdge < minThis)
270  {
271  // Edge completely envelops this. Take any this point and
272  // determine nearest on edge.
273  thisPt = minThisPt;
274  edgePt = edge.nearestDist(thisPt).rawPoint();
275  }
276  else if (minEdge < maxThis)
277  {
278  // minEdge inside this interval.
279  edgePt = minEdgePt;
280  thisPt = nearestDist(edgePt).rawPoint();
281  }
282  else
283  {
284  // minEdge outside this interval
285  edgePt = minEdgePt;
286  thisPt = maxThisPt;
287  }
288  }
289  }
290 
291  return Foam::mag(thisPt - edgePt);
292 }
293 
294 
295 // * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * //
296 
297 template<class Point, class PointRef>
298 inline Foam::Istream& Foam::operator>>
299 (
300  Istream& is,
302 )
303 {
304  is.readBegin("line");
305  is >> l.a_ >> l.b_;
306  is.readEnd("line");
307 
308  is.check(FUNCTION_NAME);
309  return is;
310 }
311 
312 
313 template<class Point, class PointRef>
314 inline Foam::Ostream& Foam::operator<<
315 (
316  Ostream& os,
317  const line<Point, PointRef>& l
318 )
319 {
321  << l.a_ << token::SPACE
322  << l.b_
323  << token::END_LIST;
324  return os;
325 }
326 
327 
328 // ************************************************************************* //
Foam::line::line
line(const Point &from, const Point &to)
Construct from two points.
Definition: lineI.H:35
Foam::line::centre
Point centre() const
Return centre (centroid)
Definition: lineI.H:98
p
volScalarField & p
Definition: createFieldRefs.H:8
IOstreams.H
Useful combination of include files which define Sin, Sout and Serr and the use of IO streams general...
s
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputSpray.H:25
Foam::PointHit
Describes the interaction of a face and a point. It carries the info of a successful hit and (if succ...
Definition: PointHit.H:53
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::edge
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:63
Foam::line::start
PointRef start() const noexcept
Return first point.
Definition: lineI.H:85
Foam::PointHit::distance
scalar distance() const noexcept
Return distance to hit.
Definition: PointHit.H:139
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::Istream::readBegin
bool readBegin(const char *funcName)
Begin read of data chunk, starts with '('.
Definition: Istream.C:111
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::edge::start
label start() const
Return start (first) vertex label.
Definition: edgeI.H:95
Foam::constant::physicoChemical::c1
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
Foam::constant::physicoChemical::b
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Foam::Istream
An Istream is an abstract base class for all input systems (streams, files, token lists etc)....
Definition: Istream.H:61
Foam::PointHit::rawPoint
const point_type & rawPoint() const noexcept
The point, no checks.
Definition: PointHit.H:172
Foam::line::end
PointRef end() const noexcept
Return second (last) point.
Definition: lineI.H:91
Foam::edge::end
label end() const
Return end (last/second) vertex label.
Definition: edgeI.H:106
Foam::line::second
PointRef second() const noexcept
Return second (last) point.
Definition: lineI.H:71
zero.H
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::line::mag
scalar mag() const
Return scalar magnitude.
Definition: lineI.H:105
Foam::line::unitVec
Point unitVec() const
Return the unit vector (start-to-end)
Definition: lineI.H:119
Foam::line::last
PointRef last() const noexcept
Return last (second) point.
Definition: lineI.H:78
os
OBJstream os(runTime.globalPath()/outputName)
stdFoam::end
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Definition: stdFoam.H:121
Foam::constant::physicoChemical::c2
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::FixedList
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: HashTable.H:104
Foam::token::SPACE
Space [isspace].
Definition: token.H:125
Foam::UList
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:103
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::line
A line primitive.
Definition: line.H:53
FUNCTION_NAME
#define FUNCTION_NAME
Definition: messageStream.H:295
Foam::token::END_LIST
End list [isseparator].
Definition: token.H:156
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Point
CGAL::Point_3< K > Point
Definition: CGALIndexedPolyhedron.H:53
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::token::BEGIN_LIST
Begin list [isseparator].
Definition: token.H:155
Foam::line::nearestDist
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:130
Foam::line::vec
Point vec() const
Return start-to-end vector.
Definition: lineI.H:112
Foam::line::first
PointRef first() const noexcept
Return first point.
Definition: lineI.H:64