55Foam::rigidBodyMeshMotionSolver::bodyMesh::bodyMesh
60 const dictionary&
dict
65 patches_(
dict.get<wordRes>(
"patches")),
66 patchSet_(
mesh.boundaryMesh().patchSet(patches_))
83 "rigidBodyMotionState",
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")
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_;
225 forcesDict.
add(
"patches", bodyMeshes_[bi].patches_);
226 forcesDict.
add(
"rhoInf", rhoInf_);
227 forcesDict.
add(
"rho", rhoName_);
231 f.calcForcesMoments();
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);
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
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,...
The IOstreamOption is a simple container for options an IOstream can normally have.
streamFormat format() const noexcept
Get the current stream format.
void append(T *ptr)
Append an element to the end of the list.
label bodyID(const word &name) const
Return the ID of the body with the given name.
scalar deltaTValue() const noexcept
Return time step value.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
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.
bool found(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Search for an entry (const access) with the given keyword.
entry * add(entry *entryPtr, bool mergeEntry=false)
Add a new entry.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
const Type & value() const
Return const reference to value.
Virtual base class for displacement motion solver.
A keyword and a list of tokens is an 'entry'.
Computes forces and moments over a given list of patches by integrating pressure and viscous forces a...
const Time & time() const
Return the top-level database.
A class for handling keywords in dictionaries.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
void movePoints()
Update for new mesh geometry.
void updateMesh()
Update for new mesh topology.
Virtual base class for mesh motion solver.
const polyMesh & mesh() const
Return reference to mesh.
const dictionary & coeffDict() const
Const access to the coefficients dictionary.
virtual bool read()
Read dynamicMeshDict dictionary.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
const Type & lookupObject(const word &name, const bool recursive=false) const
constant condensation/saturation model.
Mesh consisting of general polyhedral cells.
label nPoints() const noexcept
Number of mesh points.
Rigid-body mesh motion solver for fvMesh.
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
virtual bool writeObject(IOstreamOption streamOpt, const bool valid) const
Write state using stream options.
virtual void solve()
Solve for motion.
virtual bool read()
Read dynamicMeshDict dictionary.
splitCell * master() const
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
To & refCast(From &r)
Reference type cast template function.
Ostream & endl(Ostream &os)
Add newline and flush stream.
SpatialVector< scalar > spatialVector
SpatialVector of scalars.
static constexpr const zero Zero
Global zero (0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
errorManipArg< error, int > exit(error &err, const int errNo=1)
#define forAll(list, i)
Loop across all elements in list.
static const char *const typeName
The type name used in ensight case files.
A non-counting (dummy) refCount.