57static const scalar defaultMergeTol = 1
e-6;
60scalar getMergeDistance
67 const scalar mergeTol =
71 Foam::pow(scalar(10), -scalar(IOstream::defaultPrecision()));
73 Info<<
"Merge tolerance : " << mergeTol <<
nl
74 <<
"Write tolerance : " << writeTol <<
endl;
76 if (
runTime.writeFormat() == IOstream::ASCII && mergeTol < writeTol)
79 <<
"Your current settings specify ASCII writing with "
80 << IOstream::defaultPrecision() <<
" digits precision." <<
endl
81 <<
"Your merging tolerance (" << mergeTol <<
") is finer than this."
83 <<
"Please change your writeFormat to binary"
84 <<
" or increase the writePrecision" <<
endl
85 <<
"or adjust the merge tolerance (-mergeTol)."
89 scalar mergeDist = mergeTol * bb.
mag();
91 Info<<
"Overall meshes bounding box : " << bb <<
nl
92 <<
"Relative tolerance : " << mergeTol <<
nl
93 <<
"Absolute matching distance : " << mergeDist <<
nl
107 const labelList& pPatches =
mesh.globalData().processorPatches();
108 patchNeiProcNo[Pstream::myProcNo()].setSize(pPatches.
size());
109 patchSize[Pstream::myProcNo()].setSize(pPatches.
size());
114 mesh.boundaryMesh()[pPatches[i]]
116 patchNeiProcNo[Pstream::myProcNo()][i] = ppp.
neighbProcNo();
117 patchSize[Pstream::myProcNo()][i] = ppp.size();
119 Pstream::gatherList(patchNeiProcNo);
120 Pstream::gatherList(patchSize);
127 label maxProcCells = 0;
128 label totProcFaces = 0;
129 label maxProcPatches = 0;
130 label totProcPatches = 0;
131 label maxProcFaces = 0;
133 for (
const int proci : Pstream::allProcs())
136 <<
"Processor " << proci <<
nl
137 <<
" Number of cells = " << globalCells.localSize(proci)
140 label nProcFaces = 0;
142 const labelList& nei = patchNeiProcNo[proci];
144 forAll(patchNeiProcNo[proci], i)
146 Info<<
" Number of faces shared with processor "
147 << patchNeiProcNo[proci][i] <<
" = " << patchSize[proci][i]
150 nProcFaces += patchSize[proci][i];
153 Info<<
" Number of processor patches = " << nei.
size() <<
nl
154 <<
" Number of processor faces = " << nProcFaces <<
nl
155 <<
" Number of boundary faces = "
156 << globalBoundaryFaces.localSize(proci) <<
endl;
158 maxProcCells =
max(maxProcCells, globalCells.localSize(proci));
159 totProcFaces += nProcFaces;
160 totProcPatches += nei.
size();
161 maxProcPatches =
max(maxProcPatches, nei.
size());
162 maxProcFaces =
max(maxProcFaces, nProcFaces);
167 scalar avgProcCells = scalar(globalCells.size())/Pstream::nProcs();
168 scalar avgProcPatches = scalar(totProcPatches)/Pstream::nProcs();
169 scalar avgProcFaces = scalar(totProcFaces)/Pstream::nProcs();
172 if (totProcPatches == 0)
176 if (totProcFaces == 0)
182 <<
"Number of processor faces = " << totProcFaces/2 <<
nl
183 <<
"Max number of cells = " << maxProcCells
184 <<
" (" << 100.0*(maxProcCells-avgProcCells)/avgProcCells
185 <<
"% above average " << avgProcCells <<
")" << nl
186 <<
"Max number of processor patches = " << maxProcPatches
187 <<
" (" << 100.0*(maxProcPatches-avgProcPatches)/avgProcPatches
188 <<
"% above average " << avgProcPatches <<
")" <<
nl
189 <<
"Max number of faces between processors = " << maxProcFaces
190 <<
" (" << 100.0*(maxProcFaces-avgProcFaces)/avgProcFaces
191 <<
"% above average " << avgProcFaces <<
")" << nl
205 return i*nCells[1]*nCells[2]+j*nCells[2]+
k;
235 if (nCells[0]+nCells[1]+nCells[2] > 0)
237 points.setSize((nCells[0]+1)*(nCells[1]+1)*(nCells[2]+1));
240 for (label i = 0; i <= nCells[0]; i++)
242 for (label j = 0; j <= nCells[1]; j++)
244 for (label
k = 0;
k <= nCells[2];
k++)
247 pt.
x() += i*cellSize[0];
248 pt.
y() += j*cellSize[1];
249 pt.
z() +=
k*cellSize[2];
250 points[vtxLabel(nCells, i, j,
k)] = pt;
262 for (label i = 0; i < nCells[0]; i++)
264 for (label j = 0; j < nCells[1]; j++)
266 for (label
k = 0;
k < nCells[2];
k++)
268 hexPoints[0] = vtxLabel(nCells, i, j,
k);
269 hexPoints[1] = vtxLabel(nCells, i+1, j,
k);
270 hexPoints[2] = vtxLabel(nCells, i+1, j+1,
k);
271 hexPoints[3] = vtxLabel(nCells, i, j+1,
k);
272 hexPoints[4] = vtxLabel(nCells, i, j,
k+1);
273 hexPoints[5] = vtxLabel(nCells, i+1, j,
k+1);
274 hexPoints[6] = vtxLabel(nCells, i+1, j+1,
k+1);
275 hexPoints[7] = vtxLabel(nCells, i, j+1,
k+1);
319 searchableSurfacesQueries::findNearest
338 forAll(nearestSurfaces, i)
340 if (nearestSurfaces[i] == surfI)
342 surfPoints.append(
points[i]);
343 surfIndices.append(i);
348 label geomI = surfaces[surfI];
350 geometry[geomI].getVolumeType(surfPoints, volType);
355 label pointi = surfIndices[i];
356 scalar dist =
mag(
points[pointi] - nearest[pointi].hitPoint());
360 if (vT == volumeType::OUTSIDE)
364 else if (vT == volumeType::INSIDE)
371 <<
"getVolumeType failure, neither INSIDE or OUTSIDE"
384int main(
int argc,
char *argv[])
388 "Generate foamyHexMesh-consistent representation of surfaces"
390 argList::addBoolOption
393 "Write the resulting mesh and distance fields"
399 "The merge distance relative to the bounding box size (default 1e-6)"
402 argList::noFunctionObjects();
407 const bool writeMesh =
args.
found(
"writeMesh");
411 Info<<
"Writing resulting mesh and cellDistance, pointDistance fields."
423 IOobject::MUST_READ_IF_MODIFIED,
433 "cvSearchableSurfaces",
440 foamyHexMeshDict.subDict(
"geometry"),
441 foamyHexMeshDict.getOrDefault(
"singleRegionName",
true)
451 foamyHexMeshDict.subDict(
"surfaceConformation")
457 foamyHexMeshDict.subDict(
"motionControl"),
466 const treeBoundBox& bb = geometryToConformTo.globalBounds();
470 vector nScalarCells = span/cellShapeControls.defaultCellSize();
476 label(nScalarCells.x())+2,
477 label(nScalarCells.y())+2,
478 label(nScalarCells.z())+2
487 Info<<
"Generating initial hex mesh with" <<
nl
488 <<
" bounding box : " << bb <<
nl
489 <<
" nCells : " << nCells <<
nl
490 <<
" cellSize : " << cellSize <<
nl
499 polyMesh::defaultRegion,
512 Info<<
"Writing initial hex mesh to " <<
meshPtr().instance() <<
nl
518 if (Pstream::parRun())
528 labelList decomp = decompositionModel::New
532 ).decomposer().decompose(
mesh,
mesh.cellCentres());
545 Info<<
"Wanted distribution:"
546 << distributor.countCells(decomp) <<
nl <<
endl;
556 Info<<
"Writing redistributed mesh" <<
nl <<
endl;
561 Info<<
"Refining background mesh according to cell size specification" <<
nl
565 foamyHexMeshDict.
subDict(
"backgroundMeshDecomposition");
579 backgroundMesh.mesh().
write();
582 const scalar tolDim = getMergeDistance
586 backgroundMesh.mesh().bounds()
595 const fvMesh& fvm = backgroundMesh.mesh();
613 const labelList& surfaces = geometryToConformTo.surfaces();
619 const labelList& cellLevel = backgroundMesh.cellLevel();
624 distSqr[celli] =
magSqr(cellSize)/
pow(2, cellLevel[celli]);
629 cellDistance.primitiveFieldRef() = signedDistance
639 cellDistance.boundaryFieldRef();
646 fld.patch().patchInternalField(distSqr)
648 fld = signedDistance(patchDistSqr, cc, geometry, surfaces);
656 cellDistance.write();
681 label ownDistSqr = distSqr[own];
686 pointDistSqr[
f[fp]] =
max(pointDistSqr[
f[fp]], ownDistSqr);
689 syncTools::syncPointList
697 pointDistance.primitiveFieldRef() = signedDistance
707 pointDistance.write();
722 isoFaces[i] = iso[i];
724 isoPoints = iso.points();
731 PatchTools::gatherAndMerge
740 if (Pstream::master())
vtk::internalMeshWriter writer(topoMesh, topoCells, vtk::formatType::INLINE_ASCII, runTime.path()/"blockTopology")
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))
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
void setSize(const label n)
Alias for resize()
virtual bool write(const token &tok)=0
Write token to stream or otherwise handle it.
A list of faces which address into the list of points.
A List obtained as a section of another List.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static word timeName(const scalar t, const int precision=precision_)
void size(const label n)
Older name for setAddressableSize.
Templated 3D Vector derived from VectorSpace adding construction from 3 components,...
const Cmpt & z() const
Access to the vector z component.
const Cmpt & y() const
Access to the vector y component.
const Cmpt & x() const
Access to the vector x component.
Extract command arguments and options from the supplied argc and argv parameters.
bool found(const word &optName) const
Return true if the named option is found.
T getOrDefault(const word &optName, const T &deflt) const
Get a value from the named option if present, or return default.
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
A bounding box defined in terms of min/max extrema points.
const point & min() const
Minimum describing the bounding box.
scalar mag() const
The magnitude of the bounding box span.
vector span() const
The bounding box span (from minimum to maximum)
Maps a geometry to a set of cell primitives.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
A face is a list of labels corresponding to mesh vertices.
A class for handling file names.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Mesh data needed to do the Finite Volume discretisation.
const volVectorField & C() const
Return cell centres as volVectorField.
const Time & time() const
Return the top-level database.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
A surface formed by the iso value. After "Polygonising A Scalar Field Using Tetrahedrons",...
Mesh consisting of general polyhedral cells.
virtual const faceList & faces() const
Return raw faces.
virtual const labelList & faceOwner() const
Return face owner.
virtual const pointField & points() const
Return raw points.
label nInternalFaces() const noexcept
Number of internal faces.
label nPoints() const noexcept
Number of mesh points.
Neighbour processor patch.
int neighbProcNo() const
Return neighbour processor number.
Container for searchableSurfaces. The collection is specified as a dictionary. For example,...
A class for managing temporary objects.
Standard boundBox with extra functionality for use in octree.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
An enumeration wrapper for classification of a location as being inside/outside of a volume.
Write faces/points (optionally with fields) as a vtp file or a legacy vtk file.
A class for handling words, derived from Foam::string.
Foam::autoPtr< Foam::dynamicFvMesh > meshPtr
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
IOobject io("surfaceFilmProperties", mesh.time().constant(), mesh, IOobject::READ_IF_PRESENT, IOobject::NO_WRITE, false)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
constexpr char nl
The newline '\n' character (0x0a)
wordList patchTypes(nPatches)
wordList patchNames(nPatches)
Foam::argList args(argc, argv)
#define forAll(list, i)
Loop across all elements in list.