sixDoFRigidBodyMotionI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2015 OpenFOAM Foundation
9 Copyright (C) 2016-2019 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30
31inline Foam::tensor Foam::sixDoFRigidBodyMotion::rotationTensorX
32(
33 scalar phi
34) const
35{
36 return tensor
37 (
38 1, 0, 0,
41 );
42}
43
44
45inline Foam::tensor Foam::sixDoFRigidBodyMotion::rotationTensorY
46(
47 scalar phi
48) const
49{
50 return tensor
51 (
53 0, 1, 0,
55 );
56}
57
58
59inline Foam::tensor Foam::sixDoFRigidBodyMotion::rotationTensorZ
60(
61 scalar phi
62) const
63{
64 return tensor
65 (
68 0, 0, 1
69 );
70}
71
72
74Foam::sixDoFRigidBodyMotion::rotate
75(
76 const tensor& Q0,
77 const vector& pi0,
78 const scalar deltaT
79) const
80{
81 Tuple2<tensor, vector> Qpi(Q0, pi0);
82 tensor& Q = Qpi.first();
83 vector& pi = Qpi.second();
84
85 tensor R = rotationTensorX(0.5*deltaT*pi.x()/momentOfInertia_.xx());
86 pi = pi & R;
87 Q = Q & R;
88
89 R = rotationTensorY(0.5*deltaT*pi.y()/momentOfInertia_.yy());
90 pi = pi & R;
91 Q = Q & R;
92
93 R = rotationTensorZ(deltaT*pi.z()/momentOfInertia_.zz());
94 pi = pi & R;
95 Q = Q & R;
96
97 R = rotationTensorY(0.5*deltaT*pi.y()/momentOfInertia_.yy());
98 pi = pi & R;
99 Q = Q & R;
100
101 R = rotationTensorX(0.5*deltaT*pi.x()/momentOfInertia_.xx());
102 pi = pi & R;
103 Q = Q & R;
104
105 return Qpi;
106}
107
108
110Foam::sixDoFRigidBodyMotion::restraints() const
111{
112 return restraints_;
113}
114
115
117Foam::sixDoFRigidBodyMotion::constraints() const
118{
119 return constraints_;
120}
121
122
123inline const Foam::point&
124Foam::sixDoFRigidBodyMotion::initialCentreOfRotation() const
125{
126 return initialCentreOfRotation_;
127}
128
129
130inline const Foam::tensor&
131Foam::sixDoFRigidBodyMotion::initialQ() const
132{
133 return initialQ_;
134}
135
136
138{
139 return motionState_.Q();
140}
141
142
144{
145 return motionState_.v();
146}
147
148
149inline const Foam::vector& Foam::sixDoFRigidBodyMotion::a() const
150{
151 return motionState_.a();
152}
153
154
155inline const Foam::vector& Foam::sixDoFRigidBodyMotion::pi() const
156{
157 return motionState_.pi();
158}
159
160
161inline const Foam::vector& Foam::sixDoFRigidBodyMotion::tau() const
162{
163 return motionState_.tau();
164}
165
166
167inline Foam::point& Foam::sixDoFRigidBodyMotion::initialCentreOfRotation()
168{
169 return initialCentreOfRotation_;
170}
171
172
173inline Foam::tensor& Foam::sixDoFRigidBodyMotion::initialQ()
174{
175 return initialQ_;
176}
177
178
180{
181 return motionState_.Q();
182}
183
184
186{
187 return motionState_.v();
188}
189
190
191inline Foam::vector& Foam::sixDoFRigidBodyMotion::a()
192{
193 return motionState_.a();
194}
195
196
197inline Foam::vector& Foam::sixDoFRigidBodyMotion::pi()
198{
199 return motionState_.pi();
200}
201
202
203inline Foam::vector& Foam::sixDoFRigidBodyMotion::tau()
204{
205 return motionState_.tau();
206}
207
208
209// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
210
211inline Foam::scalar Foam::sixDoFRigidBodyMotion::mass() const
212{
213 return mass_;
214}
215
216
217inline const Foam::diagTensor&
219{
220 return momentOfInertia_;
221}
222
223
226{
227 return motionState_;
228}
229
230
232{
233 return motionState_.centreOfRotation();
234}
235
236
237inline const Foam::point&
239{
240 return initialCentreOfMass_;
241}
242
243
245{
246 return transform(initialCentreOfMass_);
247}
248
249
251{
252 return centreOfMass() - motionState_.centreOfRotation();
253}
254
255
256inline const Foam::tensor&
258{
259 return Q();
260}
261
262
264{
265 return Q() & (inv(momentOfInertia_) & pi());
266}
267
269{
270 return time_;
271}
272
274{
275 return report_;
276}
277
278
280{
281 motionState0_ = motionState_;
282}
283
284
286{
287 return motionState_.centreOfRotation();
288}
289
290
292(
293 const point& pt
294) const
295{
296 return (omega() ^ (pt - centreOfRotation())) + v();
297}
298
299
301(
302 const point& initialPoint
303) const
304{
305 return
306 (
307 centreOfRotation()
308 + (Q() & initialQ_.T() & (initialPoint - initialCentreOfRotation_))
309 );
310}
311
312// ************************************************************************* //
#define R(A, B, C, D, E, F, K, M)
surfaceScalarField & phi
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
A 2-tuple for storing two objects of dissimilar types. The container is similar in purpose to std::pa...
Definition: Tuple2.H:58
Computes the second invariant of the velocity gradient tensor .
Definition: Q.H:153
Q(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
Definition: Q.C:69
Default transformation behaviour.
Holds the motion state of sixDoF object. Wrapped up together to allow rapid scatter to other processo...
point centreOfMass() const
Return the current centre of mass.
const diagTensor & momentOfInertia() const
Return the inertia tensor.
void newTime()
Store the motion state at the beginning of the time-step.
const vector & v() const
Return the current velocity.
const point & initialCentreOfMass() const
Return the initial centre of mass.
bool report() const
Return the report Switch.
const Time & time() const
Return time.
const tensor & orientation() const
Return the orientation tensor, Q.
vector momentArm() const
Return the current momentArm.
vector omega() const
Return the angular velocity in the global frame.
scalar mass() const
Return the mass.
point velocity(const point &pt) const
Return the velocity of a position.
const point & centreOfRotation() const
Return the current centre of rotation.
const sixDoFRigidBodyMotionState & state() const
Return the motion state.
Tensor of scalars, i.e. Tensor<scalar>.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
constexpr scalar pi(M_PI)
dimensionSet transform(const dimensionSet &ds)
Return the argument; transformations do not change the dimensions.
Definition: dimensionSet.C:536
dimensionedScalar sin(const dimensionedScalar &ds)
Tensor< scalar > tensor
Definition: symmTensor.H:61
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
dimensionedScalar cos(const dimensionedScalar &ds)