36 #include "indexedVertex.H" 60 template<
class Triangulation,
class Type>
63 const Triangulation&
mesh,
68 auto& newField = tNewField.
ref();
75 typename Triangulation::Finite_vertices_iterator vit =
76 mesh.finite_vertices_begin();
77 vit !=
mesh.finite_vertices_end();
89 newField.resize(added);
108 typename T::Finite_vertices_iterator vit =
mesh.finite_vertices_begin();
109 vit !=
mesh.finite_vertices_end();
118 std::list<typename T::Vertex_handle> adjVerts;
119 mesh.finite_adjacent_vertices(vit, std::back_inserter(adjVerts));
125 typename std::list<typename T::Vertex_handle>::const_iterator
126 adjVertI = adjVerts.begin();
127 adjVertI != adjVerts.end();
131 typename T::Vertex_handle vh = *adjVertI;
137 globalIndexing.toGlobal(vh->procIndex(), vh->index())
142 pointPoints[vit->index()].
transfer(indices);
167 typename T::Finite_vertices_iterator vit =
mesh.finite_vertices_begin();
168 vit !=
mesh.finite_vertices_end();
177 alignments[vit->index()] = vit->alignment();
195 typename T::Finite_vertices_iterator vit =
mesh.finite_vertices_begin();
196 vit !=
mesh.finite_vertices_end();
216 const label maxRefinementIterations,
217 const scalar defaultCellSize
220 for (label iter = 0; iter < maxRefinementIterations; ++iter)
226 CellSizeDelaunay::Finite_cells_iterator cit =
227 mesh.finite_cells_begin();
228 cit !=
mesh.finite_cells_end();
232 const point newPoint =
237 cit->vertex(0)->point(),
238 cit->vertex(1)->point(),
239 cit->vertex(2)->point(),
240 cit->vertex(3)->point()
244 if (geometryToConformTo.
inside(newPoint))
246 ptsToInsert.
append(newPoint);
267 int main(
int argc,
char *argv[])
274 label maxRefinementIterations = 2;
275 label maxSmoothingIterations = 200;
276 scalar minResidual = 0;
277 scalar defaultCellSize = 0.001;
278 scalar nearFeatDistSqrCoeff = 1
e-8;
306 "cvSearchableSurfaces",
313 foamyHexMeshDict.subDict(
"geometry"),
314 foamyHexMeshDict.getOrDefault(
"singleRegionName",
true)
322 foamyHexMeshDict.subDict(
"surfaceConformation")
335 foamyHexMeshDict.subDict(
"backgroundMeshDecomposition")
351 const label surfI = geometryToConformTo.
surfaces()[sI];
354 geometryToConformTo.
geometry()[surfI];
356 Info<<
nl <<
"Inserting points from surface " <<
surface.name()
375 nearFeatDistSqrCoeff,
386 geometryToConformTo.
features()[infoFeature];
394 pointAlignment() += norms[nI];
397 pointAlignment().normalize();
398 pointAlignment().orthogonalize();
405 nearFeatDistSqrCoeff,
413 geometryToConformTo.
features()[infoFeature];
421 pointAlignment() += norms[nI];
424 pointAlignment().normalize();
425 pointAlignment().orthogonalize();
433 surface.findNearest(ptField, distField, infoList);
436 surface.getNormal(infoList, normals);
438 pointAlignment.
set(
new triad(normals[0]));
446 CellSizeDelaunay::Vertex_handle vh =
mesh.
insert 457 CellSizeDelaunay::Vertex_handle vh =
mesh.
insert 474 maxRefinementIterations,
500 CellSizeDelaunay::Finite_vertices_iterator vit =
501 mesh.finite_vertices_begin();
502 vit !=
mesh.finite_vertices_end();
506 if (vit->nearBoundary())
508 fixedAlignments[vit->index()] = vit->alignment();
514 for (label iter = 0; iter < maxSmoothingIterations; iter++)
516 Info<<
"Iteration " << iter;
518 meshDistributor().distribute(
points);
519 meshDistributor().distribute(alignments);
527 const labelList& pPoints = pointPoints[pI];
534 const triad& oldTriad = alignments[pI];
535 triad& newTriad = triadAv[pI];
538 const triad& fixedAlignment = fixedAlignments[pI];
540 forAll(pPoints, adjPointi)
542 const label adjPointIndex = pPoints[adjPointi];
548 triad tmpTriad = alignments[adjPointIndex];
552 if (tmpTriad.
set(dir))
554 tmpTriad[dir] *= (1.0/dist);
558 newTriad += tmpTriad;
567 forAll(fixedAlignment, dirI)
569 if (fixedAlignment.
set(dirI))
577 forAll(fixedAlignment, dirI)
579 if (fixedAlignment.
set(dirI))
581 newTriad.
align(fixedAlignment[dirI]);
585 else if (nFixed == 2)
587 forAll(fixedAlignment, dirI)
589 if (fixedAlignment.
set(dirI))
591 newTriad[dirI] = fixedAlignment[dirI];
601 else if (nFixed == 3)
603 forAll(fixedAlignment, dirI)
605 if (fixedAlignment.
set(dirI))
607 newTriad[dirI] = fixedAlignment[dirI];
621 scalar dotProd = (oldTriad[dir] & newTriad[dir]);
622 scalar
diff =
mag(dotProd) - 1.0;
631 alignments[pI] = triadAv[pI].sortxyz();
636 Info<<
", Residual = " << residual <<
endl;
638 if (residual <= minResidual)
650 const triad& tri = alignments[pI];
687 filterFarPoints(
mesh, sizes)
700 filterFarPoints(
mesh, alignments)
705 alignmentsIO.write();
const labelListList & featurePointNormals() const
Return the indices of the normals that are adjacent to the.
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
void size(const label n)
Older name for setAddressableSize.
fileName path() const
Return path = rootPath/caseName. Same as TimePaths::path()
pointFromPoint topoint(const Point &P)
const labelListList & edgeNormals() const
Return the indices of the normals that are adjacent to the.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
void transfer(List< T > &list)
Transfer the contents of the argument List into this list and annul the argument list.
void align(const vector &v)
Align this triad with the given vector v.
T & ref() const
Return non-const reference to the contents of a non-null managed pointer.
Output to file stream as an OSstream, normally using std::ofstream for the actual output...
constexpr char nl
The newline '\n' character (0x0a)
bool empty() const noexcept
True if List is empty (ie, size() is zero)
Ostream & endl(Ostream &os)
Add newline and flush stream.
PointIndexHit< point > pointIndexHit
A PointIndexHit with a 3D point.
static bool & parRun() noexcept
Test if this a parallel run.
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.
This class describes the interaction of an object (often a face) and a point. It carries the info of ...
Base class of (analytical or triangulated) surface. Encapsulates all the search routines. WIP.
Ignore writing from objectRegistry::writeObject()
static int myProcNo(const label communicator=worldComm)
Rank of this process in the communicator (starting from masterNo()). Can be negative if the process i...
void orthogonalize()
Same as orthogonalise.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Perform reduction on a copy, using specified binary operation.
bool set(const direction d) const
Is the vector in the direction d set.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
#define forAll(list, i)
Loop across all elements in list.
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
void normalize()
Same as normalise.
A list of faces which address into the list of points.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
vectorField pointField
pointField is a vectorField.
const dimensionedScalar e
Elementary charge.
void setSize(const label n)
Alias for resize()
Generic templated field type.
Ostream & printInfo(Ostream &os) const
Print information.
const word & system() const noexcept
Return system name.
virtual bool write(const bool writeOnProc=true) const
Write using setting from DB.
Container for searchableSurfaces. The collection is specified as a dictionary. For example...
void append(const T &val)
Copy append an element to the end of this list.
label index() const noexcept
Return the hit index.
const wordList surface
Standard surface field types (scalar, vector, tensor, etc)
Field< triad > triadField
Specialisation of Field<T> for triad.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
const word & constant() const noexcept
Return constant name.
Ostream & printExecutionTime(OSstream &os) const
Print the elapsed ExecutionTime (cpu-time), ClockTime.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Representation of a 3D Cartesian coordinate system as a Vector of row vectors.
bool hit() const noexcept
Is there a hit?
vector point
Point is a vector.
void set(T *p) noexcept
Deprecated(2018-02) Identical to reset().
Automatically write from objectRegistry::writeObject()
messageStream Info
Information stream (stdout output on master, null elsewhere)
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
Pointer management similar to std::unique_ptr, with some additional methods and type checking...
A class for managing temporary objects.
static autoPtr< T > New(Args &&... args)
Construct autoPtr with forwarding arguments.
Defines the attributes of an object for which implicit objectRegistry management is supported...
PrimitivePatch< List< face >, const pointField > bMesh
Holder of faceList and points. (v.s. e.g. primitivePatch which references points) ...
A primitive field of type <T> with automated input and output.
void reduce(T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce inplace (cf. MPI Allreduce) using linear/tree communication schedule.