src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2009-2010 OpenCFD Ltd.
   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 Class
  25     Foam::sixDoFRigidBodyMotion
  26
  27 Description
  28     Six degree of freedom motion for a rigid body.  Angular momentum
  29     stored in body fixed reference frame.  Reference orientation of
  30     the body (where Q = I) must align with the cartesian axes such
  31     that the Inertia tensor is in principle component form.
  32
  33     Symplectic motion as per:
  34
  35     title = {Symplectic splitting methods for rigid body molecular dynamics},
  36     publisher = {AIP},
  37     year = {1997},
  38     journal = {The Journal of Chemical Physics},
  39     volume = {107},
  40     number = {15},
  41     pages = {5840-5851},
  42     url = {http://link.aip.org/link/?JCP/107/5840/1},
  43     doi = {10.1063/1.474310}
  44
  45     Can add restraints (i.e. a spring) and constraints (i.e. motion
  46     may only be on a plane).
  47
  48 SourceFiles
  49     sixDoFRigidBodyMotionI.H
  50     sixDoFRigidBodyMotion.C
  51     sixDoFRigidBodyMotionIO.C
  52
  53 \*---------------------------------------------------------------------------*/
  54
  55 #ifndef sixDoFRigidBodyMotion_H
  56 #define sixDoFRigidBodyMotion_H
  57
  58 #include "sixDoFRigidBodyMotionState.H"
  59 #include "pointField.H"
  60 #include "sixDoFRigidBodyMotionRestraint.H"
  61 #include "sixDoFRigidBodyMotionConstraint.H"
  62
  63
  64 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  65
  66 namespace Foam
  67 {
  68
  69 // Forward declaration of classes
  70 class Istream;
  71 class Ostream;
  72
  73 // Forward declaration of friend functions and operators
  74 class sixDoFRigidBodyMotion;
  75 Istream& operator>>(Istream&, sixDoFRigidBodyMotion&);
  76 Ostream& operator<<(Ostream&, const sixDoFRigidBodyMotion&);
  77
  78
  79 /*---------------------------------------------------------------------------*\
  80                       Class sixDoFRigidBodyMotion Declaration
  81 \*---------------------------------------------------------------------------*/
  82
  83 class sixDoFRigidBodyMotion
  84 {
  85     // Private data
  86
  87         //- Motion state data object
  88         sixDoFRigidBodyMotionState motionState_;
  89
  90         //- Restraints on the motion
  91         PtrList<sixDoFRigidBodyMotionRestraint> restraints_;
  92
  93         //- Names of the restraints
  94         wordList restraintNames_;
  95
  96         //- Constaints on the motion
  97         PtrList<sixDoFRigidBodyMotionConstraint> constraints_;
  98
  99         //- Names of the constraints
 100         wordList constraintNames_;
 101
 102         //- Maximum number of iterations allowed to attempt to obey
 103         //  constraints
 104         label maxConstraintIterations_;
 105
 106         //- Centre of mass of initial state
 107         point initialCentreOfMass_;
 108
 109         //- Orientation of initial state
 110         tensor initialQ_;
 111
 112         //- Moment of inertia of the body in reference configuration
 113         //  (Q = I)
 114         diagTensor momentOfInertia_;
 115
 116         //- Mass of the body
 117         scalar mass_;
 118
 119         //- Switch to turn reporting of motion data on and off
 120         Switch report_;
 121
 122
 123     // Private Member Functions
 124
 125         //- Calculate the rotation tensor around the body reference
 126         //  frame x-axis by the given angle
 127         inline tensor rotationTensorX(scalar deltaT) const;
 128
 129         //- Calculate the rotation tensor around the body reference
 130         //  frame y-axis by the given angle
 131         inline tensor rotationTensorY(scalar deltaT) const;
 132
 133         //- Calculate the rotation tensor around the body reference
 134         //  frame z-axis by the given angle
 135         inline tensor rotationTensorZ(scalar deltaT) const;
 136
 137         //- Apply rotation tensors to Q for the given torque (pi) and deltaT
 138         inline void rotate(tensor& Q, vector& pi, scalar deltaT) const;
 139
 140         //- Apply the restraints to the object
 141         void applyRestraints();
 142
 143         //- Apply the constraints to the object
 144         void applyConstraints(scalar deltaT);
 145
 146         // Access functions retained as private because of the risk of
 147         // confusion over what is a body local frame vector and what is global
 148
 149         // Access
 150
 151             //- Return access to the motion state
 152             inline const sixDoFRigidBodyMotionState& motionState() const;
 153
 154             //- Return access to the restraints
 155             inline const PtrList<sixDoFRigidBodyMotionRestraint>&
 156                 restraints() const;
 157
 158             //- Return access to the restraintNames
 159             inline const wordList& restraintNames() const;
 160
 161             //- Return access to the constraints
 162             inline const PtrList<sixDoFRigidBodyMotionConstraint>&
 163                 constraints() const;
 164
 165             //- Return access to the constraintNames
 166             inline const wordList& constraintNames() const;
 167
 168             //- Return access to the maximum allowed number of
 169             //  constraint iterations
 170             inline label maxConstraintIterations() const;
 171
 172             //- Return access to the initial centre of mass
 173             inline const point& initialCentreOfMass() const;
 174
 175             //- Return access to the initial orientation
 176             inline const tensor& initialQ() const;
 177
 178             //- Return access to the orientation
 179             inline const tensor& Q() const;
 180
 181             //- Return access to velocity
 182             inline const vector& v() const;
 183
 184             //- Return access to acceleration
 185             inline const vector& a() const;
 186
 187             //- Return access to angular momentum
 188             inline const vector& pi() const;
 189
 190             //- Return access to torque
 191             inline const vector& tau() const;
 192
 193
 194         // Edit
 195
 196             //- Return access to the centre of mass
 197             inline point& initialCentreOfMass();
 198
 199             //- Return access to the centre of mass
 200             inline tensor& initialQ();
 201
 202             //- Return non-const access to the orientation
 203             inline tensor& Q();
 204
 205             //- Return non-const access to vector
 206             inline vector& v();
 207
 208             //- Return non-const access to acceleration
 209             inline vector& a();
 210
 211             //- Return non-const access to angular momentum
 212             inline vector& pi();
 213
 214             //- Return non-const access to torque
 215             inline vector& tau();
 216
 217
 218 public:
 219
 220     // Constructors
 221
 222         //- Construct null
 223         sixDoFRigidBodyMotion();
 224
 225         //- Construct from components
 226         sixDoFRigidBodyMotion
 227         (
 228             const point& centreOfMass,
 229             const tensor& Q,
 230             const vector& v,
 231             const vector& a,
 232             const vector& pi,
 233             const vector& tau,
 234             scalar mass,
 235             const point& initialCentreOfMass,
 236             const tensor& initialQ,
 237             const diagTensor& momentOfInertia,
 238             bool report = false
 239         );
 240
 241         //- Construct from dictionary
 242         sixDoFRigidBodyMotion(const dictionary& dict);
 243
 244         //- Construct as copy
 245         sixDoFRigidBodyMotion(const sixDoFRigidBodyMotion&);
 246
 247
 248     //- Destructor
 249     ~sixDoFRigidBodyMotion();
 250
 251
 252     // Member Functions
 253
 254         //- Add restraints to the motion, public to allow external
 255         //  addition of restraints after construction
 256         void addRestraints(const dictionary& dict);
 257
 258         //- Add restraints to the motion, public to allow external
 259         //  addition of restraints after construction
 260         void addConstraints(const dictionary& dict);
 261
 262         //- First leapfrog velocity adjust and motion part, required
 263         // before force calculation
 264         void updatePosition
 265         (
 266             scalar deltaT
 267         );
 268
 269         //- Second leapfrog velocity adjust part, required after motion and
 270         // force calculation
 271         void updateForce
 272         (
 273             const vector& fGlobal,
 274             const vector& tauGlobal,
 275             scalar deltaT
 276         );
 277
 278         //- Global forces supplied at locations, calculating net force
 279         //  and moment
 280         void updateForce
 281         (
 282             const pointField& positions,
 283             const vectorField& forces,
 284             scalar deltaT
 285         );
 286
 287         //- Transform the given initial state pointField by the current
 288         //  motion state
 289         inline tmp<pointField> currentPosition
 290         (
 291             const pointField& pInitial
 292         ) const;
 293
 294         //- Transform the given initial state point by the current motion
 295         //  state
 296         inline point currentPosition(const point& pInitial) const;
 297
 298         //- Transform the given initial state direction by the current
 299         //  motion state
 300         inline vector currentOrientation(const vector& vInitial) const;
 301
 302         //- Access the orientation tensor, Q.
 303         //  globalVector = Q & bodyLocalVector
 304         //  bodyLocalVector = Q.T() & globalVector
 305         inline const tensor& orientation() const;
 306
 307         //- Predict the position of the supplied initial state point
 308         //  after deltaT given the current motion state and the
 309         //  additional supplied force and moment
 310         point predictedPosition
 311         (
 312             const point& pInitial,
 313             const vector& deltaForce,
 314             const vector& deltaMoment,
 315             scalar deltaT
 316         ) const;
 317
 318         //- Predict the orientation of the supplied initial state
 319         //  vector after deltaT given the current motion state and the
 320         //  additional supplied moment
 321         vector predictedOrientation
 322         (
 323             const vector& vInitial,
 324             const vector& deltaMoment,
 325             scalar deltaT
 326         ) const;
 327
 328         //- Return the angular velocity in the global frame
 329         inline vector omega() const;
 330
 331         //- Return the velocity of a position given by the current
 332         //  motion state
 333         inline point currentVelocity(const point& pt) const;
 334
 335         //- Report the status of the motion
 336         void status() const;
 337
 338
 339         // Access
 340
 341             //- Return const access to the centre of mass
 342             inline const point& centreOfMass() const;
 343
 344             //- Return access to the inertia tensor
 345             inline const diagTensor& momentOfInertia() const;
 346
 347             //- Return const access to the mass
 348             inline scalar mass() const;
 349
 350             //- Return the report Switch
 351             inline bool report() const;
 352
 353
 354         // Edit
 355
 356             //- Return non-const access to the centre of mass
 357             inline point& centreOfMass();
 358
 359             //- Return non-const access to the inertia tensor
 360             inline diagTensor& momentOfInertia();
 361
 362             //- Return non-const access to the mass
 363             inline scalar& mass();
 364
 365
 366         //- Write
 367         void write(Ostream&) const;
 368
 369
 370     // IOstream Operators
 371
 372         friend Istream& operator>>(Istream&, sixDoFRigidBodyMotion&);
 373         friend Ostream& operator<<(Ostream&, const sixDoFRigidBodyMotion&);
 374 };
 375
 376
 377 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 378
 379 } // End namespace Foam
 380
 381 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 382
 383 #include "sixDoFRigidBodyMotionI.H"
 384
 385 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 386
 387 #endif
 388
 389 // ************************************************************************* //