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-------------------------------------------------------------------------------
11License
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
34namespace Foam
35{
38 (
41 dict
42 );
44 (
47 dict,
48 plate
49 );
50}
51
52
53// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
54
55Foam::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)
93Foam::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
139Foam::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
221(
222 const IOobject& io,
223 const point& origin,
224 const vector& span
225)
226:
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
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
320void Foam::searchablePlate::findNearest
321(
322 const pointField& samples,
323 const scalarField& nearestDistSqr,
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
336void Foam::searchablePlate::findLine
337(
338 const pointField& start,
339 const pointField& end,
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,
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// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addNamedToRunTimeSelectionTable(baseType, thisType, argNames, lookupName)
Add to construction table with 'lookupName' as the key.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Minimal example by using system/controlDict.functions:
void max(const dimensioned< Type > &dt)
Use the maximum of the field and specified value.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
void setSize(const label n)
Alias for resize()
Definition: List.H:218
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:116
int overlaps
Flag to control which overlap calculations are performed.
Definition: PDRparams.H:97
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:66
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:64
bool overlaps(const boundBox &bb) const
Overlaps/touches boundingBox?
Definition: boundBoxI.H:221
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
static const char *const componentNames[]
Definition: bool.H:104
static constexpr direction nComponents
Number of components in bool is 1.
Definition: bool.H:98
Searching on finite plate. Plate has to be aligned with coordinate axes. Plate defined as origin and ...
virtual void getVolumeType(const pointField &, List< volumeType > &) const
Determine type (inside/outside/mixed) for point. unknown if.
virtual void findLineAll(const pointField &start, const pointField &end, List< List< pointIndexHit > > &) const
Get all intersections in order from start to end.
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 void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
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 tmp< pointField > points() const
Get the points that define the surface.
Base class of (analytical or triangulated) surface. Encapsulates all the search routines....
virtual const boundBox & bounds() const
Return const reference to boundBox.
A class for managing temporary objects.
Definition: tmp.H:65
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
const pointField & points
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
#define DebugInFunction
Report an information message using Foam::Info.
Namespace for OpenFOAM.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
PointIndexHit< point > pointIndexHit
A PointIndexHit for 3D points.
Definition: pointIndexHit.H:46
vector point
Point is a vector.
Definition: point.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
uint8_t direction
Definition: direction.H:56
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
dictionary dict
volScalarField & e
Definition: createFields.H:11
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333
scalarField samples(nIntervals, Zero)