sixDoFRigidBodyMotionSolver.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) 2013-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 "mathematicalConstants.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39  defineTypeNameAndDebug(sixDoFRigidBodyMotionSolver, 0);
40 
42  (
43  motionSolver,
44  sixDoFRigidBodyMotionSolver,
45  dictionary
46  );
47 }
48 
49 
50 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 
53 (
54  const polyMesh& mesh,
55  const dictionary& dict
56 )
57 :
58  displacementMotionSolver(mesh, dict, typeName),
60  (
61  coeffDict(),
62  IOobject
63  (
64  "sixDoFRigidBodyMotionState",
65  mesh.time().timeName(),
66  "uniform",
67  mesh
68  ).typeHeaderOk<IOdictionary>(true)
69  ? IOdictionary
70  (
71  IOobject
72  (
73  "sixDoFRigidBodyMotionState",
74  mesh.time().timeName(),
75  "uniform",
76  mesh,
79  false
80  )
81  )
82  : coeffDict()
83  ),
84  patches_(wordReList(coeffDict().lookup("patches"))),
85  patchSet_(mesh.boundaryMesh().patchSet(patches_)),
86  di_(coeffDict().lookup<scalar>("innerDistance")),
87  do_(coeffDict().lookup<scalar>("outerDistance")),
88  test_(coeffDict().lookupOrDefault<Switch>("test", false)),
89  rhoInf_(1.0),
90  rhoName_(coeffDict().lookupOrDefault<word>("rho", "rho")),
91  scale_
92  (
93  IOobject
94  (
95  "motionScale",
96  mesh.time().timeName(),
97  mesh,
98  IOobject::NO_READ,
99  IOobject::NO_WRITE,
100  false
101  ),
102  pointMesh::New(mesh),
104  ),
105  curTimeIndex_(-1)
106 {
107  if (rhoName_ == "rhoInf")
108  {
109  rhoInf_ = coeffDict().lookup<scalar>("rhoInf");
110  }
111 
112  // Calculate scaling factor everywhere
113  // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
114 
115  {
116  const pointMesh& pMesh = pointMesh::New(mesh);
117 
118  pointPatchDist pDist(pMesh, patchSet_, points0());
119 
120  // Scaling: 1 up to di then linear down to 0 at do away from patches
121  scale_.primitiveFieldRef() =
122  min
123  (
124  max
125  (
126  (do_ - pDist.primitiveField())/(do_ - di_),
127  scalar(0)
128  ),
129  scalar(1)
130  );
131 
132  // Convert the scale function to a cosine
133  scale_.primitiveFieldRef() =
134  min
135  (
136  max
137  (
138  0.5
139  - 0.5
140  *cos(scale_.primitiveField()
142  scalar(0)
143  ),
144  scalar(1)
145  );
146 
147  pointConstraints::New(pMesh).constrain(scale_);
148  scale_.write();
149  }
150 }
151 
152 
153 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
154 
156 {}
157 
158 
159 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
160 
163 {
164  return points0() + pointDisplacement_.primitiveField();
165 }
166 
167 
169 {
170  const Time& t = mesh().time();
171 
172  if (mesh().nPoints() != points0().size())
173  {
175  << "The number of points in the mesh seems to have changed." << endl
176  << "In constant/polyMesh there are " << points0().size()
177  << " points; in the current mesh there are " << mesh().nPoints()
178  << " points." << exit(FatalError);
179  }
180 
181  // Store the motion state at the beginning of the time-stepbool
182  bool firstIter = false;
183  if (curTimeIndex_ != t.timeIndex())
184  {
185  newTime();
186  curTimeIndex_ = t.timeIndex();
187  firstIter = true;
188  }
189 
191 
192  if (mesh().foundObject<uniformDimensionedVectorField>("g"))
193  {
195  }
196  else if (coeffDict().found("g"))
197  {
198  coeffDict().lookup("g") >> g;
199  }
200 
201  // scalar ramp = min(max((t.value() - 5)/10, 0), 1);
202  scalar ramp = 1.0;
203 
204  if (test_)
205  {
206  update
207  (
208  firstIter,
209  ramp*(mass()*g.value()),
210  ramp*(mass()*(momentArm() ^ g.value())),
211  t.deltaTValue(),
212  t.deltaT0Value()
213  );
214  }
215  else
216  {
217  dictionary forcesDict;
218 
219  forcesDict.add("type", functionObjects::forces::typeName);
220  forcesDict.add("patches", patches_);
221  forcesDict.add("rhoInf", rhoInf_);
222  forcesDict.add("rho", rhoName_);
223  forcesDict.add("CofR", centreOfRotation());
224 
225  functionObjects::forces f("forces", t, forcesDict);
226 
227  f.calcForcesMoment();
228 
229  update
230  (
231  firstIter,
232  ramp*(f.forceEff() + mass()*g.value()),
233  ramp
234  *(
235  f.momentEff()
236  + mass()*(momentArm() ^ g.value())
237  ),
238  t.deltaTValue(),
239  t.deltaT0Value()
240  );
241  }
242 
243  // Update the displacements
244  pointDisplacement_.primitiveFieldRef() =
245  transform(points0(), scale_) - points0();
246 
247  // Displacement has changed. Update boundary conditions
249  (
250  pointDisplacement_.mesh()
251  ).constrainDisplacement(pointDisplacement_);
252 }
253 
254 
256 {
258  (
259  IOobject
260  (
261  "sixDoFRigidBodyMotionState",
262  mesh().time().timeName(),
263  "uniform",
264  mesh(),
267  false
268  )
269  );
270 
271  state().write(dict);
272 
273  return
274  dict.regIOobject::writeObject
275  (
278  mesh().time().writeCompression(),
279  true
280  )
282 }
283 
284 
285 // ************************************************************************* //
Six degree of freedom motion for a rigid body.
virtual vector momentEff() const
Return the total moment.
Definition: forces.C:1025
virtual bool write() const
Write motion state information for restart.
dictionary dict
Virtual base class for displacement motion solver.
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
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
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
virtual bool write() const
Optionally write motion state information for restart.
Definition: motionSolver.C:151
const Internal::FieldType & primitiveField() const
Return a const-reference to the internal field.
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
const dimensionSet dimless
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
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
Macros for easy insertion into run-time selection tables.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:1133
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
dimensionedScalar cos(const dimensionedScalar &ds)
const dimensionSet dimAcceleration
scalar deltaT0Value() const
Return old time step value.
Definition: TimeStateI.H:47
A class for handling words, derived from string.
Definition: word.H:59
label nPoints
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
Calculation of distance to nearest patch for all points.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
sixDoFRigidBodyMotionSolver(const polyMesh &, const dictionary &dict)
Construct from polyMesh and dictionary.
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
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Time & time() const
Return time.
defineTypeNameAndDebug(combustionModel, 0)
labelList f(nPoints)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual void solve()
Solve for motion.
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
static const versionNumber currentVersion
Current version number.
Definition: IOstream.H:203
label nPoints() const
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
Calculates the forces and moments by integrating the pressure and skin-friction forces over a given l...
Definition: forces.H:204
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedVector & g
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void constrain(GeometricField< Type, pointPatchField, pointMesh > &pf, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints) and.
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:477