searchablePlate.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) 2020 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 "searchablePlate.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36  defineTypeNameAndDebug(searchablePlate, 0);
38  (
39  searchableSurface,
40  searchablePlate,
41  dict
42  );
44  (
45  searchableSurface,
46  searchablePlate,
47  dict,
48  plate
49  );
50 }
51 
52 
53 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
54 
55 Foam::direction Foam::searchablePlate::calcNormal(const point& span)
56 {
57  direction normalDir = 3;
58 
59  for (direction dir = 0; dir < vector::nComponents; ++dir)
60  {
61  if (span[dir] < 0)
62  {
63  // Negative entry. Flag and exit.
64  normalDir = 3;
65  break;
66  }
67  else if (span[dir] < VSMALL)
68  {
69  if (normalDir != 3)
70  {
71  // Multiple zero entries. Flag and exit.
72  normalDir = 3;
73  break;
74  }
75 
76  normalDir = dir;
77  }
78  }
79 
80  if (normalDir == 3)
81  {
83  << "Span should have two positive and one zero entry: "
84  << span << nl
85  << exit(FatalError);
86  }
87 
88  return normalDir;
89 }
90 
91 
92 // Returns miss or hit with face (always 0)
93 Foam::pointIndexHit Foam::searchablePlate::findNearest
94 (
95  const point& sample,
96  const scalar nearestDistSqr
97 ) const
98 {
99  // For every component direction can be
100  // left of min, right of max or inbetween.
101  // - outside points: project first one x plane (either min().x()
102  // or max().x()), then onto y plane and finally z. You should be left
103  // with intersection point
104  // - inside point: find nearest side (compare to mid point). Project onto
105  // that.
106 
107  // Project point on plane.
108  pointIndexHit info(true, sample, 0);
109  info.rawPoint()[normalDir_] = origin_[normalDir_];
110 
111  // Clip to edges if outside
112  for (direction dir = 0; dir < vector::nComponents; ++dir)
113  {
114  if (dir != normalDir_)
115  {
116  if (info.rawPoint()[dir] < origin_[dir])
117  {
118  info.rawPoint()[dir] = origin_[dir];
119  }
120  else if (info.rawPoint()[dir] > origin_[dir]+span_[dir])
121  {
122  info.rawPoint()[dir] = origin_[dir]+span_[dir];
123  }
124  }
125  }
126 
127  // Check if outside. Optimisation: could do some checks on distance already
128  // on components above
129  if (magSqr(info.rawPoint() - sample) > nearestDistSqr)
130  {
131  info.setMiss();
132  info.setIndex(-1);
133  }
134 
135  return info;
136 }
137 
138 
139 Foam::pointIndexHit Foam::searchablePlate::findLine
140 (
141  const point& start,
142  const point& end
143 ) const
144 {
145  pointIndexHit info
146  (
147  true,
148  Zero,
149  0
150  );
151 
152  const vector dir(end-start);
153 
154  if (mag(dir[normalDir_]) < VSMALL)
155  {
156  info.setMiss();
157  info.setIndex(-1);
158  }
159  else
160  {
161  scalar t = (origin_[normalDir_]-start[normalDir_]) / dir[normalDir_];
162 
163  if (t < 0 || t > 1)
164  {
165  info.setMiss();
166  info.setIndex(-1);
167  }
168  else
169  {
170  info.rawPoint() = start+t*dir;
171  info.rawPoint()[normalDir_] = origin_[normalDir_];
172 
173  // Clip to edges
174  for (direction dir = 0; dir < vector::nComponents; ++dir)
175  {
176  if (dir != normalDir_)
177  {
178  if (info.rawPoint()[dir] < origin_[dir])
179  {
180  info.setMiss();
181  info.setIndex(-1);
182  break;
183  }
184  else if (info.rawPoint()[dir] > origin_[dir]+span_[dir])
185  {
186  info.setMiss();
187  info.setIndex(-1);
188  break;
189  }
190  }
191  }
192  }
193  }
194 
195  // Debug
196  if (info.hit())
197  {
198  treeBoundBox bb(origin_, origin_+span_);
199  bb.min()[normalDir_] -= 1e-6;
200  bb.max()[normalDir_] += 1e-6;
201 
202  if (!bb.contains(info.hitPoint()))
203  {
205  << "bb:" << bb << endl
206  << "origin_:" << origin_ << endl
207  << "span_:" << span_ << endl
208  << "normalDir_:" << normalDir_ << endl
209  << "hitPoint:" << info.hitPoint()
210  << abort(FatalError);
211  }
212  }
213 
214  return info;
215 }
216 
217 
218 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
219 
220 Foam::searchablePlate::searchablePlate
221 (
222  const IOobject& io,
223  const point& origin,
224  const vector& span
225 )
226 :
227  searchableSurface(io),
228  origin_(origin),
229  span_(span),
230  normalDir_(calcNormal(span_))
231 {
233  << " origin:" << origin_
234  << " origin+span:" << origin_+span_
235  << " normal:" << vector::componentNames[normalDir_] << nl;
236 
237  bounds() = boundBox(origin_, origin_ + span_);
238 }
239 
240 
241 Foam::searchablePlate::searchablePlate
242 (
243  const IOobject& io,
244  const dictionary& dict
245 )
246 :
247  searchablePlate(io, dict.get<point>("origin"), dict.get<vector>("span"))
248 {}
249 
250 
251 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
252 
254 {
255  if (regions_.empty())
256  {
257  regions_.resize(1);
258  regions_.first() = "region0";
259  }
260  return regions_;
261 }
262 
263 
265 {
266  return tmp<pointField>::New(1, origin_ + 0.5*span_);
267 }
268 
269 
271 (
272  pointField& centres,
273  scalarField& radiusSqr
274 ) const
275 {
276  centres.resize(1);
277  radiusSqr.resize(1);
278 
279  centres[0] = origin_ + 0.5*span_;
280  radiusSqr[0] = Foam::magSqr(0.5*span_);
281 
282  // Add a bit to make sure all points are tested inside
283  radiusSqr += Foam::sqr(SMALL);
284 }
285 
286 
288 {
289  auto tpts = tmp<pointField>::New(4, origin_);
290  auto& pts = tpts.ref();
291 
292  pts[2] += span_;
293 
294  if (span_.x() < span_.y() && span_.x() < span_.z())
295  {
296  pts[1].y() += span_.y();
297  pts[3].z() += span_.z();
298  }
299  else if (span_.y() < span_.z())
300  {
301  pts[1].x() += span_.x();
302  pts[3].z() += span_.z();
303  }
304  else
305  {
306  pts[1].x() += span_.x();
307  pts[3].y() += span_.y();
308  }
309 
310  return tpts;
311 }
312 
313 
315 {
316  return bb.overlaps(bounds());
317 }
318 
319 
320 void Foam::searchablePlate::findNearest
321 (
322  const pointField& samples,
323  const scalarField& nearestDistSqr,
324  List<pointIndexHit>& info
325 ) const
326 {
327  info.setSize(samples.size());
328 
329  forAll(samples, i)
330  {
331  info[i] = findNearest(samples[i], nearestDistSqr[i]);
332  }
333 }
334 
335 
336 void Foam::searchablePlate::findLine
337 (
338  const pointField& start,
339  const pointField& end,
340  List<pointIndexHit>& info
341 ) const
342 {
343  info.setSize(start.size());
344 
345  forAll(start, i)
346  {
347  info[i] = findLine(start[i], end[i]);
348  }
349 }
350 
351 
353 (
354  const pointField& start,
355  const pointField& end,
356  List<pointIndexHit>& info
357 ) const
358 {
359  findLine(start, end, info);
360 }
361 
362 
364 (
365  const pointField& start,
366  const pointField& end,
368 ) const
369 {
370  List<pointIndexHit> nearestInfo;
371  findLine(start, end, nearestInfo);
372 
373  info.setSize(start.size());
374  forAll(info, pointi)
375  {
376  if (nearestInfo[pointi].hit())
377  {
378  info[pointi].setSize(1);
379  info[pointi][0] = nearestInfo[pointi];
380  }
381  else
382  {
383  info[pointi].clear();
384  }
385  }
386 }
387 
388 
390 (
391  const List<pointIndexHit>& info,
392  labelList& region
393 ) const
394 {
395  region.setSize(info.size());
396  region = 0;
397 }
398 
399 
401 (
402  const List<pointIndexHit>& info,
403  vectorField& normal
404 ) const
405 {
406  normal.setSize(info.size());
407  normal = Zero;
408  forAll(normal, i)
409  {
410  normal[i][normalDir_] = 1.0;
411  }
412 }
413 
414 
416 (
417  const pointField& points,
418  List<volumeType>& volType
419 ) const
420 {
422  << "Volume type not supported for plate."
423  << exit(FatalError);
424 }
425 
426 
427 // ************************************************************************* //
Foam::addToRunTimeSelectionTable
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::VectorSpace< Vector< scalar >, scalar, 3 >::componentNames
static const char *const componentNames[]
Definition: VectorSpace.H:114
Foam::List::resize
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
Foam::searchablePlate::findLineAll
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit >> &) const
Get all intersections in order from start to end.
Definition: searchablePlate.C:364
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::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::dictionary::get
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:107
Foam::searchablePlate::overlaps
virtual bool overlaps(const boundBox &bb) const
Does any part of the surface overlap the supplied bound box?
Definition: searchablePlate.C:314
Foam::searchablePlate
Searching on finite plate. Plate has to be aligned with coordinate axes. Plate defined as origin and ...
Definition: searchablePlate.H:93
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::searchablePlate::points
virtual tmp< pointField > points() const
Get the points that define the surface.
Definition: searchablePlate.C:287
Foam::magSqr
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Foam::boundBox::overlaps
bool overlaps(const boundBox &bb) const
Overlaps/touches boundingBox?
Definition: boundBoxI.H:221
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::Field< vector >
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
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)
searchablePlate.H
samples
scalarField samples(nIntervals, Zero)
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
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::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::vector
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:51
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
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::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Vector< scalar >
Foam::List< word >
Foam::searchablePlate::coordinates
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
Definition: searchablePlate.C:264
Foam::searchablePlate::findLineAny
virtual void findLineAny(const pointField &start, const pointField &end, List< pointIndexHit > &) const
Return any intersection on segment from start to end.
Definition: searchablePlate.C:353
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
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.
Foam::direction
uint8_t direction
Definition: direction.H:52
Foam::boundBox
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:63
Foam::searchablePlate::regions
virtual const wordList & regions() const
Names of regions.
Definition: searchablePlate.C:253
Foam::point
vector point
Point is a vector.
Definition: point.H:43
Foam::searchablePlate::getNormal
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
Definition: searchablePlate.C:401
Foam::searchablePlate::boundingSpheres
virtual void boundingSpheres(pointField &centres, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
Definition: searchablePlate.C:271
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
Foam::searchablePlate::getVolumeType
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
Definition: searchablePlate.C:416
Foam::VectorSpace< Vector< scalar >, scalar, 3 >::nComponents
static constexpr direction nComponents
Number of components in this vector space.
Definition: VectorSpace.H:101
Foam::searchablePlate::getRegion
virtual void getRegion(const List< pointIndexHit > &, labelList &region) const
From a set of points and indices get the region.
Definition: searchablePlate.C:390
sample
Minimal example by using system/controlDict.functions: