src/dynamicMesh/meshMotion/solidBodyMotion/harmonicOscillation/harmonicOscillation.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   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 \*---------------------------------------------------------------------------*/
  26
  27 #include "harmonicOscillation.H"
  28 #include "addToRunTimeSelectionTable.H"
  29 #include "mathematicalConstants.H"
  30
  31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  32
  33 namespace Foam
  34 {
  35 namespace solidBodyMotionFunctions
  36 {
  37
  38 defineTypeNameAndDebug(harmonicOscillation, 0);
  39 addToRunTimeSelectionTable
  40 (
  41     solidBodyMotionFunction,
  42     harmonicOscillation,
  43     dictionary
  44 );
  45
  46 }
  47 }
  48
  49 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  50
  51 Foam::septernion
  52 Foam::solidBodyMotionFunctions::harmonicOscillation::calcTransformation
  53 (
  54     const scalar t
  55 ) const
  56 {
  57     const vector translation =
  58     cmptMultiply
  59     (
  60         transAmplitudes_,
  61         vector
  62         (
  63             sin(transAngularFreq_.x()*t + transPhaseAngles_.x()),
  64             sin(transAngularFreq_.y()*t + transPhaseAngles_.y()),
  65             sin(transAngularFreq_.z()*t + transPhaseAngles_.z())
  66         )
  67     );
  68
  69
  70     const vector rotation =
  71     cmptMultiply
  72     (
  73         rotAmplitudes_,
  74         vector
  75         (
  76             sin(rotAngularFreq_.x()*t + rotPhaseAngles_.x()),
  77             sin(rotAngularFreq_.y()*t + rotPhaseAngles_.y()),
  78             sin(rotAngularFreq_.z()*t + rotPhaseAngles_.z())
  79         )
  80     );
  81
  82     const   quaternion R(rotation.x(), rotation.y(), rotation.z());
  83     const   septernion TR
  84         (
  85             septernion(origin_ + translation)*R*septernion(-origin_)
  86         );
  87
  88     return TR;
  89 }
  90
  91 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  92
  93 Foam::solidBodyMotionFunctions::harmonicOscillation::harmonicOscillation
  94 (
  95     const dictionary& SBMFCoeffs,
  96     const Time& runTime
  97 )
  98 :
  99     solidBodyMotionFunction(SBMFCoeffs, runTime)
 100 {
 101     read(SBMFCoeffs);
 102 }
 103
 104 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 105
 106 Foam::septernion
 107 Foam::solidBodyMotionFunctions::harmonicOscillation::transformation() const
 108 {
 109    const scalar t = time_.value();
 110
 111    const septernion TR = calcTransformation(t);
 112
 113    Info<< "solidBodyMotionFunctions::harmonicOscillation::transformation(): "
 114    << "Time = " << t << " transformation: " << TR << endl;
 115
 116    return TR;
 117 }
 118
 119
 120 Foam::septernion
 121 Foam::solidBodyMotionFunctions::harmonicOscillation::velocity() const
 122 {
 123     const scalar t = time_.value();
 124     const scalar dt = time_.deltaT().value();
 125
 126     const septernion velocity
 127     (
 128         (calcTransformation(t).t() - calcTransformation(t - dt).t())/dt,
 129         (calcTransformation(t).r()/calcTransformation(t - dt).r())/dt
 130     );
 131
 132     return velocity;
 133 }
 134
 135
 136 bool Foam::solidBodyMotionFunctions::harmonicOscillation::read
 137 (
 138     const dictionary& SBMFCoeffs
 139 )
 140 {
 141     solidBodyMotionFunction::read(SBMFCoeffs);
 142
 143     SBMFCoeffs_.lookup("origin") >> origin_;
 144     SBMFCoeffs_.lookup("translationalAmplitudes") >> transAmplitudes_;
 145     SBMFCoeffs_.lookup("translationalAngularFrequencies") >> transAngularFreq_;
 146     SBMFCoeffs_.lookup("translationalPhaseAngles") >> transPhaseAngles_;
 147     SBMFCoeffs_.lookup("rotationalAmplitudes") >> rotAmplitudes_;
 148     SBMFCoeffs_.lookup("rotationalAngularFrequencies") >> rotAngularFreq_;
 149     SBMFCoeffs_.lookup("rotationalPhaseAngles") >> rotPhaseAngles_;
 150     SBMFCoeffs_.lookup("inDegrees") >> inDegrees_;
 151
 152     // Sanity check for negative frequencies
 153     if
 154     (
 155         transAngularFreq_.x() < -SMALL
 156      || transAngularFreq_.y() < -SMALL
 157      || transAngularFreq_.z() < -SMALL
 158      || rotAngularFreq_.x() < -SMALL
 159      || rotAngularFreq_.y() < -SMALL
 160      || rotAngularFreq_.z() < -SMALL
 161     )
 162     {
 163         FatalErrorIn
 164         (
 165            "harmonicOscillation::read"
 166              "(\n"
 167              "    const fvMesh& mesh\n"
 168              "    const dictionary& dict\n"
 169              ")"
 170         )   << "Negative angular frequency detected."
 171             << abort(FatalError);
 172     }
 173
 174     // Convert to radians if necessary, printing data
 175     if (inDegrees_)
 176     {
 177         const scalar piBy180 = mathematicalConstant::pi/180.0;
 178
 179         transPhaseAngles_ *= piBy180;
 180         rotAmplitudes_ *= piBy180;
 181         rotPhaseAngles_ *= piBy180;
 182
 183         Info<< "Data in degrees converted:" << nl
 184             << "translationalPhaseAngles = " << transPhaseAngles_ << nl
 185             << "rotationalAmplitudes = " << rotAmplitudes_ << nl
 186             << "rotationalPhaseAngles = " << rotPhaseAngles_ << endl;
 187     }
 188
 189     // Count prescribed rotations: only two harmonic rotations can be prescribed
 190     label nRot = 0;
 191
 192     // Rotation
 193     mag(rotAmplitudes_.x()) > SMALL ? ++nRot : /* null expression */ 0 ;
 194     mag(rotAmplitudes_.y()) > SMALL ? ++nRot : /* null expression */ 0 ;
 195     mag(rotAmplitudes_.z()) > SMALL ? ++nRot : /* null expression */ 0 ;
 196
 197     // Check if more than two rotations prescribed
 198     if (nRot > 2)
 199     {
 200         FatalErrorIn
 201         (
 202             "harmonicOscillation::read"
 203              "(\n"
 204              "    const fvMesh& mesh\n"
 205              "    const dictionary& dict\n"
 206              ")"
 207         )   << "Only two harmonic rotations can be prescribed. Prescribing "
 208             << "more than two yields unphysical results. "
 209             << abort(FatalError);
 210     }
 211
 212     return true;
 213 }
 214
 215
 216 // ************************************************************************* //