39 namespace PatchFunction1Types
50 Foam::PatchFunction1Types::FilterField::RBF_typeNames_
52 { RBF_type::RBF_linear,
"linear" },
53 { RBF_type::RBF_quadratic,
"quadratic" },
54 { RBF_type::RBF_linear,
"default" },
68 const scalar radiusSqr
71 return (1 - ::
sqrt(
p.distSqr(
p0) / radiusSqr));
80 const scalar radiusSqr
83 return (1 -
p.distSqr(
p0) / radiusSqr);
117 const auto&
tree = treePtr();
122 Info<<
"have " <<
tree.nodes().size() <<
" nodes" <<
nl;
138 class BasisFunctionOp,
139 class PointWeightFieldType
141 void Foam::PatchFunction1Types::FilterField::buildWeightsImpl
143 const TreeType&
tree,
144 const RadiusSqrOp& radiusSqrOp,
145 const BasisFunctionOp& basisFuncOp,
146 const PointWeightFieldType& pointWeightFld
149 tmp<pointField> tpoints =
tree.shapes().centres();
151 const auto&
points = tpoints();
164 bool usesNeighbours =
false;
166 for (label pointi = 0; pointi <
nPoints; ++pointi)
169 auto& currAddr = addressing_[pointi];
170 auto& currWeight = weights_[pointi];
173 const scalar radiusSqr = radiusSqrOp(pointi);
175 tree.findSphere(
p0, radiusSqr, neighbours);
178 if (neighbours.size() < 2)
180 currAddr.resize(1, pointi);
181 currWeight.resize(1, scalar(1));
185 usesNeighbours =
true;
187 currAddr = neighbours.sortedToc();
188 currWeight.resize_nocopy(currAddr.size());
190 scalar sumWeight = 0;
196 currWeight[i] = basisFuncOp(
p,
p0, radiusSqr);
200 currWeight[i] *=
static_cast<scalar
>(pointWeightFld[i]);
201 sumWeight += currWeight[i];
204 for (
auto& w : currWeight)
216 if (
debug && !addressing_.empty())
222 for (
const labelList& addr : addressing_)
224 const label
n = addr.size();
230 Pout<<
"Weight neighbours: min=" << limits.min()
231 <<
" avg=" << (total / addressing_.size())
232 <<
" max=" << limits.max() <<
endl;
246 reset(
points, radius, interp);
258 reset(geom, radius,
relative, interp);
282 Pout<<
"Apply " << RBF_typeNames_[interp] <<
" filter," 283 <<
" radius=" << radius <<
nl 284 <<
"Create tree..." <<
endl;
290 const scalar radiusSqr =
sqr(radius);
294 if (interp == RBF_type::RBF_quadratic)
304 [=](label index) -> scalar {
return radiusSqr; },
322 reset(geom.
points(), radius, interp);
328 Pout<<
"Apply " << RBF_typeNames_[interp] <<
" filter," 329 << (
relative ?
" relative" :
"") <<
" radius=" << radius <<
nl 330 <<
"Create tree..." <<
endl;
336 const scalar radiusSqr =
sqr(radius);
340 if (interp == RBF_type::RBF_quadratic)
352 [&](label index) -> scalar
354 return (radiusSqr * geom.
sphere(index));
366 [=](label index) -> scalar {
return radiusSqr; },
void size(const label n)
Older name for setAddressableSize.
Holds (reference to) pointField. Encapsulation of data needed for octree searches. Used for searching for nearest point. No bounding boxes around points. Only overlaps and calcNearest are implemented, rest makes little sense.
void resize(const label len)
Adjust allocated size of list.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Output to file stream, using an OSstream.
constexpr char nl
The newline '\n' character (0x0a)
FilterField() noexcept=default
Default construct.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tf1, const word &name, const dimensionSet &dimensions, const bool initCopy=false)
Global function forwards to reuseTmpDimensionedField::New.
defineTypeNameAndDebug(FilterField, 0)
#define forAll(list, i)
Loop across all elements in list.
HashSet< label, Hash< label > > labelHashSet
A HashSet of labels, uses label hasher.
const dimensionedScalar e
Elementary charge.
static autoPtr< indexedOctree< treeDataPoint > > createTree(const pointField &points)
Construct search tree for points.
void clear()
Clear the list, i.e. set size to zero.
const scalarField & magSf() const
Face area magnitudes.
const Field< point_type > & faceCentres() const
Return face centres for patch.
Tree tree(triangles.begin(), triangles.end())
MinMax< label > labelMinMax
A label min/max range.
const Field< point_type > & points() const noexcept
Return reference to global points.
void reset()
Reset to unweighted (pass-through)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
int debug
Static debugging option.
OBJstream os(runTime.globalPath()/outputName)
static scalar quadraticWeight(const point &p, const point &p0, const scalar radiusSqr)
Quadratic basis function: 1 - (r/r0)^2.
A class representing the concept of a field of 1 used to avoid unnecessary manipulations for objects ...
vector point
Point is a vector.
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
scalar sphere(const label facei) const
The enclosing (bounding) sphere radius^2 for specified face.
Standard boundBox with extra functionality for use in octree.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
List< label > labelList
A List of labels.
static scalar linearWeight(const point &p, const point &p0, const scalar radiusSqr)
Linear basis function: 1 - (r/r0)
label nFaces() const noexcept
Number of faces in the patch.
void grow(const scalar delta)
Expand box by adjusting min/max by specified amount in each dimension.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
const volScalarField & p0
tmp< surfaceScalarField > relative(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given absolute flux in relative form.