64Foam::displacementLayeredMotionMotionSolver::getCylindrical
66 const label cellZoneI,
67 const dictionary& zoneDict
70 auto iter = cylSystems_.cfind(cellZoneI);
77 cylSystems_.set(cellZoneI,
new coordSystem::cylindrical(zoneDict));
79 return *cylSystems_[cellZoneI];
83void Foam::displacementLayeredMotionMotionSolver::calcZoneMask
85 const label cellZoneI,
107 for (
const label celli : cz)
109 isZonePoint.
set(
mesh().cellPoints(celli));
110 isZoneEdge.set(
mesh().cellEdges(celli));
117 orEqOp<unsigned int>(),
125 orEqOp<unsigned int>(),
130 <<
"On cellZone " << cz.name() <<
" marked "
131 <<
returnReduce(isZonePoint.count(), sumOp<label>()) <<
" points and "
132 <<
returnReduce(isZoneEdge.count(), sumOp<label>()) <<
" edges" <<
nl;
137void Foam::displacementLayeredMotionMotionSolver::walkStructured
139 const label cellZoneI,
140 const bitSet& isZonePoint,
141 const bitSet& isZoneEdge,
149 List<pointEdgeStructuredWalk> seedInfo(seedPoints.size());
153 const label pointi = seedPoints[i];
155 seedInfo[i] = pointEdgeStructuredWalk
166 List<pointEdgeStructuredWalk> allPointInfo(
mesh().
nPoints());
171 for (
const label pointi : isZonePoint)
173 allPointInfo[pointi] = pointEdgeStructuredWalk
184 List<pointEdgeStructuredWalk> allEdgeInfo(
mesh().
nEdges());
187 for (
const label edgei : isZoneEdge)
189 allEdgeInfo[edgei] = pointEdgeStructuredWalk
199 PointEdgeWave<pointEdgeStructuredWalk> wallCalc
207 mesh().globalData().nTotalPoints()
211 for (
const label pointi : isZonePoint)
213 const auto& pointInfo = allPointInfo[pointi];
215 distance[pointi] = pointInfo.dist();
216 data[pointi] = pointInfo.data();
219 if (patchPoints.size())
221 patchPoints[pointi] = pointInfo.index();
229Foam::displacementLayeredMotionMotionSolver::faceZoneEvaluate
233 const dictionary&
dict,
234 const PtrList<pointVectorField>& patchDisp,
239 auto&
fld = tfld.ref();
243 if (
type ==
"fixedValue")
247 else if (
type ==
"timeVaryingUniformFixedValue")
249 interpolationTable<vector> timeSeries(
dict);
251 fld = timeSeries(
mesh().time().timeOutputValue());
253 else if (
type ==
"slip")
255 if ((patchi % 2) != 1)
258 <<
"FaceZone:" << fz.name()
264 else if (
type ==
"follow")
269 else if (
type ==
"uniformFollow")
273 const word patchName(
dict.
get<word>(
"patch"));
277 pointDisplacement_.boundaryField()[
patchID].patchInternalField()
284 <<
"Unknown faceZonePatch type " <<
type
285 <<
" for faceZone " << fz.name() <<
nl
293void Foam::displacementLayeredMotionMotionSolver::cellZoneSolve
295 const label cellZoneI,
296 const dictionary& zoneDict
299 bitSet isZonePoint, isZoneEdge;
300 calcZoneMask(cellZoneI, isZonePoint, isZoneEdge);
302 const dictionary& patchesDict = zoneDict.subDict(
"boundaryField");
304 if (patchesDict.size() != 2)
307 <<
"Two faceZones (patches) must be specified per cellZone. "
308 <<
" cellZone:" << cellZoneI
309 <<
" patches:" << patchesDict.toc()
313 PtrList<labelField> patchPoints(patchesDict.size());
314 PtrList<scalarField> patchDist(patchesDict.size());
315 PtrList<pointVectorField> patchDisp(patchesDict.size());
319 for (
const entry& dEntry : patchesDict)
321 const word& faceZoneName = dEntry.keyword();
326 <<
"Cannot find faceZone " << faceZoneName
343 mesh().cellZones()[cellZoneI].
name() +
"_" + fz.name(),
364 pointDisplacement_.correctBoundaryConditions();
367 for (
const entry& dEntry : patchesDict)
369 const word& faceZoneName = dEntry.keyword();
370 const dictionary& faceZoneDict = dEntry.dict();
374 const labelList& fzMeshPoints = fz().meshPoints();
375 DynamicList<label> meshPoints(fzMeshPoints.size());
378 if (isZonePoint[fzMeshPoints[i]])
380 meshPoints.append(fzMeshPoints[i]);
385 tmp<vectorField> tseed = faceZoneEvaluate
395 <<
"For cellZone:" << cellZoneI
396 <<
" for faceZone:" << fz.name()
397 <<
" nPoints:" << tseed().size()
398 <<
" have patchField:"
399 <<
" max:" <<
gMax(tseed())
400 <<
" min:" <<
gMin(tseed())
420 patchDisp[patchi].correctBoundaryConditions();
436 mesh().cellZones()[cellZoneI].
name() +
":distance",
447 for (
const label pointi : isZonePoint)
449 const scalar d1 = patchDist[0][pointi];
450 const scalar d2 = patchDist[1][pointi];
464 const word interpolationScheme(zoneDict.get<word>(
"interpolationScheme"));
466 if (interpolationScheme ==
"oneSided")
468 for (
const label pointi : isZonePoint)
470 pointDisplacement_[pointi] = patchDisp[0][pointi];
473 else if (interpolationScheme ==
"linear")
475 for (
const label pointi : isZonePoint)
477 const scalar d1 = patchDist[0][pointi];
478 const scalar d2 = patchDist[1][pointi];
479 const scalar
s = d1/(d1 + d2 + VSMALL);
481 const vector& pd1 = patchDisp[0][pointi];
482 const vector& pd2 = patchDisp[1][pointi];
484 pointDisplacement_[pointi] = (1 -
s)*pd1 +
s*pd2;
487 else if (interpolationScheme ==
"cylindrical")
489 const coordSystem::cylindrical& cs =
490 this->getCylindrical(cellZoneI, zoneDict);
496 <<
"cylindrical interpolation not yet available" <<
nl
497 <<
"coordinate system " << cs <<
nl
503 <<
"Invalid interpolationScheme: " << interpolationScheme
504 <<
". Valid schemes: (oneSided linear cylindrical)" <<
nl
512Foam::displacementLayeredMotionMotionSolver::
513displacementLayeredMotionMotionSolver
523Foam::displacementLayeredMotionMotionSolver::
524displacementLayeredMotionMotionSolver
543 points0() + pointDisplacement_.primitiveField()
556 pointDisplacement_.boundaryFieldRef().updateCoeffs();
559 for (
const entry& dEntry : coeffDict().subDict(
"regions"))
561 const word& cellZoneName = dEntry.keyword();
566 Info<<
"solving for zone: " << cellZoneName <<
endl;
571 <<
"Cannot find cellZone " << cellZoneName
576 cellZoneSolve(zoneI, regionDict);
593 const vectorField displacement(this->newPoints() - points0_);
597 const label oldPointi = mpm.
pointMap()[pointi];
603 if ((masterPointi != pointi))
609 points0_[pointi] -= displacement[pointi];
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
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 primitive field of type <T> with automated input and output.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
const T * set(const label i) const
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
wordList names() const
A list of the zone names.
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Mesh motion solver for an (multi-block) extruded fvMesh. Gets given the structure of the mesh blocks ...
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
virtual void solve()
Solve for motion.
Virtual base class for displacement motion solver.
A keyword and a list of tokens is an 'entry'.
const Time & time() const
Return the top-level database.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const labelList & pointMap() const
Old point map.
const labelList & reversePointMap() const
Reverse point map.
void updateMesh()
Update for new mesh topology.
Virtual base class for mesh motion solver.
static const complex max
complex (VGREAT,VGREAT)
Application of (multi-)patch point constraints.
void constrainDisplacement(pointVectorField &displacement, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints),.
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Mesh consisting of general polyhedral cells.
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
A class for managing temporary objects.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const labelList nEdges(UPstream::listGatherValues< label >(aMesh.nEdges()))
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))
#define DebugInfo
Report an information message using Foam::Info.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
scalar distance(const vector &p1, const vector &p2)
GeometricField< scalar, pointPatchField, pointMesh > pointScalarField
List< label > labelList
A List of labels.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
messageStream Info
Information stream (stdout output on master, null elsewhere)
GeometricField< vector, pointPatchField, pointMesh > pointVectorField
vectorField pointField
pointField is a vectorField.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Field< vector > vectorField
Specialisation of Field<T> for vector.
Field< label > labelField
Specialisation of Field<T> for label.
Type gAverage(const FieldField< Field, Type > &f)
static constexpr const zero Zero
Global zero (0)
Type gMin(const FieldField< Field, Type > &f)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
T returnReduce(const T &value, const BinaryOp &bop, const int tag=UPstream::msgType(), const label comm=UPstream::worldComm)
Reduce (copy) and return value.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Type gMax(const FieldField< Field, Type > &f)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, false)))