src/ODE/translationODE/translationODE.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original authors
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 Class
  26     translationODE
  27
  28 Description
  29     Ordinary differential equation for three degrees of freedom
  30     solid body motion
  31
  32 Author
  33     Hrvoje Jasak
  34     Dubravko Matijasevic
  35
  36 \*---------------------------------------------------------------------------*/
  37
  38 #include "translationODE.H"
  39 #include "Time.H"
  40
  41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  42
  43 // namespace Foam
  44 // {
  45 //     defineTypeNameAndDebug(translationODE, 0);
  46 // }
  47
  48 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  49
  50 void Foam::translationODE::setCoeffs()
  51 {
  52     // Set ODE coefficients from position and rotation
  53
  54     // Linear displacement in relative coordinate system
  55     {
  56         const vector& xVal = Xrel_.value();
  57         coeffs_[0] = xVal.x();
  58         coeffs_[1] = xVal.y();
  59         coeffs_[2] = xVal.z();
  60     }
  61
  62     // Linear velocity in relative coordinate system
  63     {
  64         const vector& uVal = U_.value();
  65         coeffs_[3] = uVal.x();
  66         coeffs_[4] = uVal.y();
  67         coeffs_[5] = uVal.z();
  68     }
  69 }
  70
  71
  72 Foam::dimensionedVector Foam::translationODE::A
  73 (
  74     const dimensionedVector& xR,
  75     const dimensionedVector& uR
  76 ) const
  77 {
  78     return
  79     (
  80        - (linSpringCoeffs_ & xR)    // spring
  81        - (linDampingCoeffs_ & uR)   // damping
  82        + force()
  83     )/mass_;                        // gravity
  84 }
  85
  86
  87 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  88
  89 // Construct from components
  90 Foam::translationODE::translationODE
  91 (
  92     const IOobject& io
  93 )
  94 :
  95     IOdictionary(io),
  96     mass_(lookup("mass")),
  97     xEquilibrium_(lookup("equilibriumPosition")),
  98     linSpringCoeffs_(lookup("linearSpring")),
  99     linDampingCoeffs_(lookup("linearDamping")),
 100     Xrel_(lookup("Xrel")),
 101     U_(lookup("U")),
 102     Uold_(lookup("Uold")),
 103     force_(lookup("force")),
 104     coeffs_(6, 0.0)
 105 {
 106     setCoeffs();
 107 }
 108
 109
 110 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 111
 112 Foam::translationODE::~translationODE()
 113 {}
 114
 115
 116 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 117
 118 void Foam::translationODE::derivatives
 119 (
 120     const scalar x,
 121     const scalarField& y,
 122     scalarField& dydx
 123 ) const
 124 {
 125     // Set the derivatives for displacement
 126     dydx[0] = y[3];
 127     dydx[1] = y[4];
 128     dydx[2] = y[5];
 129
 130     dimensionedVector curX("curX", dimLength, vector(y[0], y[1], y[2]));
 131     dimensionedVector curU("curU", dimLength/dimTime, vector(y[3], y[4], y[5]));
 132
 133     const vector& accel = A(curX, curU).value();
 134
 135     dydx[3] = accel.x();
 136     dydx[4] = accel.y();
 137     dydx[5] = accel.z();
 138 }
 139
 140
 141 void Foam::translationODE::jacobian
 142 (
 143     const scalar x,
 144     const scalarField& y,
 145     scalarField& dfdx,
 146     scalarSquareMatrix& dfdy
 147 ) const
 148 {
 149     notImplemented("translationODE::jacobian(...) const");
 150 }
 151
 152
 153 void Foam::translationODE::update(const scalar delta)
 154 {
 155     // Update position
 156     vector& Xval = Xrel_.value();
 157
 158     Xval.x() = coeffs_[0];
 159     Xval.y() = coeffs_[1];
 160     Xval.z() = coeffs_[2];
 161
 162     // Update velocity
 163     Uold_ = U_;
 164
 165     vector& Uval = U_.value();
 166
 167     Uval.x() = coeffs_[3];
 168     Uval.y() = coeffs_[4];
 169     Uval.z() = coeffs_[5];
 170 }
 171
 172
 173 bool Foam::translationODE::writeData(Ostream& os) const
 174 {
 175     os << *this;
 176     return os.good();
 177 }
 178
 179
 180 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 181
 182 Foam::Ostream& Foam::operator<<(Ostream& os, const translationODE& sds)
 183 {
 184     os.writeKeyword("mass") << sds.mass_ << token::END_STATEMENT << endl;
 185     os.writeKeyword("equilibriumPosition") << sds.xEquilibrium_
 186         << token::END_STATEMENT << endl;
 187     os.writeKeyword("linearSpring") << sds.linSpringCoeffs_
 188         << token::END_STATEMENT << endl;
 189     os.writeKeyword("linearDamping") << sds.linDampingCoeffs_
 190         << token::END_STATEMENT << endl;
 191
 192     os.writeKeyword("Xrel") << sds.Xrel() << token::END_STATEMENT << endl;
 193     os.writeKeyword("U") << sds.U() << token::END_STATEMENT << endl;
 194     os.writeKeyword("Uold") << tab << sds.Uold() << token::END_STATEMENT << nl;
 195
 196     os.writeKeyword("force") << sds.force() << token::END_STATEMENT << endl;
 197
 198     return os;
 199 }
 200
 201
 202 // ************************************************************************* //