sixDoFRigidBodyDisplacementPointPatchVectorField.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) 2011-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 
27 #include "pointPatchFields.H"
29 #include "Time.H"
30 #include "fvMesh.H"
31 #include "volFields.H"
33 #include "forces.H"
34 
35 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
36 
37 namespace Foam
38 {
39 
40 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
41 
44 (
45  const pointPatch& p,
47 )
48 :
50  motion_(),
51  initialPoints_(p.localPoints()),
52  rhoInf_(1.0),
53  rhoName_("rho"),
54  lookupGravity_(-1),
55  g_(vector::zero),
56  curTimeIndex_(-1)
57 {}
58 
59 
62 (
63  const pointPatch& p,
65  const dictionary& dict
66 )
67 :
69  motion_(dict, dict),
70  rhoInf_(1.0),
71  rhoName_(dict.lookupOrDefault<word>("rhoName", "rho")),
72  lookupGravity_(-1),
73  g_(vector::zero),
74  curTimeIndex_(-1)
75 {
76  if (rhoName_ == "rhoInf")
77  {
78  rhoInf_ = readScalar(dict.lookup("rhoInf"));
79  }
80 
81  if (dict.readIfPresent("g", g_))
82  {
83  lookupGravity_ = -2;
84  }
85 
86  if (!dict.found("value"))
87  {
88  updateCoeffs();
89  }
90 
91  if (dict.found("initialPoints"))
92  {
93  initialPoints_ = vectorField("initialPoints", dict , p.size());
94  }
95  else
96  {
97  initialPoints_ = p.localPoints();
98  }
99 }
100 
101 
104 (
106  const pointPatch& p,
108  const pointPatchFieldMapper& mapper
109 )
110 :
111  fixedValuePointPatchField<vector>(ptf, p, iF, mapper),
112  motion_(ptf.motion_),
113  initialPoints_(ptf.initialPoints_, mapper),
114  rhoInf_(ptf.rhoInf_),
115  rhoName_(ptf.rhoName_),
116  lookupGravity_(ptf.lookupGravity_),
117  g_(ptf.g_),
118  curTimeIndex_(-1)
119 {}
120 
121 
124 (
127 )
128 :
130  motion_(ptf.motion_),
131  initialPoints_(ptf.initialPoints_),
132  rhoInf_(ptf.rhoInf_),
133  rhoName_(ptf.rhoName_),
134  lookupGravity_(ptf.lookupGravity_),
135  g_(ptf.g_),
136  curTimeIndex_(-1)
137 {}
138 
139 
140 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
141 
143 (
144  const pointPatchFieldMapper& m
145 )
146 {
148 
149  initialPoints_.autoMap(m);
150 }
151 
152 
154 (
155  const pointPatchField<vector>& ptf,
156  const labelList& addr
157 )
158 {
160  refCast<const sixDoFRigidBodyDisplacementPointPatchVectorField>(ptf);
161 
163 
164  initialPoints_.rmap(sDoFptf.initialPoints_, addr);
165 }
166 
167 
169 {
170  if (this->updated())
171  {
172  return;
173  }
174 
175  if (lookupGravity_ < 0)
176  {
177  if (db().foundObject<uniformDimensionedVectorField>("g"))
178  {
179  if (lookupGravity_ == -2)
180  {
182  (
183  "void sixDoFRigidBodyDisplacementPointPatchVectorField"
184  "::updateCoeffs()"
185  )
186  << "Specifying the value of g in this boundary condition "
187  << "when g is available from the database is considered "
188  << "a fatal error to avoid the possibility of inconsistency"
189  << exit(FatalError);
190  }
191  else
192  {
193  lookupGravity_ = 1;
194  }
195  }
196  else
197  {
198  lookupGravity_ = 0;
199  }
200  }
201 
202  const polyMesh& mesh = this->dimensionedInternalField().mesh()();
203  const Time& t = mesh.time();
204  const pointPatch& ptPatch = this->patch();
205 
206  // Store the motion state at the beginning of the time-step
207  bool firstIter = false;
208  if (curTimeIndex_ != t.timeIndex())
209  {
210  motion_.newTime();
211  curTimeIndex_ = t.timeIndex();
212  firstIter = true;
213  }
214 
215  dictionary forcesDict;
216 
217  forcesDict.add("type", forces::typeName);
218  forcesDict.add("patches", wordList(1, ptPatch.name()));
219  forcesDict.add("rhoInf", rhoInf_);
220  forcesDict.add("rhoName", rhoName_);
221  forcesDict.add("CofR", motion_.centreOfRotation());
222 
223  forces f("forces", db(), forcesDict);
224 
225  f.calcForcesMoment();
226 
227  // Get the forces on the patch faces at the current positions
228 
229  if (lookupGravity_ == 1)
230  {
233 
234  g_ = g.value();
235  }
236 
237  // scalar ramp = min(max((t.value() - 5)/10, 0), 1);
238  scalar ramp = 1.0;
239 
240  motion_.update
241  (
242  firstIter,
243  ramp*(f.forceEff() + motion_.mass()*g_),
244  ramp*(f.momentEff() + motion_.mass()*(motion_.momentArm() ^ g_)),
245  t.deltaTValue(),
246  t.deltaT0Value()
247  );
248 
250  (
251  motion_.transform(initialPoints_) - initialPoints_
252  );
253 
255 }
256 
257 
259 {
261 
262  os.writeKeyword("rhoName") << rhoName_ << token::END_STATEMENT << nl;
263 
264  if (rhoName_ == "rhoInf")
265  {
266  os.writeKeyword("rhoInf") << rhoInf_ << token::END_STATEMENT << nl;
267  }
268 
269  if (lookupGravity_ == 0 || lookupGravity_ == -2)
270  {
271  os.writeKeyword("g") << g_ << token::END_STATEMENT << nl;
272  }
273 
274  motion_.write(os);
275 
276  initialPoints_.writeEntry("initialPoints", os);
277 
278  writeEntry("value", os);
279 }
280 
281 
282 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
283 
285 (
288 );
289 
290 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
291 
292 } // End namespace Foam
293 
294 // ************************************************************************* //
const objectRegistry & db() const
Return local objectRegistry.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:306
Basic pointPatch represents a set of points from the mesh.
Definition: pointPatch.H:56
virtual void calcForcesMoment()
Calculate the forces and moments.
Definition: forces.C:809
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
labelList f(nPoints)
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
A class for handling words, derived from string.
Definition: word.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dynamicFvMesh & mesh
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar mass() const
Return the mass.
Namespace for OpenFOAM.
virtual void rmap(const pointPatchField< vector > &, const labelList &)
Reverse map the given pointPatchField onto this pointPatchField.
point transform(const point &initialPoints) const
Transform the given initial state point by the current motion.
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
virtual void autoMap(const pointPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
makePointPatchTypeField(pointPatchVectorField, solidBodyMotionDisplacementPointPatchVectorField)
dictionary dict
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
static const char nl
Definition: Ostream.H:260
void newTime()
Store the motion state at the beginning of the time-step.
scalar deltaTValue() const
Return time step value.
Definition: TimeState.H:100
label timeIndex() const
Return current time index.
Definition: TimeState.C:73
volScalarField & p
Definition: createFields.H:51
const DimensionedField< vector, pointMesh > & dimensionedInternalField() const
Return dimensioned internal field reference.
Field< vector > vectorField
Specialisation of Field<T> for vector.
void update(bool firstIter, const vector &fGlobal, const vector &tauGlobal, scalar deltaT, scalar deltaT0)
Symplectic integration of velocities, orientation and position.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
void writeEntry(const word &keyword, Ostream &os) const
Write the field as a dictionary entry.
Definition: Field.C:634
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Macros for easy insertion into run-time selection tables.
virtual void rmap(const pointPatchField< Type > &, const labelList &)
Reverse map the given PointPatchField onto.
bool add(entry *, bool mergeEntry=false)
Add a new entry.
Definition: dictionary.C:739
Pre-declare SubField and related Field type.
Definition: Field.H:57
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:59
Foam::pointPatchFieldMapper.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
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
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const dimensionedVector & g
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:452
List< word > wordList
A List of words.
Definition: fileName.H:54
error FatalError
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
static const Vector zero
Definition: Vector.H:80
virtual const word & name() const =0
Return name.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
void write(Ostream &) const
Write.
sixDoFRigidBodyDisplacementPointPatchVectorField(const pointPatch &, const DimensionedField< vector, pointMesh > &)
Construct from patch and internal field.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
const Time & time() const
Return time.
const pointPatch & patch() const
Return patch.
bool updated() const
Return true if the boundary condition has already been updated.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const point & centreOfRotation() const
Return the current centre of rotation.
virtual vector forceEff() const
Return the total force.
Definition: forces.C:954
label size() const
Return size.
const Type & value() const
Return const reference to value.
virtual void write(Ostream &) const
Write.