searchableCylinder.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 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 "searchableCylinder.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(searchableCylinder, 0);
38  (
39  searchableSurface,
40  searchableCylinder,
41  dict
42  );
44  (
45  searchableSurface,
46  searchableCylinder,
47  dict,
48  cylinder
49  );
50 }
51 
52 
53 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
54 
56 {
57  return tmp<pointField>::New(1, 0.5*(point1_ + point2_));
58 }
59 
60 
62 (
63  pointField& centres,
64  scalarField& radiusSqr
65 ) const
66 {
67  centres.setSize(1);
68  centres[0] = 0.5*(point1_ + point2_);
69 
70  radiusSqr.setSize(1);
71  radiusSqr[0] = Foam::magSqr(point1_-centres[0]) + Foam::sqr(radius_);
72 
73  // Add a bit to make sure all points are tested inside
74  radiusSqr += Foam::sqr(SMALL);
75 }
76 
77 
79 {
80  auto tpts = tmp<pointField>::New(2);
81  auto& pts = tpts.ref();
82 
83  pts[0] = point1_;
84  pts[1] = point2_;
85 
86  return tpts;
87 }
88 
89 
90 Foam::pointIndexHit Foam::searchableCylinder::findNearest
91 (
92  const point& sample,
93  const scalar nearestDistSqr
94 ) const
95 {
96  pointIndexHit info(false, sample, -1);
97 
98  vector v(sample - point1_);
99 
100  // Decompose sample-point1 into normal and parallel component
101  scalar parallel = (v & unitDir_);
102 
103  // Remove the parallel component and normalise
104  v -= parallel*unitDir_;
105  scalar magV = mag(v);
106 
107  if (magV < ROOTVSMALL)
108  {
109  v = Zero;
110  }
111  else
112  {
113  v /= magV;
114  }
115 
116  if (parallel <= 0)
117  {
118  // nearest is at point1 end cap. Clip to radius.
119  info.setPoint(point1_ + min(magV, radius_)*v);
120  }
121  else if (parallel >= magDir_)
122  {
123  // nearest is at point2 end cap. Clip to radius.
124  info.setPoint(point2_ + min(magV, radius_)*v);
125  }
126  else
127  {
128  // inbetween endcaps. Might either be nearer endcaps or cylinder wall
129 
130  // distance to endpoint: parallel or parallel-magDir
131  // distance to cylinder wall: magV-radius_
132 
133  // Nearest cylinder point
134  point cylPt;
135  if (magV < ROOTVSMALL)
136  {
137  // Point exactly on centre line. Take any point on wall.
138  vector e1 = point(1,0,0) ^ unitDir_;
139  scalar magE1 = mag(e1);
140  if (magE1 < SMALL)
141  {
142  e1 = point(0,1,0) ^ unitDir_;
143  magE1 = mag(e1);
144  }
145  e1 /= magE1;
146  cylPt = sample + radius_*e1;
147  }
148  else
149  {
150  cylPt = sample + (radius_-magV)*v;
151  }
152 
153  if (parallel < 0.5*magDir_)
154  {
155  // Project onto p1 endcap
156  point end1Pt = point1_ + min(magV, radius_)*v;
157 
158  if (magSqr(sample-cylPt) < magSqr(sample-end1Pt))
159  {
160  info.setPoint(cylPt);
161  }
162  else
163  {
164  info.setPoint(end1Pt);
165  }
166  }
167  else
168  {
169  // Project onto p2 endcap
170  point end2Pt = point2_ + min(magV, radius_)*v;
171 
172  if (magSqr(sample-cylPt) < magSqr(sample-end2Pt))
173  {
174  info.setPoint(cylPt);
175  }
176  else
177  {
178  info.setPoint(end2Pt);
179  }
180  }
181  }
182 
183  if (magSqr(sample - info.rawPoint()) < nearestDistSqr)
184  {
185  info.setHit();
186  info.setIndex(0);
187  }
188 
189  return info;
190 }
191 
192 
193 Foam::scalar Foam::searchableCylinder::radius2(const point& pt) const
194 {
195  const vector x = (pt-point1_) ^ unitDir_;
196  return (x & x);
197 }
198 
199 
200 // From http://www.gamedev.net/community/forums/topic.asp?topic_id=467789 -
201 // intersection of cylinder with ray
202 void Foam::searchableCylinder::findLineAll
203 (
204  const point& start,
205  const point& end,
206  pointIndexHit& near,
207  pointIndexHit& far
208 ) const
209 {
210  near.setMiss();
211  far.setMiss();
212 
213  vector point1Start(start-point1_);
214  vector point2Start(start-point2_);
215  vector point1End(end-point1_);
216 
217  // Quick rejection of complete vector outside endcaps
218  scalar s1 = point1Start&unitDir_;
219  scalar s2 = point1End&unitDir_;
220 
221  if ((s1 < 0 && s2 < 0) || (s1 > magDir_ && s2 > magDir_))
222  {
223  return;
224  }
225 
226  // Line as P = start+t*V where V is unit vector and t=[0..mag(end-start)]
227  vector V(end-start);
228  scalar magV = mag(V);
229  if (magV < ROOTVSMALL)
230  {
231  return;
232  }
233  V /= magV;
234 
235 
236  // We now get the nearest intersections to start. This can either be
237  // the intersection with the end plane or with the cylinder side.
238 
239  // Get the two points (expressed in t) on the end planes. This is to
240  // clip any cylinder intersection against.
241  scalar tPoint1;
242  scalar tPoint2;
243 
244  // Maintain the two intersections with the endcaps
245  scalar tNear = VGREAT;
246  scalar tFar = VGREAT;
247 
248  {
249  scalar s = (V&unitDir_);
250  if (mag(s) > VSMALL)
251  {
252  tPoint1 = -s1/s;
253  tPoint2 = -(point2Start&unitDir_)/s;
254  if (tPoint2 < tPoint1)
255  {
256  std::swap(tPoint1, tPoint2);
257  }
258  if (tPoint1 > magV || tPoint2 < 0)
259  {
260  return;
261  }
262 
263  // See if the points on the endcaps are actually inside the cylinder
264  if (tPoint1 >= 0 && tPoint1 <= magV)
265  {
266  if (radius2(start+tPoint1*V) <= sqr(radius_))
267  {
268  tNear = tPoint1;
269  }
270  }
271  if (tPoint2 >= 0 && tPoint2 <= magV)
272  {
273  if (radius2(start+tPoint2*V) <= sqr(radius_))
274  {
275  // Check if already have a near hit from point1
276  if (tNear <= 1)
277  {
278  tFar = tPoint2;
279  }
280  else
281  {
282  tNear = tPoint2;
283  }
284  }
285  }
286  }
287  else
288  {
289  // Vector perpendicular to cylinder. Check for outside already done
290  // above so just set tpoint to allow all.
291  tPoint1 = -VGREAT;
292  tPoint2 = VGREAT;
293  }
294  }
295 
296 
297  const vector x = point1Start ^ unitDir_;
298  const vector y = V ^ unitDir_;
299  const scalar d = sqr(radius_);
300 
301  // Second order equation of the form a*t^2 + b*t + c
302  const scalar a = (y&y);
303  const scalar b = 2*(x&y);
304  const scalar c = (x&x)-d;
305 
306  const scalar disc = b*b-4*a*c;
307 
308  scalar t1 = -VGREAT;
309  scalar t2 = VGREAT;
310 
311  if (disc < 0)
312  {
313  // Fully outside
314  return;
315  }
316  else if (disc < ROOTVSMALL)
317  {
318  // Single solution
319  if (mag(a) > ROOTVSMALL)
320  {
321  t1 = -b/(2*a);
322 
323  //Pout<< "single solution t:" << t1
324  // << " for start:" << start << " end:" << end
325  // << " c:" << c << endl;
326 
327  if (t1 >= 0 && t1 <= magV && t1 >= tPoint1 && t1 <= tPoint2)
328  {
329  // valid. Insert sorted.
330  if (t1 < tNear)
331  {
332  tFar = tNear;
333  tNear = t1;
334  }
335  else if (t1 < tFar)
336  {
337  tFar = t1;
338  }
339  }
340  else
341  {
342  return;
343  }
344  }
345  else
346  {
347  // Aligned with axis. Check if outside radius
348  //Pout<< "small discriminant:" << disc
349  // << " for start:" << start << " end:" << end
350  // << " magV:" << magV
351  // << " c:" << c << endl;
352  if (c > 0)
353  {
354  return;
355  }
356  }
357  }
358  else
359  {
360  if (mag(a) > ROOTVSMALL)
361  {
362  scalar sqrtDisc = sqrt(disc);
363 
364  t1 = (-b - sqrtDisc)/(2*a);
365  t2 = (-b + sqrtDisc)/(2*a);
366  if (t2 < t1)
367  {
368  std::swap(t1, t2);
369  }
370 
371  if (t1 >= 0 && t1 <= magV && t1 >= tPoint1 && t1 <= tPoint2)
372  {
373  // valid. Insert sorted.
374  if (t1 < tNear)
375  {
376  tFar = tNear;
377  tNear = t1;
378  }
379  else if (t1 < tFar)
380  {
381  tFar = t1;
382  }
383  }
384  if (t2 >= 0 && t2 <= magV && t2 >= tPoint1 && t2 <= tPoint2)
385  {
386  // valid. Insert sorted.
387  if (t2 < tNear)
388  {
389  tFar = tNear;
390  tNear = t2;
391  }
392  else if (t2 < tFar)
393  {
394  tFar = t2;
395  }
396  }
397  //Pout<< "two solutions t1:" << t1 << " t2:" << t2
398  // << " for start:" << start << " end:" << end
399  // << " magV:" << magV
400  // << " c:" << c << endl;
401  }
402  else
403  {
404  // Aligned with axis. Check if outside radius
405  //Pout<< "large discriminant:" << disc
406  // << " small a:" << a
407  // << " for start:" << start << " end:" << end
408  // << " magV:" << magV
409  // << " c:" << c << endl;
410  if (c > 0)
411  {
412  return;
413  }
414  }
415  }
416 
417  // Check tNear, tFar
418  if (tNear >= 0 && tNear <= magV)
419  {
420  near.setPoint(start+tNear*V);
421  near.setHit();
422  near.setIndex(0);
423 
424  if (tFar <= magV)
425  {
426  far.setPoint(start+tFar*V);
427  far.setHit();
428  far.setIndex(0);
429  }
430  }
431  else if (tFar >= 0 && tFar <= magV)
432  {
433  near.setPoint(start+tFar*V);
434  near.setHit();
435  near.setIndex(0);
436  }
437 }
438 
439 
440 Foam::boundBox Foam::searchableCylinder::calcBounds() const
441 {
442 
443  // Adapted from
444  // http://www.gamedev.net/community/forums
445  // /topic.asp?topic_id=338522&forum_id=20&gforum_id=0
446 
447  // Let cylinder have end points A,B and radius r,
448 
449  // Bounds in direction X (same for Y and Z) can be found as:
450  // Let A.X<B.X (otherwise swap points)
451  // Good approximate lowest bound is A.X-r and highest is B.X+r (precise for
452  // capsule). At worst, in one direction it can be larger than needed by 2*r.
453 
454  // Accurate bounds for cylinder is
455  // A.X-kx*r, B.X+kx*r
456  // where
457  // kx=sqrt(((A.Y-B.Y)^2+(A.Z-B.Z)^2)/((A.X-B.X)^2+(A.Y-B.Y)^2+(A.Z-B.Z)^2))
458 
459  // similar thing for Y and Z
460  // (i.e.
461  // ky=sqrt(((A.X-B.X)^2+(A.Z-B.Z)^2)/((A.X-B.X)^2+(A.Y-B.Y)^2+(A.Z-B.Z)^2))
462  // kz=sqrt(((A.X-B.X)^2+(A.Y-B.Y)^2)/((A.X-B.X)^2+(A.Y-B.Y)^2+(A.Z-B.Z)^2))
463  // )
464 
465  // How derived: geometric reasoning. Bounds of cylinder is same as for 2
466  // circles centered on A and B. This sqrt thingy gives sine of angle between
467  // axis and direction, used to find projection of radius.
468 
469  vector kr
470  (
471  sqrt(sqr(unitDir_.y()) + sqr(unitDir_.z())),
472  sqrt(sqr(unitDir_.x()) + sqr(unitDir_.z())),
473  sqrt(sqr(unitDir_.x()) + sqr(unitDir_.y()))
474  );
475 
476  kr *= radius_;
477 
478  point min = point1_ - kr;
479  point max = point1_ + kr;
480 
481  min = ::Foam::min(min, point2_ - kr);
482  max = ::Foam::max(max, point2_ + kr);
483 
484  return boundBox(min, max);
485 }
486 
487 
488 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
489 
490 Foam::searchableCylinder::searchableCylinder
491 (
492  const IOobject& io,
493  const point& point1,
494  const point& point2,
495  const scalar radius
496 )
497 :
498  searchableSurface(io),
499  point1_(point1),
500  point2_(point2),
501  magDir_(mag(point2_-point1_)),
502  unitDir_((point2_-point1_)/magDir_),
503  radius_(radius)
504 {
505  bounds() = calcBounds();
506 }
507 
508 
509 Foam::searchableCylinder::searchableCylinder
510 (
511  const IOobject& io,
512  const dictionary& dict
513 )
514 :
515  searchableSurface(io),
516  point1_(dict.get<point>("point1")),
517  point2_(dict.get<point>("point2")),
518  magDir_(mag(point2_-point1_)),
519  unitDir_((point2_-point1_)/magDir_),
520  radius_(dict.get<scalar>("radius"))
521 {
522  bounds() = calcBounds();
523 }
524 
525 
526 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
527 
529 {
530  if (regions_.empty())
531  {
532  regions_.setSize(1);
533  regions_[0] = "region0";
534  }
535  return regions_;
536 }
537 
538 
539 void Foam::searchableCylinder::findNearest
540 (
541  const pointField& samples,
542  const scalarField& nearestDistSqr,
543  List<pointIndexHit>& info
544 ) const
545 {
546  info.setSize(samples.size());
547 
548  forAll(samples, i)
549  {
550  info[i] = findNearest(samples[i], nearestDistSqr[i]);
551  }
552 }
553 
554 
556 (
557  const pointField& start,
558  const pointField& end,
559  List<pointIndexHit>& info
560 ) const
561 {
562  info.setSize(start.size());
563 
564  forAll(start, i)
565  {
566  // Pick nearest intersection. If none intersected take second one.
568  findLineAll(start[i], end[i], info[i], b);
569  if (!info[i].hit() && b.hit())
570  {
571  info[i] = b;
572  }
573  }
574 }
575 
576 
578 (
579  const pointField& start,
580  const pointField& end,
581  List<pointIndexHit>& info
582 ) const
583 {
584  info.setSize(start.size());
585 
586  forAll(start, i)
587  {
588  // Discard far intersection
590  findLineAll(start[i], end[i], info[i], b);
591  if (!info[i].hit() && b.hit())
592  {
593  info[i] = b;
594  }
595  }
596 }
597 
598 
599 void Foam::searchableCylinder::findLineAll
600 (
601  const pointField& start,
602  const pointField& end,
604 ) const
605 {
606  info.setSize(start.size());
607 
608  forAll(start, i)
609  {
610  pointIndexHit near, far;
611  findLineAll(start[i], end[i], near, far);
612 
613  if (near.hit())
614  {
615  if (far.hit())
616  {
617  info[i].setSize(2);
618  info[i][0] = near;
619  info[i][1] = far;
620  }
621  else
622  {
623  info[i].setSize(1);
624  info[i][0] = near;
625  }
626  }
627  else
628  {
629  if (far.hit())
630  {
631  info[i].setSize(1);
632  info[i][0] = far;
633  }
634  else
635  {
636  info[i].clear();
637  }
638  }
639  }
640 }
641 
642 
644 (
645  const List<pointIndexHit>& info,
646  labelList& region
647 ) const
648 {
649  region.setSize(info.size());
650  region = 0;
651 }
652 
653 
655 (
656  const List<pointIndexHit>& info,
657  vectorField& normal
658 ) const
659 {
660  normal.setSize(info.size());
661  normal = Zero;
662 
663  forAll(info, i)
664  {
665  if (info[i].hit())
666  {
667  vector v(info[i].hitPoint() - point1_);
668 
669  // Decompose sample-point1 into normal and parallel component
670  const scalar parallel = (v & unitDir_);
671 
672  // Remove the parallel component and normalise
673  v -= parallel*unitDir_;
674  scalar magV = mag(v);
675 
676  if (parallel <= 0)
677  {
678  if ((magV-radius_) < mag(parallel))
679  {
680  // either above endcap (magV<radius) or outside but closer
681  normal[i] = -unitDir_;
682  }
683  else
684  {
685  normal[i] = v/magV;
686  }
687  }
688  else if (parallel <= 0.5*magDir_)
689  {
690  // See if endcap closer or sidewall
691  if (magV >= radius_ || (radius_-magV) < parallel)
692  {
693  normal[i] = v/magV;
694  }
695  else
696  {
697  // closer to endcap
698  normal[i] = -unitDir_;
699  }
700  }
701  else if (parallel <= magDir_)
702  {
703  // See if endcap closer or sidewall
704  if (magV >= radius_ || (radius_-magV) < (magDir_-parallel))
705  {
706  normal[i] = v/magV;
707  }
708  else
709  {
710  // closer to endcap
711  normal[i] = unitDir_;
712  }
713  }
714  else // beyond cylinder
715  {
716  if ((magV-radius_) < (parallel-magDir_))
717  {
718  // above endcap
719  normal[i] = unitDir_;
720  }
721  else
722  {
723  normal[i] = v/magV;
724  }
725  }
726  }
727  }
728 }
729 
730 
732 (
733  const pointField& points,
734  List<volumeType>& volType
735 ) const
736 {
737  volType.setSize(points.size());
738 
739  forAll(points, pointi)
740  {
741  const point& pt = points[pointi];
742 
743  volType[pointi] = volumeType::OUTSIDE;
744 
745  vector v(pt - point1_);
746 
747  // Decompose sample-point1 into normal and parallel component
748  const scalar parallel = (v & unitDir_);
749 
750  // Quick rejection. Left of point1 endcap, or right of point2 endcap
751  if (parallel < 0 || parallel > magDir_)
752  {
753  volType[pointi] = volumeType::OUTSIDE;
754  }
755  else
756  {
757  // Remove the parallel component
758  v -= parallel*unitDir_;
759 
760  volType[pointi] =
761  (
762  mag(v) <= radius_
765  );
766  }
767  }
768 }
769 
770 
771 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Foam::searchableCylinder::boundingSpheres
virtual void boundingSpheres(pointField &centres, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
Definition: searchableCylinder.C:62
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::searchableCylinder::regions
virtual const wordList & regions() const
Names of regions.
Definition: searchableCylinder.C:528
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::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::searchableCylinder::findLine
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
Definition: searchableCylinder.C:556
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
Foam::searchableCylinder::coordinates
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
Definition: searchableCylinder.C:55
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
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::searchableCylinder::getRegion
virtual void getRegion(const List< pointIndexHit > &, labelList &region) const
From a set of points and indices get the region.
Definition: searchableCylinder.C:644
Foam::PointIndexHit
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:52
Foam::PointIndexHit::hit
bool hit() const noexcept
Is there a hit?
Definition: PointIndexHit.H:130
searchableCylinder.H
Foam::constant::physicoChemical::b
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
Foam::Field< vector >
Foam::List::setSize
void setSize(const label n)
Alias for resize()
Definition: List.H:222
Foam::searchableSurface
Base class of (analytical or triangulated) surface. Encapsulates all the search routines....
Definition: searchableSurface.H:69
Foam::addNamedToRunTimeSelectionTable
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
samples
scalarField samples(nIntervals, Zero)
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
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
stdFoam::end
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
Definition: stdFoam.H:121
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::searchableCylinder::points
virtual tmp< pointField > points() const
Get the points that define the surface.
Definition: searchableCylinder.C:78
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::pointIndexHit
PointIndexHit< point > pointIndexHit
A PointIndexHit for 3D points.
Definition: pointIndexHit.H:46
Foam::Vector< scalar >
Foam::List< word >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::searchableCylinder::findLineAny
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
Definition: searchableCylinder.C:578
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::List::clear
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:116
Foam::tmp::New
static tmp< T > New(Args &&... args)
Construct tmp of T with forwarding arguments.
x
x
Definition: LISASMDCalcMethod2.H:52
Foam::volumeType::INSIDE
A location inside the volume.
Definition: volumeType.H:68
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::searchableCylinder::getNormal
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
Definition: searchableCylinder.C:655
Foam::volumeType::OUTSIDE
A location outside the volume.
Definition: volumeType.H:69
Foam::searchableCylinder::getVolumeType
virtual void getVolumeType(const pointField &points, List< volumeType > &volType) const
Determine type (inside/outside/mixed) for point.
Definition: searchableCylinder.C:732
y
scalar y
Definition: LISASMDCalcMethod1.H:14
sample
Minimal example by using system/controlDict.functions: