line.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011 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 "line.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 
33 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
34 
35 template<>
37 (
39  point2D& thisPt,
40  point2D& edgePt
41 ) const
42 {
43  vector2D u = end()-start();
44  vector2D v = e.end()-e.start();
45  vector2D w = start()-e.start();
46 
47  scalar d = u.perp(v);
48 
49  if (Foam::mag(d) > VSMALL)
50  {
51  scalar s = v.perp(w) / d;
52 
53  if (s <= SMALL)
54  {
55  thisPt = start();
56  }
57  else if (s >= (1-SMALL))
58  {
59  thisPt = end();
60  }
61  else
62  {
63  thisPt = start()+s*u;
64  }
65 
66 
67  scalar t = u.perp(w) / d;
68 
69  if (t <= SMALL)
70  {
71  edgePt = e.start();
72  }
73  else if (t >= (1-SMALL))
74  {
75  edgePt = e.end();
76  }
77  else
78  {
79  edgePt = e.start()+t*v;
80  }
81  }
82  else
83  {
84  // Parallel lines. Find overlap of both lines by projecting onto
85  // direction vector (now equal for both lines).
86 
87  scalar edge0 = e.start() & u;
88  scalar edge1 = e.end() & u;
89  bool edgeOrder = edge0 < edge1;
90 
91  scalar minEdge = (edgeOrder ? edge0 : edge1);
92  scalar maxEdge = (edgeOrder ? edge1 : edge0);
93  const point2D& minEdgePt = (edgeOrder ? e.start() : e.end());
94  const point2D& maxEdgePt = (edgeOrder ? e.end() : e.start());
95 
96  scalar this0 = start() & u;
97  scalar this1 = end() & u;
98  bool thisOrder = this0 < this1;
99 
100  scalar minThis = min(this0, this1);
101  scalar maxThis = max(this1, this0);
102  const point2D& minThisPt = (thisOrder ? start() : end());
103  const point2D& maxThisPt = (thisOrder ? end() : start());
104 
105  if (maxEdge < minThis)
106  {
107  // edge completely below *this
108  edgePt = maxEdgePt;
109  thisPt = minThisPt;
110  }
111  else if (maxEdge < maxThis)
112  {
113  // maxEdge inside interval of *this
114  edgePt = maxEdgePt;
115  thisPt = nearestDist(edgePt).rawPoint();
116  }
117  else
118  {
119  // maxEdge outside. Check if minEdge inside.
120  if (minEdge < minThis)
121  {
122  // Edge completely envelops this. Take any this point and
123  // determine nearest on edge.
124  thisPt = minThisPt;
125  edgePt = e.nearestDist(thisPt).rawPoint();
126  }
127  else if (minEdge < maxThis)
128  {
129  // minEdge inside this interval.
130  edgePt = minEdgePt;
131  thisPt = nearestDist(edgePt).rawPoint();
132  }
133  else
134  {
135  // minEdge outside this interval
136  edgePt = minEdgePt;
137  thisPt = maxThisPt;
138  }
139  }
140  }
141 
142  return Foam::mag(thisPt - edgePt);
143 }
144 
145 
146 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
147 
148 } // End namespace Foam
149 
150 // ************************************************************************* //
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Vector2D< scalar > vector2D
vector2D obtained from generic Vector2D
Definition: vector2D.H:49
const double e
Definition: doubleFloat.H:78
dimensioned< scalar > mag(const dimensioned< Type > &)
vector2D point2D
Point2D is a vector.
Definition: point2D.H:41
PointRef start() const
Return first vertex.
Definition: lineI.H:60
volScalarField::GeometricBoundaryField d(nearWallDist(mesh).y())
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
PointHit< Point > nearestDist(const Point &p) const
Return nearest distance to line from a given point.
Definition: lineI.H:95
PointRef end() const
Return second vertex.
Definition: lineI.H:66
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().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))
A line primitive.
Definition: line.H:57