Go to the documentation of this file.
47 rigidBodyMeshMotionSolver,
55 Foam::rigidBodyMeshMotionSolver::bodyMesh::bodyMesh
60 const dictionary&
dict
65 patches_(
dict.
get<wordRes>(
"patches")),
66 patchSet_(
mesh.boundaryMesh().patchSet(patches_))
70 Foam::rigidBodyMeshMotionSolver::rigidBodyMeshMotionSolver
83 "rigidBodyMotionState",
87 ).typeHeaderOk<IOdictionary>(
true)
92 "rigidBodyMotionState",
103 test_(coeffDict().getOrDefault(
"test", false)),
105 rhoName_(coeffDict().getOrDefault<
word>(
"rho",
"rho")),
116 "rigidBodyMotionSolver:meshSolver",
120 coeffDict().subDict(
"meshSolver")
126 if (rhoName_ ==
"rhoInf")
128 coeffDict().
readEntry(
"rhoInf", rhoInf_);
133 for (
const entry& dEntry : bodiesDict)
135 const keyType& bodyName = dEntry.keyword();
138 if (bodyDict.
found(
"patches"))
140 const label bodyID = model_.bodyID(bodyName);
145 <<
"Body " << bodyName
146 <<
" has been merged with another body"
147 " and cannot be assigned a set of patches"
171 return meshSolverPtr_->curPoints();
179 if (
mesh().
nPoints() != meshSolver_.points0().size())
182 <<
"The number of points in the mesh seems to have changed." <<
endl
183 <<
"In constant/polyMesh there are " << meshSolver_.points0().size()
184 <<
" points; in the current mesh there are " <<
mesh().
nPoints()
189 if (curTimeIndex_ != this->db().time().
timeIndex())
192 curTimeIndex_ = this->db().time().timeIndex();
202 const label nIter(coeffDict().get<label>(
"nIter"));
204 for (label i=0; i<nIter; i++)
221 const label bodyID = bodyMeshes_[bi].bodyID_;
224 forcesDict.
add(
"type", functionObjects::forces::typeName);
225 forcesDict.
add(
"patches", bodyMeshes_[bi].patches_);
226 forcesDict.
add(
"rhoInf", rhoInf_);
227 forcesDict.
add(
"rho", rhoName_);
231 f.calcForcesMoment();
249 model_.status(bodyMeshes_[bi].bodyID_);
256 for (
const label patchi : bodyMeshes_[bi].patchSet_)
260 meshSolver_.pointDisplacement().boundaryField()[patchi]
261 .patchInternalField(meshSolver_.points0())
264 meshSolver_.pointDisplacement().boundaryFieldRef()[patchi] ==
266 model_.transformPoints
268 bodyMeshes_[bi].bodyID_,
275 meshSolverPtr_->solve();
292 "rigidBodyMotionState",
302 model_.state().write(
dict);
303 return dict.regIOobject::writeObject(streamOpt, valid);
311 model_.read(coeffDict());
322 meshSolverPtr_->movePoints(
points);
328 meshSolverPtr_->updateMesh(mpm);
addToRunTimeSelectionTable(decompositionMethod, kahipDecomp, dictionary)
A keyword and a list of tokens is an 'entry'.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
A class for handling words, derived from Foam::string.
Calculates the forces and moments by integrating the pressure and skin-friction forces over a given l...
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
virtual void solve()
Solve for motion.
static word timeName(const scalar t, const int precision=precision_)
streamFormat format() const noexcept
Get the current stream format.
static bool master(const label communicator=worldComm)
Am I the master process.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const Type & value() const
Return const reference to value.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
Mesh consisting of general polyhedral cells.
label nPoints() const noexcept
Number of mesh points.
scalar deltaTValue() const noexcept
Return time step value.
#define forAll(list, i)
Loop across all elements in list.
A class for handling keywords in dictionaries.
virtual void movePoints(const pointField &)
Update local data for geometry changes.
virtual bool read()
Read dynamicMeshDict dictionary.
Generic templated field type.
const dictionary & subDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary.
The IOstreamOption is a simple container for options an IOstream can normally have.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
const Type & lookupObject(const word &name, const bool recursive=false) const
virtual bool read()
Read dynamicMeshDict dictionary.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Virtual base class for displacement motion solver.
Macros for easy insertion into run-time selection tables.
To & refCast(From &r)
Reference type cast template function.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
errorManipArg< error, int > exit(error &err, const int errNo=1)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write state using stream options.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Virtual base class for mesh motion solver.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const Time & time() const
Return the top-level database.
static const Vector< scalar > zero
entry * add(entry *entryPtr, bool mergeEntry=false)
Add a new entry.
SpatialVector< scalar > spatialVector
SpatialVector of scalars.
constant condensation/saturation model.
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
defineTypeNameAndDebug(combustionModel, 0)