distanceSurface.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-2016 OpenFOAM Foundation
9 Copyright (C) 2016-2021 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 "distanceSurface.H"
30#include "dictionary.H"
31#include "volFields.H"
34#include "fvMesh.H"
35
36// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37
38namespace Foam
39{
41}
42
43
44// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
45
46const Foam::Enum
47<
48 Foam::distanceSurface::topologyFilterType
49>
50Foam::distanceSurface::topoFilterNames_
51({
52 { topologyFilterType::NONE, "none" },
53 { topologyFilterType::LARGEST_REGION, "largestRegion" },
54 { topologyFilterType::NEAREST_POINTS, "nearestPoints" },
55 { topologyFilterType::PROXIMITY_REGIONS, "proximityRegions" },
56 { topologyFilterType::PROXIMITY_FACES, "proximityFaces" },
57 { topologyFilterType::PROXIMITY_FACES, "proximity" },
58});
59
60
61// * * * * * * * * * * * * * * * Local Functions * * * * * * * * * * * * * * //
62
63namespace Foam
64{
65
66// Check that all point hits are valid
67static inline void checkAllHits(const UList<pointIndexHit>& nearest)
68{
69 label notHit = 0;
70 for (const pointIndexHit& pHit : nearest)
71 {
72 if (!pHit.hit())
73 {
74 ++notHit;
75 }
76 }
77
78 if (notHit)
79 {
81 << "Had " << notHit << " faces/cells from "
82 << nearest.size() << " without a point hit." << nl
83 << "May be caused by a severely degenerate input surface" << nl
84 << abort(FatalError);
85 }
86}
87
88
89// Normal distance from surface hit point to a point in the mesh
90static inline scalar normalDistance_zero
91(
92 const point& pt,
93 const pointIndexHit& pHit,
94 const vector& norm
95)
96{
97 const vector diff(pt - pHit.point());
98
99 return (diff & norm);
100}
101
102
103// Signed distance from surface hit point to a point in the mesh,
104// the sign is dictated by the normal
105static inline scalar normalDistance_nonzero
106(
107 const point& pt,
108 const pointIndexHit& pHit,
109 const vector& norm
110)
111{
112 const vector diff(pt - pHit.point());
113 const scalar normDist = (diff & norm);
114
115 return Foam::sign(normDist) * Foam::mag(diff);
116}
117
118
119// Normal distance from surface hit point to a point in the mesh
120static inline void calcNormalDistance_zero
121(
123 const pointField& points,
124 const List<pointIndexHit>& nearest,
125 const vectorField& normals
126)
127{
128 forAll(nearest, i)
129 {
130 distance[i] =
131 normalDistance_zero(points[i], nearest[i], normals[i]);
132 }
133}
134
135
136// Signed distance from surface hit point -> point in the mesh,
137// the sign is dictated by the normal
139(
141 const pointField& points,
142 const List<pointIndexHit>& nearest,
143 const vectorField& normals
144)
145{
146 forAll(nearest, i)
147 {
148 distance[i] =
149 normalDistance_nonzero(points[i], nearest[i], normals[i]);
150 }
151}
152
153
154// Close to the surface: normal distance from surface hit point
155// Far from surface: distance from surface hit point
156//
157// Note
158// This switch may be helpful when working directly with
159// distance/gradient fields. Has low overhead otherwise.
160// May be replaced in the future (2020-11)
162(
164 const bitSet& ignoreLocation,
165 const pointField& points,
166 const List<pointIndexHit>& nearest,
167 const vectorField& normals
168)
169{
170 forAll(nearest, i)
171 {
172 if (ignoreLocation.test(i))
173 {
174 distance[i] =
175 normalDistance_nonzero(points[i], nearest[i], normals[i]);
176 }
177 else
178 {
179 distance[i] =
180 normalDistance_zero(points[i], nearest[i], normals[i]);
181 }
182 }
183}
184
185
186// Flat surfaces (eg, a plane) have an extreme change in
187// the normal at the edge, which creates a zero-crossing
188// extending to infinity.
189//
190// Ad hoc treatment: require that the surface hit point is within a
191// somewhat generous bounding box for the cell (+10%)
192//
193// Provisioning for filtering based on the cell points,
194// but its usefulness remains to be seen (2020-12-09)
195template<bool WantPointFilter = false>
197(
198 bitSet& ignoreCells,
199 const List<pointIndexHit>& nearest,
200 const polyMesh& mesh,
201 const scalar boundBoxInflate = 0.1 // 10% to catch corners
202)
203{
204 // A deny filter. Initially false (accept everything)
205 ignoreCells.resize(mesh.nCells());
206
207 bitSet pointFilter;
208 if (WantPointFilter)
209 {
210 // Create as accept filter. Initially false (deny everything)
211 pointFilter.resize(mesh.nPoints());
212 }
213
214 boundBox cellBb;
215
216 forAll(nearest, celli)
217 {
218 const point& pt = nearest[celli].point();
219
220 const labelList& cPoints = mesh.cellPoints(celli);
221
222 cellBb.clear();
223 cellBb.add(mesh.points(), cPoints);
224 cellBb.inflate(boundBoxInflate);
225
226 if (!cellBb.contains(pt))
227 {
228 ignoreCells.set(celli);
229 }
230 else if (WantPointFilter)
231 {
232 // Good cell candidate, accept its points
233 pointFilter.set(cPoints);
234 }
235 }
236
237 // Flip from accept to deny filter
238 pointFilter.flip();
239
240 return pointFilter;
241}
242
243
244} // End namespace Foam
245
246
247// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
248
250(
251 const word& defaultSurfaceName,
252 const polyMesh& mesh,
253 const dictionary& dict
254)
255:
256 mesh_(mesh),
257 geometryPtr_
258 (
260 (
261 dict.get<word>("surfaceType"),
263 (
264 dict.getOrDefault("surfaceName", defaultSurfaceName),
265 mesh.time().constant(), // directory
266 "triSurface", // instance
267 mesh.time(), // registry
268 IOobject::MUST_READ,
269 IOobject::NO_WRITE
270 ),
271 dict
272 )
273 ),
274 distance_(dict.getOrDefault<scalar>("distance", 0)),
275 withZeroDistance_(equal(distance_, 0)),
276 withSignDistance_
277 (
278 withZeroDistance_
279 || (distance_ < 0)
280 || dict.getOrDefault<bool>("signed", true)
281 ),
282
283 isoParams_
284 (
285 dict,
286 isoSurfaceParams::ALGO_TOPO,
287 isoSurfaceParams::filterType::DIAGCELL
288 ),
289 topoFilter_
290 (
291 topoFilterNames_.getOrDefault
292 (
293 "topology",
294 dict,
295 topologyFilterType::NONE
296 )
297 ),
298 nearestPoints_(),
299 maxDistanceSqr_(Foam::sqr(GREAT)),
300 absProximity_(dict.getOrDefault<scalar>("absProximity", 1e-5)),
301 cellDistancePtr_(nullptr),
302 pointDistance_(),
303 surface_(),
304 meshCells_(),
305 isoSurfacePtr_(nullptr)
306{
307 if (topologyFilterType::NEAREST_POINTS == topoFilter_)
308 {
309 dict.readEntry("nearestPoints", nearestPoints_);
310 }
311
312 if (dict.readIfPresent("maxDistance", maxDistanceSqr_))
313 {
314 maxDistanceSqr_ = Foam::sqr(maxDistanceSqr_);
315 }
316}
317
318
320(
321 const polyMesh& mesh,
322 const word& surfaceType,
323 const word& surfaceName,
324 const isoSurfaceParams& params,
325 const bool interpolate
326)
327:
329 (
330 mesh,
332 surfaceType,
333 surfaceName,
334 scalar(0),
335 true, // redundant - must be signed
336 params
337 )
338{}
339
340
342(
343 const polyMesh& mesh,
344 const bool interpolate,
345 const word& surfaceType,
346 const word& surfaceName,
347 const scalar distance,
348 const bool useSignedDistance,
349 const isoSurfaceParams& params
350)
351:
353 (
354 mesh,
357 (
358 surfaceType,
360 (
361 surfaceName, // name
362 mesh.time().constant(), // directory
363 "triSurface", // instance
364 mesh.time(), // registry
365 IOobject::MUST_READ,
366 IOobject::NO_WRITE
367 ),
368 dictionary()
369 ),
370 distance,
371 useSignedDistance,
372 params
373 )
374{}
375
376
378(
379 const polyMesh& mesh,
380 const bool interpolate,
382 const scalar distance,
383 const bool useSignedDistance,
384 const isoSurfaceParams& params
385)
386:
387 mesh_(mesh),
388 geometryPtr_(surface),
389 distance_(distance),
390 withZeroDistance_(equal(distance_, 0)),
391 withSignDistance_
392 (
393 withZeroDistance_
394 || (distance_ < 0)
395 || useSignedDistance
396 ),
397
398 isoParams_(params),
399 topoFilter_(topologyFilterType::NONE),
400 nearestPoints_(),
401 maxDistanceSqr_(Foam::sqr(GREAT)),
402 absProximity_(1e-5),
403 cellDistancePtr_(nullptr),
404 pointDistance_(),
405 surface_(),
406 meshCells_(),
407 isoSurfacePtr_(nullptr)
408{}
409
410
411// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
412
414{
415 if (debug)
416 {
417 Pout<< "distanceSurface::createGeometry updating geometry." << endl;
418 }
419
420 // Clear any previously stored topologies
421 isoSurfacePtr_.reset(nullptr);
422 surface_.clear();
423 meshCells_.clear();
424
425 // Doing searches on this surface
426 const searchableSurface& geom = geometryPtr_();
427
428 const fvMesh& fvmesh = static_cast<const fvMesh&>(mesh_);
429
430 // Distance to cell centres
431 // ~~~~~~~~~~~~~~~~~~~~~~~~
432
433 cellDistancePtr_.reset
434 (
436 (
438 (
439 "distanceSurface.cellDistance",
440 fvmesh.time().timeName(),
441 fvmesh.time(),
444 false
445 ),
446 fvmesh,
448 )
449 );
450 auto& cellDistance = *cellDistancePtr_;
451
452
453 // For distance = 0 we apply additional geometric filtering
454 // to limit the extent of open edges.
455 //
456 // Does not work with ALGO_POINT
457
458 bitSet ignoreCells, ignoreCellPoints;
459
460 const bool filterCells =
461 (
462 withZeroDistance_
463 && isoParams_.algorithm() != isoSurfaceParams::ALGO_POINT
464 );
465
466
467 // Internal field
468 {
469 const pointField& cc = fvmesh.C();
470 scalarField& fld = cellDistance.primitiveFieldRef();
471
472 List<pointIndexHit> nearest;
473 geom.findNearest
474 (
475 cc,
476 // Use initialized field (GREAT) to limit search too
477 fld,
478 nearest
479 );
480 checkAllHits(nearest);
481
482 // Geometric pre-filtering when distance == 0
483
484 // NOTE (2021-05-31)
485 // Can skip the prefilter if we use proximity-regions filter anyhow
486 // but it makes the iso algorithm more expensive and doesn't help
487 // unless we start relying on area-based weighting for rejecting regions.
488
489 if (filterCells)
490 {
491 ignoreCellPoints = simpleGeometricFilter<false>
492 (
493 ignoreCells,
494 nearest,
495 fvmesh,
496
497 // Inflate bound box.
498 // - To catch corners: approx. 10%
499 // - Extra generous for PROXIMITY_REGIONS
500 // (extra weighting for 'bad' faces)
501 (
502 topologyFilterType::PROXIMITY_REGIONS == topoFilter_
503 ? 1
504 : 0.1
505 )
506 );
507 }
508
509 if (withSignDistance_)
510 {
511 vectorField norms;
512 geom.getNormal(nearest, norms);
513
514 if (filterCells)
515 {
516 // With inside/outside switching (see note above)
518 (
519 fld,
520 ignoreCells,
521 cc,
522 nearest,
523 norms
524 );
525 }
526 else if (withZeroDistance_)
527 {
528 calcNormalDistance_zero(fld, cc, nearest, norms);
529 }
530 else
531 {
532 calcNormalDistance_nonzero(fld, cc, nearest, norms);
533 }
534 }
535 else
536 {
537 calcAbsoluteDistance(fld, cc, nearest);
538 }
539 }
540
541
542 // Patch fields
543 {
544 forAll(fvmesh.C().boundaryField(), patchi)
545 {
546 const pointField& cc = fvmesh.C().boundaryField()[patchi];
547 scalarField& fld = cellDistance.boundaryFieldRef()[patchi];
548
549 List<pointIndexHit> nearest;
550 geom.findNearest
551 (
552 cc,
553 scalarField(cc.size(), GREAT),
554 nearest
555 );
556 checkAllHits(nearest);
557
558 if (withSignDistance_)
559 {
560 vectorField norms;
561 geom.getNormal(nearest, norms);
562
563 if (withZeroDistance_)
564 {
565 // Slight inconsistency in boundary vs interior when
566 // cells are filtered, but the patch fields are only
567 // used by isoSurfacePoint, and filtering is disabled
568 // for that anyhow.
569
570 calcNormalDistance_zero(fld, cc, nearest, norms);
571 }
572 else
573 {
574 calcNormalDistance_nonzero(fld, cc, nearest, norms);
575 }
576 }
577 else
578 {
579 calcAbsoluteDistance(fld, cc, nearest);
580 }
581 }
582 }
583
584
585 // On processor patches the mesh.C() will already be the cell centre
586 // on the opposite side so no need to swap cellDistance.
587
588
589 // Distance to points
590 pointDistance_.resize(fvmesh.nPoints());
591 pointDistance_ = GREAT;
592 {
593 const pointField& pts = fvmesh.points();
594 scalarField& fld = pointDistance_;
595
596 List<pointIndexHit> nearest;
597 geom.findNearest
598 (
599 pts,
600 // Use initialized field (GREAT) to limit search too
601 pointDistance_,
602 nearest
603 );
604 checkAllHits(nearest);
605
606 if (withSignDistance_)
607 {
608 vectorField norms;
609 geom.getNormal(nearest, norms);
610
611 if (filterCells)
612 {
613 // With inside/outside switching (see note above)
615 (
616 fld,
617 ignoreCellPoints,
618 pts,
619 nearest,
620 norms
621 );
622 }
623 else if (withZeroDistance_)
624 {
625 calcNormalDistance_zero(fld, pts, nearest, norms);
626 }
627 else
628 {
629 calcNormalDistance_nonzero(fld, pts, nearest, norms);
630 }
631 }
632 else
633 {
634 calcAbsoluteDistance(fld, pts, nearest);
635 }
636 }
637
638
639 // Don't need ignoreCells if there is nothing to ignore.
640 if (ignoreCells.none())
641 {
642 ignoreCells.clearStorage();
643 }
644 else if (filterCells && topologyFilterType::NONE != topoFilter_)
645 {
646 // For refine blocked cells (eg, checking actual cells cut)
647 isoSurfaceBase isoCutter
648 (
649 mesh_,
650 cellDistance,
651 pointDistance_,
652 distance_
653 );
654
655 if (topologyFilterType::LARGEST_REGION == topoFilter_)
656 {
657 refineBlockedCells(ignoreCells, isoCutter);
658 filterKeepLargestRegion(ignoreCells);
659 }
660 else if (topologyFilterType::NEAREST_POINTS == topoFilter_)
661 {
662 refineBlockedCells(ignoreCells, isoCutter);
663 filterKeepNearestRegions(ignoreCells);
664 }
665
666 // Note: apply similar filtering for PROXIMITY_REGIONS later instead
667 }
668
669 // Don't need point filter beyond this point
670 ignoreCellPoints.clearStorage();
671
672
673 if (debug)
674 {
675 Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
676 cellDistance.write();
677
678 pointScalarField pDist
679 (
681 (
682 "distanceSurface.pointDistance",
683 fvmesh.time().timeName(),
684 fvmesh.time(),
687 false
688 ),
689 pointMesh::New(fvmesh),
691 );
692 pDist.primitiveFieldRef() = pointDistance_;
693
694 Pout<< "Writing point distance:" << pDist.objectPath() << endl;
695 pDist.write();
696 }
697
698 isoSurfacePtr_.reset
699 (
701 (
702 isoParams_,
703 cellDistance,
704 pointDistance_,
705 distance_,
706 ignoreCells
707 )
708 );
709
710
711 // Restrict ignored cells to those actually cut
712 if (filterCells && topologyFilterType::PROXIMITY_REGIONS == topoFilter_)
713 {
714 isoSurfaceBase isoCutter
715 (
716 mesh_,
717 cellDistance,
718 pointDistance_,
719 distance_
720 );
721
722 refineBlockedCells(ignoreCells, isoCutter);
723 }
724
725
726 // ALGO_POINT still needs cell, point fields (for interpolate)
727 // The others can do straight transfer
728
729 // But also flatten into a straight transfer for proximity filtering
730 if
731 (
732 isoParams_.algorithm() != isoSurfaceParams::ALGO_POINT
733 || topologyFilterType::PROXIMITY_FACES == topoFilter_
734 || topologyFilterType::PROXIMITY_REGIONS == topoFilter_
735 )
736 {
737 surface_.transfer(static_cast<meshedSurface&>(*isoSurfacePtr_));
738 meshCells_.transfer(isoSurfacePtr_->meshCells());
739
740 isoSurfacePtr_.reset(nullptr);
741 cellDistancePtr_.reset(nullptr);
742 pointDistance_.clear();
743 }
744
745 if (topologyFilterType::PROXIMITY_FACES == topoFilter_)
746 {
747 filterFaceProximity();
748 }
749 else if (topologyFilterType::PROXIMITY_REGIONS == topoFilter_)
750 {
751 filterRegionProximity(ignoreCells);
752 }
753
754 if (debug)
755 {
756 print(Pout, debug);
757 Pout<< endl;
758 }
759}
760
761
763{
764 os << " surface:" << surfaceName()
765 << " distance:" << distance()
766 << " topology:" << topoFilterNames_[topoFilter_];
767
768 isoParams_.print(os);
769
770 if (level)
771 {
772 os << " faces:" << surface().surfFaces().size()
773 << " points:" << surface().points().size();
774 }
775}
776
777
778// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
fileName objectPath() const
The complete path + object name.
Definition: IOobjectI.H:214
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
void resize(const label numElem, const unsigned int val=0u)
Reset addressable list size, does not shrink the allocated size.
Definition: PackedListI.H:409
void clearStorage()
Clear the list and delete storage.
Definition: PackedListI.H:520
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:66
const point_type & point() const noexcept
Return point, no checks.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: UList.H:94
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:66
void flip()
Invert all bits in the addressable region.
Definition: bitSetI.H:634
void set(const bitSet &bitset)
Set specified bits from another bitset.
Definition: bitSetI.H:590
bool none() const
True if no bits in this bitset are set.
Definition: bitSetI.H:488
bool test(const label pos) const
Test value at specified position, never auto-vivify entries.
Definition: bitSetI.H:521
A bounding box defined in terms of min/max extrema points.
Definition: boundBox.H:64
bool contains(const point &pt) const
Contains point? (inside or on edge)
Definition: boundBoxI.H:271
void add(const boundBox &bb)
Extend to include the second box.
Definition: boundBoxI.H:191
void inflate(const scalar s)
Inflate box by factor*mag(span) in all dimensions.
Definition: boundBox.C:175
void clear()
Clear bounding box and make it an inverted box.
Definition: boundBoxI.H:184
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
A surface defined by a distance from an input searchable surface. Uses an iso-surface algorithm (cell...
void createGeometry()
Create/recreate the distance surface.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const volVectorField & C() const
Return cell centres as volVectorField.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Low-level components common to various iso-surface algorithms.
Preferences for controlling iso-surface algorithms.
scalar print()
Print to screen.
constant condensation/saturation model.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1083
label nPoints() const noexcept
Number of mesh points.
const labelListList & cellPoints() const
label nCells() const noexcept
Number of mesh cells.
void reset(const label nPoints, const label nInternalFaces, const label nFaces, const label nCells)
Reset this primitiveMesh given the primitive array sizes.
virtual bool write(const bool valid=true) const
Write using setting from DB.
Base class of (analytical or triangulated) surface. Encapsulates all the search routines....
virtual void findNearest(const pointField &sample, const scalarField &nearestDistSqr, List< pointIndexHit > &) const =0
virtual void getNormal(const List< pointIndexHit > &, vectorField &normal) const =0
From a set of points and indices get the normal.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
bool
Definition: EEqn.H:20
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
OBJstream os(runTime.globalPath()/outputName)
const pointField & points
Namespace for OpenFOAM.
static void calcNormalDistance_zero(scalarField &distance, const pointField &points, const List< pointIndexHit > &nearest, const vectorField &normals)
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
static scalar normalDistance_zero(const point &pt, const pointIndexHit &pHit, const vector &norm)
dimensionedScalar sign(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
static void calcNormalDistance_filtered(scalarField &distance, const bitSet &ignoreLocation, const pointField &points, const List< pointIndexHit > &nearest, const vectorField &normals)
static scalar normalDistance_nonzero(const point &pt, const pointIndexHit &pHit, const vector &norm)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
bool equal(const T &s1, const T &s2)
Compare two values for equality.
Definition: doubleFloat.H:46
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:378
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
bool interpolate(const vector &p1, const vector &p2, const vector &o, vector &n, scalar l)
Definition: curveTools.C:75
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
error FatalError
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
static void calcNormalDistance_nonzero(scalarField &distance, const pointField &points, const List< pointIndexHit > &nearest, const vectorField &normals)
static void checkAllHits(const UList< pointIndexHit > &nearest)
static bitSet simpleGeometricFilter(bitSet &ignoreCells, const List< pointIndexHit > &nearest, const polyMesh &mesh, const scalar boundBoxInflate=0.1)
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