sixDoFRigidBodyMotionSolver.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2015 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 
52 Foam::sixDoFRigidBodyMotionSolver::sixDoFRigidBodyMotionSolver
53 (
54  const polyMesh& mesh,
55  const IOdictionary& dict
56 )
57 :
58  displacementMotionSolver(mesh, dict, typeName),
59  motion_
60  (
61  coeffDict(),
62  IOobject
63  (
64  "sixDoFRigidBodyMotionState",
65  mesh.time().timeName(),
66  "uniform",
67  mesh
68  ).headerOk()
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_(readScalar(coeffDict().lookup("innerDistance"))),
87  do_(readScalar(coeffDict().lookup("outerDistance"))),
88  test_(coeffDict().lookupOrDefault<Switch>("test", false)),
89  rhoInf_(1.0),
90  rhoName_(coeffDict().lookupOrDefault<word>("rhoName", "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),
103  dimensionedScalar("zero", dimless, 0.0)
104  ),
105  curTimeIndex_(-1)
106 {
107  if (rhoName_ == "rhoInf")
108  {
109  rhoInf_ = readScalar(coeffDict().lookup("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_.internalField() =
122  min
123  (
124  max
125  (
126  (do_ - pDist.internalField())/(do_ - di_),
127  scalar(0)
128  ),
129  scalar(1)
130  );
131 
132  // Convert the scale function to a cosine
133  scale_.internalField() =
134  min
135  (
136  max
137  (
138  0.5
139  - 0.5
140  *cos(scale_.internalField()
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_.internalField();
165 }
166 
167 
169 {
170  const Time& t = mesh().time();
171 
172  if (mesh().nPoints() != points0().size())
173  {
175  (
176  "sixDoFRigidBodyMotionSolver::curPoints() const"
177  ) << "The number of points in the mesh seems to have changed." << endl
178  << "In constant/polyMesh there are " << points0().size()
179  << " points; in the current mesh there are " << mesh().nPoints()
180  << " points." << exit(FatalError);
181  }
182 
183  // Store the motion state at the beginning of the time-stepbool
184  bool firstIter = false;
185  if (curTimeIndex_ != this->db().time().timeIndex())
186  {
187  motion_.newTime();
188  curTimeIndex_ = this->db().time().timeIndex();
189  firstIter = true;
190  }
191 
193 
194  if (db().foundObject<uniformDimensionedVectorField>("g"))
195  {
196  g = db().lookupObject<uniformDimensionedVectorField>("g");
197  }
198  else if (coeffDict().found("g"))
199  {
200  coeffDict().lookup("g") >> g;
201  }
202 
203  // scalar ramp = min(max((this->db().time().value() - 5)/10, 0), 1);
204  scalar ramp = 1.0;
205 
206  if (test_)
207  {
208  motion_.update
209  (
210  firstIter,
211  ramp*(motion_.mass()*g.value()),
212  ramp*(motion_.mass()*(motion_.momentArm() ^ g.value())),
213  t.deltaTValue(),
214  t.deltaT0Value()
215  );
216  }
217  else
218  {
219  dictionary forcesDict;
220 
221  forcesDict.add("type", forces::typeName);
222  forcesDict.add("patches", patches_);
223  forcesDict.add("rhoInf", rhoInf_);
224  forcesDict.add("rhoName", rhoName_);
225  forcesDict.add("CofR", motion_.centreOfRotation());
226 
227  forces f("forces", db(), forcesDict);
228 
229  f.calcForcesMoment();
230 
231  motion_.update
232  (
233  firstIter,
234  ramp*(f.forceEff() + motion_.mass()*g.value()),
235  ramp
236  *(
237  f.momentEff()
238  + motion_.mass()*(motion_.momentArm() ^ g.value())
239  ),
240  t.deltaTValue(),
241  t.deltaT0Value()
242  );
243  }
244 
245  // Update the displacements
246  pointDisplacement_.internalField() =
247  motion_.transform(points0(), scale_) - points0();
248 
249  // Displacement has changed. Update boundary conditions
251  (
252  pointDisplacement_.mesh()
253  ).constrainDisplacement(pointDisplacement_);
254 }
255 
256 
258 (
262 ) const
263 {
265  (
266  IOobject
267  (
268  "sixDoFRigidBodyMotionState",
269  mesh().time().timeName(),
270  "uniform",
271  mesh(),
274  false
275  )
276  );
277 
278  motion_.state().write(dict);
279  return dict.regIOobject::write();
280 }
281 
282 
284 {
286  {
287  motion_.read(coeffDict());
288 
289  return true;
290  }
291  else
292  {
293  return false;
294  }
295 }
296 
297 
298 // ************************************************************************* //
void constrain(GeometricField< Type, pointPatchField, pointMesh > &pf, const bool overrideValue=false) const
Apply boundary conditions (single-patch constraints) and.
virtual void calcForcesMoment()
Calculate the forces and moments.
Definition: forces.C:809
labelList f(nPoints)
static const pointMesh & New(const polyMesh &mesh)
This function object calculates the forces and moments by integrating the pressure and skin-friction ...
Definition: forces.H:206
scalar deltaT0Value() const
Return old time step value.
Definition: TimeState.H:106
virtual bool read()
Read dynamicMeshDict dictionary.
label timeIndex
Definition: getTimeIndex.H:4
A class for handling words, derived from string.
Definition: word.H:59
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:741
Addressing for all faces on surface of mesh. Can either be read from polyMesh or from triSurface...
Definition: boundaryMesh.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
InternalField & internalField()
Return internal field.
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.
Definition: Switch.H:60
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
Namespace for OpenFOAM.
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
virtual tmp< pointField > curPoints() const
Return point location obtained from the current motion field.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp) const
Write state using given format, version and compression.
dictionary dict
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
stressControl lookup("compactNormalStress") >> compactNormalStress
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
scalar deltaTValue() const
Return time step value.
Definition: TimeState.H:100
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
dimensionedScalar cos(const dimensionedScalar &ds)
Macros for easy insertion into run-time selection tables.
virtual void solve()
Solve for motion.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:739
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Calculation of distance to nearest patch for all points.
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
Dimensioned<Type> registered with the database as a registered IOobject which has the functionality o...
virtual vector momentEff() const
Return the total moment.
Definition: forces.C:960
const dimensionSet dimAcceleration
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const dimensionedVector & g
error FatalError
static const Vector zero
Definition: Vector.H:80
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
const Time & time() const
Return time.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
Version number type.
Definition: IOstream.H:96
label nPoints
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
virtual vector forceEff() const
Return the total force.
Definition: forces.C:954
A class for managing temporary objects.
Definition: PtrList.H:118
virtual bool read()
Read dynamicMeshDict dictionary.
Definition: motionSolver.C:193
label nPoints() const
const Type & value() const
Return const reference to value.
defineTypeNameAndDebug(combustionModel, 0)
word timeName
Definition: getTimeIndex.H:3
Virtual base class for displacement motion solver.