48 displacementLayeredMotionMotionSolver,
54 displacementMotionSolver,
55 displacementLayeredMotionMotionSolver,
64 Foam::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];
83 void Foam::displacementLayeredMotionMotionSolver::calcZoneMask
85 const label cellZoneI,
91 isZoneEdge.resize(
mesh().nEdges());
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;
137 void 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();
229 Foam::displacementLayeredMotionMotionSolver::faceZoneEvaluate
233 const dictionary&
dict,
234 const PtrList<pointVectorField>& patchDisp,
239 auto&
fld = tfld.ref();
243 if (
type ==
"fixedValue")
245 fld.assign(
"value",
dict, meshPoints.size());
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 293 void 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()
359 pointDisplacement_.correctBoundaryConditions();
362 for (
const entry& dEntry : patchesDict)
364 const word& faceZoneName = dEntry.keyword();
365 const dictionary& faceZoneDict = dEntry.dict();
369 const labelList& fzMeshPoints = fz().meshPoints();
370 DynamicList<label> meshPoints(fzMeshPoints.size());
373 if (isZonePoint[fzMeshPoints[i]])
375 meshPoints.append(fzMeshPoints[i]);
380 tmp<vectorField> tseed = faceZoneEvaluate
390 <<
"For cellZone:" << cellZoneI
391 <<
" for faceZone:" << fz.name()
392 <<
" nPoints:" << tseed().size()
393 <<
" have patchField:" 394 <<
" max:" <<
gMax(tseed())
395 <<
" min:" <<
gMin(tseed())
415 patchDisp[patchi].correctBoundaryConditions();
431 mesh().cellZones()[cellZoneI].
name(),
440 for (
const label pointi : isZonePoint)
442 const scalar d1 = patchDist[0][pointi];
443 const scalar d2 = patchDist[1][pointi];
457 const word interpolationScheme(zoneDict.get<word>(
"interpolationScheme"));
459 if (interpolationScheme ==
"oneSided")
461 for (
const label pointi : isZonePoint)
463 pointDisplacement_[pointi] = patchDisp[0][pointi];
466 else if (interpolationScheme ==
"linear")
468 for (
const label pointi : isZonePoint)
470 const scalar d1 = patchDist[0][pointi];
471 const scalar d2 = patchDist[1][pointi];
472 const scalar
s = d1/(d1 + d2 + VSMALL);
474 const vector& pd1 = patchDisp[0][pointi];
475 const vector& pd2 = patchDisp[1][pointi];
477 pointDisplacement_[pointi] = (1 -
s)*pd1 +
s*pd2;
480 else if (interpolationScheme ==
"cylindrical")
482 const coordSystem::cylindrical& cs =
483 this->getCylindrical(cellZoneI, zoneDict);
489 <<
"cylindrical interpolation not yet available" <<
nl 490 <<
"coordinate system " << cs <<
nl 496 <<
"Invalid interpolationScheme: " << interpolationScheme
497 <<
". Valid schemes: (oneSided linear cylindrical)" <<
nl 505 Foam::displacementLayeredMotionMotionSolver::
506 displacementLayeredMotionMotionSolver
516 Foam::displacementLayeredMotionMotionSolver::
517 displacementLayeredMotionMotionSolver
536 points0() + pointDisplacement_.primitiveField()
549 pointDisplacement_.boundaryFieldRef().updateCoeffs();
552 for (
const entry& dEntry : coeffDict().subDict(
"regions"))
554 const word& cellZoneName = dEntry.keyword();
559 Info<<
"solving for zone: " << cellZoneName <<
endl;
564 <<
"Cannot find cellZone " << cellZoneName
569 cellZoneSolve(zoneI, regionDict);
586 const vectorField displacement(this->newPoints() - points0_);
590 const label oldPointi = mpm.
pointMap()[pointi];
596 if ((masterPointi != pointi))
602 points0_[pointi] -= displacement[pointi];
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
Field< label > labelField
Specialisation of Field<T> for label.
Virtual base class for displacement motion solver.
A cylindrical coordinate system (r-theta-z). The coordinate system angle theta is always in radians...
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Type gMin(const FieldField< Field, Type > &f)
static const pointMesh & New(const polyMesh &mesh, Args &&... args)
Get existing or create MeshObject registered with typeName.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
std::enable_if< std::is_same< bool, TypeT >::value, bool >::type set(const label i, bool val=true)
A bitSet::set() method for a list of bool.
constexpr char nl
The newline '\n' character (0x0a)
void constrainDisplacement(pointVectorField &displacement, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints),.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
GeometricField< vector, pointPatchField, pointMesh > pointVectorField
scalar distance(const vector &p1, const vector &p2)
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a T. FatalIOError if not found, or if the number of tokens is incorrect.
const Time & time() const
Return the top-level database.
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.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Macros for easy insertion into run-time selection tables.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
#define forAll(list, i)
Loop across all elements in list.
Application of (multi-)patch point constraints.
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
static word scopedName(const std::string &scope, const word &name)
Create scope:name or scope_name string.
vectorField pointField
pointField is a vectorField.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
label findZoneID(const word &zoneName) const
Find zone index by name, return -1 if not found.
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
const labelList & pointMap() const noexcept
Old point map.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
const labelList & reversePointMap() const noexcept
Reverse point map.
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
virtual void updateMesh(const mapPolyMesh &mpm)
Update topology.
#define DebugInfo
Report an information message using Foam::Info.
static tmp< GeometricField< scalar, pointPatchField, pointMesh > > New(const word &name, IOobjectOption::registerOption regOpt, const Mesh &mesh, const dimensionSet &dims, const word &patchFieldType=pointPatchField< scalar >::calculatedType())
Return tmp field (NO_READ, NO_WRITE) from name, mesh, dimensions and patch type. [Takes current timeN...
static word timeName(const scalar t, const int precision=precision_)
Return a time name for the given scalar time value formatted with the given precision.
int debug
Static debugging option.
Type gMax(const FieldField< Field, Type > &f)
const faceZoneMesh & faceZones() const noexcept
Return face zone mesh.
defineTypeNameAndDebug(combustionModel, 0)
const word typeName("volScalarField::Internal")
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;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
pointField points0(pointIOField(IOobject("points", mesh.time().constant(), polyMesh::meshSubDir, mesh, IOobject::MUST_READ, IOobject::NO_WRITE, IOobject::NO_REGISTER)))
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Type gAverage(const FieldField< Field, Type > &f)
virtual void solve()
Solve for motion.
const cellZoneMesh & cellZones() const noexcept
Return cell zone mesh.
wordList names() const
A list of the zone names.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< vector > vectorField
Specialisation of Field<T> for vector.
Mesh consisting of general polyhedral cells.
List< label > labelList
A List of labels.
A class for managing temporary objects.
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))
A primitive field of type <T> with automated input and output.
Do not request registration (bool: false)
A keyword and a list of tokens is an 'entry'.
addToRunTimeSelectionTable(functionObject, pointHistory, dictionary)
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)