Go to the documentation of this file.
46 searchableExtrudedCircle,
52 searchableExtrudedCircle,
61 Foam::searchableExtrudedCircle::searchableExtrudedCircle
84 radius_(
dict.get<scalar>(
"radius"))
92 vector halfSpan(0.5*bounds().span());
93 point ctr(bounds().centre());
95 bounds().min() = ctr -
mag(halfSpan) * vector::one;
96 bounds().max() = ctr +
mag(halfSpan) * vector::one;
103 bb.
min() -= point::uniform(ROOTVSMALL);
104 bb.
max() += point::uniform(ROOTVSMALL);
136 if (regions_.empty())
139 regions_[0] =
"region0";
147 return eMeshPtr_().points().size();
153 return eMeshPtr_().points();
163 centres = eMeshPtr_().points();
164 radiusSqr.setSize(centres.size());
184 const scalar nearestDist =
Foam::sqrt(nearestDistSqr[i]);
185 const scalar searchDistSqr =
Foam::sqr(nearestDist+radius_);
188 info[i] = tree.findNearest(
samples[i], searchDistSqr);
194 const scalar
s(
mag(d));
203 const scalar distToSurface = radius_-
s;
204 if (
mag(distToSurface) > nearestDist)
210 info[i].setPoint(info[i].hitPoint() + d/
s*radius_);
238 (rawLambdas-rawLambdas[0])
239 /(rawLambdas.last()-rawLambdas[0])
247 const pointIndexHit startInfo = tree.findNearest(start, maxDistSqr);
248 curvePoints[0] = startInfo.hitPoint();
249 axialVecs[0] = edges[startInfo.index()].vec(
points);
252 curvePoints.last() = endInfo.
hitPoint();
253 axialVecs.last() = edges[endInfo.
index()].vec(
points);
260 scalar endDistance = -1.0;
264 const point& start = curvePoints[0];
265 const point&
end = curvePoints.last();
267 label edgei = startInfo.index();
268 const edge& startE = edges[edgei];
270 label pointi = startE[0];
276 curveLambdas[pointi] =
mag(
points[pointi]-curvePoints[0]);
278 curveLambdas[otherPointi] = -
mag(
points[otherPointi]-curvePoints[0]);
286 const labelList& pEdges = pointEdges[pointi];
287 if (pEdges.size() == 1)
291 else if (pEdges.size() != 2)
294 <<
" is not a single path. This is not supported"
299 label oldEdgei = edgei;
300 if (pEdges[0] == oldEdgei)
309 if (edgei == endInfo.
index())
311 endDistance = curveLambdas[pointi] +
mag(
end-
points[pointi]);
320 label oldPointi = pointi;
321 pointi = edges[edgei].otherVertex(oldPointi);
323 if (curveLambdas[pointi] >= 0)
328 curveLambdas[pointi] =
329 curveLambdas[oldPointi] + edges[edgei].mag(
points);
335 if (curveLambdas[i] >= 0)
337 curveLambdas[i] /= endDistance;
350 for (label i = 1; i < curvePoints.size()-1; i++)
352 interpolator.
valueWeights(lambdas[i], indices, weights);
354 if (indices.size() == 1)
357 label pointi = indices[0];
359 label edge0 = pointEdges[pointi][0];
360 const edge& e0 = edges[edge0];
362 curvePoints[i] = weights[0]*
p0;
364 else if (indices.size() == 2)
368 axialVecs[i] = p1-
p0;
369 curvePoints[i] = weights[0]*
p0+weights[1]*p1;
372 axialVecs /=
mag(axialVecs);
387 radialStart = start-curvePoints[0];
388 radialStart -= (radialStart&axialVecs[0])*axialVecs[0];
398 vector radialEnd(
end-curvePoints.last());
399 radialEnd -= (radialEnd&axialVecs.last())*axialVecs.last();
402 vector projectedEnd = radialEnd;
403 projectedEnd -= (projectedEnd&axialVecs[0])*axialVecs[0];
406 qProjectedEnd =
quaternion(projectedEnd, 0.0);
411 for (label i = 1; i < lambdas.size()-1; i++)
416 radialDir -= (radialDir & axialVecs[i]) * axialVecs.last();
419 info[i] =
pointIndexHit(
true, curvePoints[i]+radius_*radialDir, 0);
446 normal.setSize(info.size());
460 normal[i] = info[i].hitPoint()-curvePt.
hitPoint();
465 normal[i] -= (normal[i] & axialVec) * axialVec;
466 normal[i].normalise();
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
virtual const pointField & points() const
Return raw points.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual tmp< pointField > coordinates() const
Get representative set of element coordinates.
label index() const
Return index.
virtual void findParametricNearest(const point &start, const point &end, const scalarField &lambdas, const scalarField &nearestDistSqr, List< pointIndexHit > &) const
Unique to parametric geometry: given points find.
A class for handling words, derived from Foam::string.
virtual void findNearest(const pointField &sample, const scalarField &nearestDistSqr, List< pointIndexHit > &) const
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))
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
quaternion slerp(const quaternion &qa, const quaternion &qb, const scalar t)
Spherical linear interpolation of quaternions.
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const
From a set of points and indices get the normal.
Standard boundBox with extra functionality for use in octree.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
const point & max() const
Maximum describing the bounding box.
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
const labelListList & pointEdges() const
Vector< Cmpt > & normalise()
Normalise the vector by its magnitude.
const Time & time() const
Return time.
const point & min() const
Minimum describing the bounding box.
#define forAll(list, i)
Loop across all elements in list.
Quaternion class used to perform rotations in 3D space.
virtual const wordList & regions() const
Names of regions.
virtual ~searchableExtrudedCircle()
Destructor.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
const Point & hitPoint() const
Return hit point.
word name(const complex &c)
Return string representation of complex.
Base class of (analytical or triangulated) surface. Encapsulates all the search routines....
const pointField & points() const
Return points.
addNamedToRunTimeSelectionTable(topoSetCellSource, badQualityToCell, word, badQuality)
Non-pointer based hierarchical recursive searching.
virtual void getRegion(const List< pointIndexHit > &, labelList ®ion) const
From a set of points and indices get the region.
Holds data for octree to work on an edges subset.
scalarField samples(nIntervals, Zero)
const edgeList & edges() const
Return edges.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
virtual bool valueWeights(const scalar t, labelList &indices, scalarField &weights) const
Calculate weights and indices to calculate t from samples.
Macros for easy insertion into run-time selection tables.
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
virtual void boundingSpheres(pointField ¢res, scalarField &radiusSqr) const
Get bounding spheres (centre and radius squared), one per element.
errorManipArg< error, int > exit(error &err, const int errNo=1)
label otherVertex(const label pointLabel) const
Given the point label for one vertex, return the other one.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
PointIndexHit< point > pointIndexHit
vector vec(const UList< point > &pts) const
Return the vector (end - start)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
labelList identity(const label len, label start=0)
Create identity map of the given length with (map[i] == i)
A bounding box defined in terms of min/max extrema points.
const volScalarField & p0
void setSize(const label newSize)
Alias for resize(const label)
const word & constant() const
Return constant name.
defineTypeNameAndDebug(combustionModel, 0)
Mesh data needed to do the Finite Area discretisation.
virtual label size() const
Range of local indices that can be returned.
vector transform(const vector &v) const
Rotate the given vector.