searchableSphere.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-2022 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 Class
25  Foam::searchableSphere
26 
27 Description
28  Surface geometry with a sphere shape, which can be used with
29  snappyHexMesh.
30 
31 Usage
32  \table
33  Property | Description | Required
34  centre | Sphere centre | yes
35  radius | Sphere radius | yes
36  \endtable
37 
38  Example specification in snappyHexMeshDict/geometry:
39  \verbatim
40  type searchableSphere;
41  centre (10 10 10);
42  radius 5;
43  \endverbatim
44 
45 SourceFiles
46  searchableSphere.C
47 
48 \*---------------------------------------------------------------------------*/
49 
50 #ifndef searchableSphere_H
51 #define searchableSphere_H
52 
53 #include "treeBoundBox.H"
54 #include "searchableSurface.H"
55 
56 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57 
58 namespace Foam
59 {
60 
61 // Forward declaration of classes
62 
63 /*---------------------------------------------------------------------------*\
64  Class searchableSphere Declaration
65 \*---------------------------------------------------------------------------*/
66 
67 class searchableSphere
68 :
69  public searchableSurface
70 {
71  // Private Member Data
72 
73  //- Centre point
74  const point centre_;
75 
76  //- Radius
77  const scalar radius_;
78 
79  //- Names of regions
80  mutable wordList regions_;
81 
82 
83  // Private Member Functions
84 
85  //- Inherit findNearest from searchableSurface
87 
88  //- Find nearest point on sphere.
89  pointIndexHit findNearest
90  (
91  const point& sample,
92  const scalar nearestDistSqr
93  ) const;
94 
95  //- Find intersection with sphere
96  void findLineAll
97  (
98  const point& start,
99  const point& end,
100  pointIndexHit& near,
101  pointIndexHit& far
102  ) const;
103 
104 
105 public:
106 
107  //- Runtime type information
108  TypeName("searchableSphere");
109 
110 
111  // Constructors
112 
113  //- Construct from components
114  searchableSphere(const IOobject& io, const point&, const scalar radius);
115 
116  //- Construct from dictionary (used by searchableSurface)
118  (
119  const IOobject& io,
120  const dictionary& dict
121  );
122 
123  //- Disallow default bitwise copy construction
124  searchableSphere(const searchableSphere&) = delete;
125 
126 
127  //- Destructor
128  virtual ~searchableSphere();
129 
130 
131  // Member Functions
132 
133  virtual const wordList& regions() const;
134 
135  //- Whether supports volume type below
136  virtual bool hasVolumeType() const
137  {
138  return true;
139  }
140 
141  //- Range of local indices that can be returned.
142  virtual label size() const
143  {
144  return 1;
145  }
146 
147  //- Get representative set of element coordinates
148  // Usually the element centres (should be of length size()).
149  virtual tmp<pointField> coordinates() const
150  {
151  tmp<pointField> tCtrs(new pointField(1, centre_));
152  return tCtrs;
153  }
154 
155  //- Get bounding spheres (centre and radius squared), one per element.
156  // Any point on element is guaranteed to be inside.
157  virtual void boundingSpheres
158  (
159  pointField& centres,
160  scalarField& radiusSqr
161  ) const;
162 
163  //- Get the points that define the surface.
164  virtual tmp<pointField> points() const
165  {
166  return coordinates();
167  }
168 
169  //- Does any part of the surface overlap the supplied bound box?
170  virtual bool overlaps(const boundBox& bb) const;
171 
172 
173  // Multiple point queries.
174 
175  virtual void findNearest
176  (
177  const pointField& sample,
178  const scalarField& nearestDistSqr,
180  ) const;
181 
182  virtual void findLine
183  (
184  const pointField& start,
185  const pointField& end,
187  ) const;
188 
189  virtual void findLineAny
190  (
191  const pointField& start,
192  const pointField& end,
194  ) const;
195 
196  //- Get all intersections in order from start to end.
197  virtual void findLineAll
198  (
199  const pointField& start,
200  const pointField& end,
202  ) const;
203 
204  //- From a set of points and indices get the region
205  virtual void getRegion
206  (
207  const List<pointIndexHit>&,
208  labelList& region
209  ) const;
210 
211  //- From a set of points and indices get the normal
212  virtual void getNormal
213  (
214  const List<pointIndexHit>&,
215  vectorField& normal
216  ) const;
217 
218  //- Determine type (inside/outside/mixed) for point. unknown if
219  // cannot be determined (e.g. non-manifold surface)
220  virtual void getVolumeType
221  (
222  const pointField&,
224  ) const;
225 
226 
227  // regIOobject implementation
228 
229  bool writeData(Ostream&) const
230  {
232  return false;
233  }
234 
235 
236  // Member Operators
237 
238  //- Disallow default bitwise assignment
239  void operator=(const searchableSphere&) = delete;
240 };
241 
242 
243 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
244 
245 } // End namespace Foam
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 #endif
250 
251 // ************************************************************************* //
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:54
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Surface geometry with a sphere shape, which can be used with snappyHexMesh.
virtual label size() const
Range of local indices that can be returned.
TypeName("searchableSphere")
Runtime type information.
searchableSphere(const IOobject &io, const point &, const scalar radius)
Construct from components.
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
bool writeData(Ostream &) const
Pure virtual writaData function.
virtual void findLine(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Find first intersection on segment from start to end.
virtual bool overlaps(const boundBox &bb) const
Does any part of the surface overlap the supplied bound box?
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
virtual void boundingSpheres(pointField &centres, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
virtual tmp< pointField > points() const
Get the points that define the surface.
virtual ~searchableSphere()
Destructor.
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
void operator=(const searchableSphere &)=delete
Disallow default bitwise assignment.
virtual void getRegion(const List< pointIndexHit > &, labelList &region) const
From a set of points and indices get the region.
virtual const wordList & regions() const
Names of regions.
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
virtual bool hasVolumeType() const
Whether supports volume type below.
virtual void findNearest(const pointField &sample, const scalarField &nearestDistSqr, List< pointIndexHit > &) const =0
A class for managing temporary objects.
Definition: tmp.H:55
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
Namespace for OpenFOAM.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
vectorField pointField
pointField is a vectorField.
Definition: pointFieldFwd.H:42
vector point
Point is a vector.
Definition: point.H:41
dictionary dict