rigidBodyMeshMotionSolver.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2016-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 #include "polyMesh.H"
29 #include "pointPatchDist.H"
30 #include "pointConstraints.H"
32 #include "forces.H"
33 #include "OneConstant.H"
34 #include "mathematicalConstants.H"
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40  defineTypeNameAndDebug(rigidBodyMeshMotionSolver, 0);
41 
43  (
44  motionSolver,
45  rigidBodyMeshMotionSolver,
46  dictionary
47  );
48 }
49 
50 
51 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
52 
53 Foam::rigidBodyMeshMotionSolver::bodyMesh::bodyMesh
54 (
55  const polyMesh& mesh,
56  const word& name,
57  const label bodyID,
58  const dictionary& dict
59 )
60 :
61  name_(name),
62  bodyID_(bodyID),
63  patches_(wordReList(dict.lookup("patches"))),
64  patchSet_(mesh.boundaryMesh().patchSet(patches_))
65 {}
66 
67 
69 (
70  const polyMesh& mesh,
71  const dictionary& dict
72 )
73 :
74  motionSolver(mesh, dict, typeName),
76  (
77  coeffDict(),
78  IOobject
79  (
80  "rigidBodyMotionState",
81  mesh.time().timeName(),
82  "uniform",
83  mesh
84  ).typeHeaderOk<IOdictionary>(true)
85  ? IOdictionary
86  (
87  IOobject
88  (
89  "rigidBodyMotionState",
90  mesh.time().timeName(),
91  "uniform",
92  mesh,
95  false
96  )
97  )
98  : coeffDict()
99  ),
100  test_(coeffDict().lookupOrDefault<Switch>("test", false)),
101  rhoInf_(1.0),
102  rhoName_(coeffDict().lookupOrDefault<word>("rho", "rho")),
103  ramp_(nullptr),
104  curTimeIndex_(-1),
105  meshSolverPtr_
106  (
108  (
109  mesh,
111  (
112  IOobject
113  (
114  "rigidBodyMotionSolver:meshSolver",
115  mesh.time().constant(),
116  mesh
117  ),
118  coeffDict().subDict("meshSolver")
119  )
120  )
121  ),
122  meshSolver_(refCast<displacementMotionSolver>(meshSolverPtr_()))
123 {
124  if (rhoName_ == "rhoInf")
125  {
126  rhoInf_ = coeffDict().lookup<scalar>("rhoInf");
127  }
128 
129  if (coeffDict().found("ramp"))
130  {
131  ramp_ = Function1<scalar>::New("ramp", coeffDict());
132  }
133  else
134  {
135  ramp_ = new Function1s::OneConstant<scalar>("ramp");
136  }
137 
138  const dictionary& bodiesDict = coeffDict().subDict("bodies");
139 
140  forAllConstIter(IDLList<entry>, bodiesDict, iter)
141  {
142  const dictionary& bodyDict = iter().dict();
143 
144  if (bodyDict.found("patches"))
145  {
146  const label bodyID = this->bodyID(iter().keyword());
147 
148  if (bodyID == -1)
149  {
151  << "Body " << iter().keyword()
152  << " has been merged with another body"
153  " and cannot be assigned a set of patches"
154  << exit(FatalError);
155  }
156 
157  bodyMeshes_.append
158  (
159  new bodyMesh
160  (
161  mesh,
162  iter().keyword(),
163  bodyID,
164  bodyDict
165  )
166  );
167  }
168  }
169 }
170 
171 
172 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
173 
175 {}
176 
177 
178 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
179 
182 {
183  return meshSolverPtr_->curPoints();
184 }
185 
186 
188 {
189  const Time& t = mesh().time();
190 
191  if (mesh().nPoints() != meshSolver_.points0().size())
192  {
194  << "The number of points in the mesh seems to have changed." << endl
195  << "In constant/polyMesh there are " << meshSolver_.points0().size()
196  << " points; in the current mesh there are " << mesh().nPoints()
197  << " points." << exit(FatalError);
198  }
199 
200  // Store the motion state at the beginning of the time-step
201  if (curTimeIndex_ != t.timeIndex())
202  {
203  newTime();
204  curTimeIndex_ = t.timeIndex();
205  }
206 
207  const scalar ramp = ramp_->value(t.value());
208 
210  {
211  g() =
212  ramp
214  }
215 
216  if (test_)
217  {
218  label nIter(coeffDict().lookup<label>("nIter"));
219 
220  for (label i=0; i<nIter; i++)
221  {
223  (
224  t.value(),
225  t.deltaTValue(),
226  scalarField(nDoF(), Zero),
228  );
229  }
230  }
231  else
232  {
234 
235  forAll(bodyMeshes_, bi)
236  {
237  const label bodyID = bodyMeshes_[bi].bodyID_;
238 
239  dictionary forcesDict;
240  forcesDict.add("type", functionObjects::forces::typeName);
241  forcesDict.add("patches", bodyMeshes_[bi].patches_);
242  forcesDict.add("rhoInf", rhoInf_);
243  forcesDict.add("rho", rhoName_);
244  forcesDict.add("CofR", vector::zero);
245 
246  functionObjects::forces f("forces", t, forcesDict);
247  f.calcForcesMoment();
248 
249  fx[bodyID] = ramp*spatialVector(f.momentEff(), f.forceEff());
250  }
251 
253  (
254  t.value(),
255  t.deltaTValue(),
256  scalarField(nDoF(), Zero),
257  fx
258  );
259  }
260 
261  if (Pstream::master() && report())
262  {
263  forAll(bodyMeshes_, bi)
264  {
265  status(bodyMeshes_[bi].bodyID_);
266  }
267  }
268 
269  // Update the displacements
270  forAll(bodyMeshes_, bi)
271  {
272  forAllConstIter(labelHashSet, bodyMeshes_[bi].patchSet_, iter)
273  {
274  label patchi = iter.key();
275 
276  pointField patchPoints0
277  (
278  meshSolver_.pointDisplacement().boundaryField()[patchi]
279  .patchInternalField(meshSolver_.points0())
280  );
281 
282  meshSolver_.pointDisplacement().boundaryFieldRef()[patchi] ==
283  (
285  (
286  bodyMeshes_[bi].bodyID_,
287  patchPoints0
288  ) - patchPoints0
289  )();
290  }
291  }
292 
293  meshSolverPtr_->solve();
294 }
295 
296 
298 {
299  meshSolverPtr_->movePoints(points);
300 }
301 
302 
304 {
305  meshSolverPtr_->updateMesh(mpm);
306 }
307 
308 
310 {
312  (
313  IOobject
314  (
315  "rigidBodyMotionState",
316  mesh().time().timeName(),
317  "uniform",
318  mesh(),
321  false
322  )
323  );
324 
325  state().write(dict);
326 
327  return
328  dict.regIOobject::writeObject
329  (
332  mesh().time().writeCompression(),
333  true
334  )
335  && motionSolver::write();
336 }
337 
338 
339 // ************************************************************************* //
void newTime()
Store the motion state at the beginning of the time-step.
Template class for intrusive linked lists.
Definition: ILList.H:50
virtual vector momentEff() const
Return the total moment.
Definition: forces.C:1025
void append(T *)
Append an element at the end of the list.
Definition: PtrListI.H:39
dictionary dict
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:643
Virtual base class for displacement motion solver.
virtual bool write() const
Write motion state information for restart.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const Boundary & boundaryField() const
Return const-reference to the boundary field.
rigidBodyMeshMotionSolver(const polyMesh &, const dictionary &dict)
Construct from polyMesh and dictionary.
tmp< pointField > transformPoints(const label bodyID, const pointField &initialPoints) const
Transform the given initial pointField of the specified body.
To & refCast(From &r)
Reference type cast template function.
Definition: typeInfo.H:106
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
void status(const label bodyID) const
Report the status of the motion of the given body.
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Definition: Switch.H:60
static bool master(const label communicator=0)
Am I the master process.
Definition: UPstream.H:423
bool foundObject(const word &name) const
Is the named Type found?
virtual bool write() const
Optionally write motion state information for restart.
Definition: motionSolver.C:151
static autoPtr< motionSolver > New(const polyMesh &, const dictionary &)
Select constructed from polyMesh and dictionary.
Definition: motionSolver.C:62
Virtual base class for mesh motion solver.
Definition: motionSolver.H:55
virtual vector forceEff() const
Return the total force.
Definition: forces.C:1019
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:636
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
label nBodies() const
Return the number of bodies in the model (bodies().size())
const dictionary & coeffDict() const
Const access to the coefficients dictionary.
Definition: motionSolver.H:129
const rigidBodyModelState & state() const
Return the motion state.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual void movePoints(const pointField &)
Update local data for geometry changes.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
Macros for easy insertion into run-time selection tables.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:1133
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:982
Templated function that returns the corresponding 1 (one).
Definition: OneConstant.H:56
pointVectorField & pointDisplacement()
Return reference to the point motion displacement field.
Six degree of freedom motion for a rigid body.
stressControl lookup("compactNormalStress") >> compactNormalStress
const pointField & points
virtual Type value(const scalar x) const =0
Return value as a function of scalar x.
Pre-declare SubField and related Field type.
Definition: Field.H:56
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
SpatialVector< scalar > spatialVector
SpatialVector of scalars.
Definition: spatialVector.H:47
void write(dictionary &dict) const
Write to dictionary.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
const Type & value() const
Return const reference to value.
word timeName
Definition: getTimeIndex.H:3
virtual void calcForcesMoment()
Calculate the forces and moments.
Definition: forces.C:884
static const zero Zero
Definition: zero.H:97
virtual void updateMesh(const mapPolyMesh &)
Update local data for topology changes.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
void solve(const scalar t, const scalar deltaT, const scalarField &tau, const Field< spatialVector > &fx)
Integrate velocities, orientation and position.
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
label nDoF() const
Return the number of degrees of freedom of the model.
const vector & g() const
Return the acceleration due to gravity.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
label patchi
static const versionNumber currentVersion
Current version number.
Definition: IOstream.H:203
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
label nPoints() const
label bodyID(const word &name) const
Return the ID of the body with the given name.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
List< wordRe > wordReList
A List of wordRe (word or regular expression)
Definition: wordReList.H:50
Calculates the forces and moments by integrating the pressure and skin-friction forces over a given l...
Definition: forces.H:204
const polyMesh & mesh() const
Return reference to mesh.
Definition: motionSolver.H:123
A class for managing temporary objects.
Definition: PtrList.H:53
pointField & points0()
Return reference to the reference field.
virtual void solve()
Solve for motion.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
motionSolver(const polyMesh &mesh, const dictionary &, const word &type)
Construct from polyMesh and dictionary and type.
Definition: motionSolver.C:41
bool report() const
Return the report Switch.
Namespace for OpenFOAM.
static autoPtr< Function1< Type > > New(const word &name, const dictionary &dict)
Selector.
Definition: Function1New.C:32
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:844