Abstract base class sixDOFODE

This commit is contained in:
Vuko Vukcevic 2017-02-24 11:13:54 +01:00 committed by Hrvoje Jasak
parent 0642e99a97
commit 6436f1a496
7 changed files with 784 additions and 487 deletions

View file

@ -9,12 +9,12 @@ $(ODESolvers)/SIBS/SIBS.C
$(ODESolvers)/SIBS/SIMPR.C
$(ODESolvers)/SIBS/polyExtrapolate.C
translationODE = translationODE
$(translationODE)/translationODE.C
sixDOF = sixDOF
$(sixDOF)/finiteRotation/finiteRotation.C
$(sixDOF)/sixDOFqODE/sixDOFqODE.C
$(sixDOF)/sixDOFbodies/sixDOFbodies.C
$(sixDOF)/sixDOFODE/sixDOFODE.C
$(sixDOF)/sixDOFODE/newSixDOFODE.C
LIB = $(FOAM_LIBBIN)/libODE

View file

@ -0,0 +1,73 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright held by original author
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software; you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation; either version 2 of the License, or (at your
option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM; if not, write to the Free Software Foundation,
Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
\*---------------------------------------------------------------------------*/
#include "error.H"
#include "objectRegistry.H"
#include "sixDOFODE.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Foam::autoPtr<Foam::sixDOFODE> Foam::sixDOFODE::New(const IOobject& io)
{
word sixDOFODETypeName;
// Get object registry
const objectRegistry& database = io.db();
// Check whether the dictionary is in the database
if (database.foundObject<IOdictionary>(io.name()))
{
sixDOFODETypeName =
word
(
database.lookupObject<IOdictionary>(io.name()).lookup("type")
);
}
else
{
sixDOFODETypeName = word(IOdictionary(io).lookup("type"));
}
dictionaryConstructorTable::iterator cstrIter =
dictionaryConstructorTablePtr_->find(sixDOFODETypeName);
if (cstrIter == dictionaryConstructorTablePtr_->end())
{
FatalErrorIn
(
"sixDOFODE::New(const IOobject& io)"
) << "Unknown sixDOFODE " << sixDOFODETypeName
<< endl << endl
<< "Valid sixDOFODE types are:" << endl
<< dictionaryConstructorTablePtr_->toc()
<< exit(FatalError);
}
return autoPtr<sixDOFODE>(cstrIter()(io));
}
// ************************************************************************* //

View file

@ -0,0 +1,252 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation, either version 3 of the License, or (at your
option) any later version.
foam-extend is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
General Public License for more details.
You should have received a copy of the GNU General Public License
along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
sixDOFODE
Description
Abstract base class for six-degrees-of-freedom (6DOF) ordinary differential
equations
Author
Dubravko Matijasevic, FSB Zagreb. All rights reserved.
Hrvoje Jasak, FSB Zagreb. All rights reserved.
Vuko Vukcevic, FSB Zagreb. All rights reserved.
\*---------------------------------------------------------------------------*/
#include "objectRegistry.H"
#include "sixDOFODE.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
defineTypeNameAndDebug(sixDOFODE, 0);
defineRunTimeSelectionTable(sixDOFODE, dictionary);
}
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
void Foam::sixDOFODE::aitkensRelaxation
(
const scalar min,
const scalar max
)
{
// Calculate translational relax factor
const scalar saveOldRelFacT = oldRelaxFactorT_;
oldRelaxFactorT_ = relaxFactorT_;
if(magSqr(A_[0] - An_[1] - A_[1] - An_[2]) > SMALL)
{
relaxFactorT_ = saveOldRelFacT + (saveOldRelFacT - 1)*
(
(A_[1] - An_[2])
& (A_[0] - An_[1] - A_[1] - An_[2])
)/
magSqr(A_[0] - An_[1] - A_[1] - An_[2]);
}
else
{
relaxFactorT_ = min;
}
// Bound the relaxation factor for stability
if (relaxFactorT_ > max)
{
relaxFactorT_ = max;
}
else if (relaxFactorT_ < min)
{
relaxFactorT_ = min;
}
// Calculate rotational relax factor
const scalar saveOldRelFacR = oldRelaxFactorR_;
oldRelaxFactorR_ = relaxFactorR_;
if
(
magSqr(OmegaDot_[0] - OmegaDotn_[1] - OmegaDot_[1] - OmegaDotn_[2])
> SMALL
)
{
relaxFactorR_ = saveOldRelFacR
+ (saveOldRelFacR - 1)*
(
(OmegaDot_[1] - OmegaDotn_[2])
& (OmegaDot_[0] - OmegaDotn_[1] - OmegaDot_[1] - OmegaDotn_[2])
)/
magSqr(OmegaDot_[0] - OmegaDotn_[1] - OmegaDot_[1] - OmegaDotn_[2]);
}
else
{
relaxFactorR_ = min;
}
// Bound the relaxation factor for stability
if(relaxFactorR_ > max)
{
relaxFactorR_ = max;
}
else if(relaxFactorR_ < min)
{
relaxFactorR_ = min;
}
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
Foam::sixDOFODE::sixDOFODE(const IOobject& io)
:
IOdictionary(io),
mass_(lookup("mass")),
momentOfInertia_(lookup("momentOfInertia")),
Xequilibrium_(lookup("equilibriumPosition")),
linSpringCoeffs_(lookup("linearSpring")),
linDampingCoeffs_(lookup("linearDamping")),
relaxFactorT_(1.0),
relaxFactorR_(1.0),
oldRelaxFactorT_(1.0),
oldRelaxFactorR_(1.0),
A_(3, vector::zero),
OmegaDot_(3, vector::zero),
An_(3, vector::zero),
OmegaDotn_(3, vector::zero),
force_(lookup("force")),
moment_(lookup("moment")),
forceRelative_(lookup("forceRelative")),
momentRelative_(lookup("momentRelative"))
{}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::sixDOFODE::~sixDOFODE()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::sixDOFODE::relaxAcceleration
(
const scalar minRelFactor,
const scalar maxRelFactor
)
{
if (mag(minRelFactor - maxRelFactor) < SMALL)
{
// Fixed relaxation
relaxFactorT_ = minRelFactor;
relaxFactorR_ = minRelFactor;
}
else
{
// Use Aitkens relaxation
// Update Aitkens relaxation factor
aitkensRelaxation(minRelFactor, maxRelFactor);
// Update non relaxed accelerations
An_[1] = An_[2];
An_[2] = (force_/mass_).value();
OmegaDotn_[1] = OmegaDotn_[2];
OmegaDotn_[2] = (inv(momentOfInertia_) & moment_).value();
}
const vector Aold = A_[2];
const vector OmegaDotold = OmegaDot_[2];
force_.value() =
Aold*mass_.value()
+ relaxFactorT_*(force_.value() - Aold*mass_.value());
moment_.value() =
(momentOfInertia_.value() & OmegaDotold)
+ relaxFactorR_*
(
moment_.value()
- (momentOfInertia_.value() & OmegaDotold)
);
// Update relaxed old accelerations
A_[0] = A_[1];
A_[1] = A_[2];
A_[2] = (force_/mass_).value();
OmegaDot_[0] = OmegaDot_[1];
OmegaDot_[1] = OmegaDot_[2];
OmegaDot_[2] = (inv(momentOfInertia_) & moment_).value();
}
// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
bool Foam::sixDOFODE::writeData(Ostream& os) const
{
os.writeKeyword("mass") << tab << mass_ << token::END_STATEMENT << nl;
os.writeKeyword("momentOfInertia") << tab << momentOfInertia_
<< token::END_STATEMENT << nl << nl;
os.writeKeyword("equilibriumPosition") << tab << Xequilibrium_
<< token::END_STATEMENT << nl;
os.writeKeyword("linearSpring") << tab << linSpringCoeffs_
<< token::END_STATEMENT << nl;
os.writeKeyword("linearDamping") << tab << linDampingCoeffs_
<< token::END_STATEMENT << nl << nl;
os.writeKeyword("force") << tab << force_
<< token::END_STATEMENT << nl;
os.writeKeyword("moment") << tab << moment_
<< token::END_STATEMENT << nl;
os.writeKeyword("forceRelative") << tab << forceRelative_
<< token::END_STATEMENT << nl;
os.writeKeyword("momentRelative") << tab << momentRelative_
<< token::END_STATEMENT << endl;
return os.good();
}
// * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * * //
Foam::Ostream& Foam::operator<<(Ostream& os, const sixDOFODE& sds)
{
sds.writeData(os);
os.check("Ostream& operator<<(Ostream&, const sixDOFODE");
return os;
}
// ************************************************************************* //

View file

@ -0,0 +1,395 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation, either version 3 of the License, or (at your
option) any later version.
foam-extend is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
General Public License for more details.
You should have received a copy of the GNU General Public License
along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
sixDOFODE
Description
Abstract base class for six-degrees-of-freedom (6DOF) ordinary differential
equations with a fat interface in order to provide representation of
rotations in quaternions or rotation tensors.
Author
Dubravko Matijasevic, FSB Zagreb. All rights reserved.
Hrvoje Jasak, FSB Zagreb. All rights reserved.
Vuko Vukcevic, FSB Zagreb. All rights reserved.
SourceFiles
sixDOFODEI.H
sixDOFODE.C
newSixDOFODE.C
\*---------------------------------------------------------------------------*/
#ifndef sixDOFODE_H
#define sixDOFODE_H
#include "ODE.H"
#include "IOdictionary.H"
#include "dimensionedTypes.H"
#include "autoPtr.H"
#include "runTimeSelectionTables.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
class finiteRotation;
class HamiltonRodriguezRot;
/*---------------------------------------------------------------------------*\
Class sixDOFODE Declaration
\*---------------------------------------------------------------------------*/
class sixDOFODE
:
public IOdictionary,
public ODE
{
// Private data
// Body data
//- Mass
dimensionedScalar mass_;
//- Rotational moment of inertia around centre of mass
// in body (relative) coordinates - aligned with main axes
dimensionedDiagTensor momentOfInertia_;
// Equilibrium position and spring/damping coefficients
//- Spring equilibrium position for translation
dimensionedVector Xequilibrium_;
//- Linear spring coeffs
dimensionedDiagTensor linSpringCoeffs_;
//- Linear damping coeffs
dimensionedDiagTensor linDampingCoeffs_;
// Aitkens relaxation data members
//- Translational relaxation factor
scalar relaxFactorT_;
//- Rotational relaxation factor
scalar relaxFactorR_;
//- Old translational relaxation factor
scalar oldRelaxFactorT_;
//- Old rotational relaxation factor
scalar oldRelaxFactorR_;
//- Accelerations from previous iterations
// A_[2] is the old value, A_[1] old old, A_[0] old old old
List<vector> A_;
List<vector> OmegaDot_;
//- Previos iteration non relaxed accelerations
List<vector> An_;
List<vector> OmegaDotn_;
// External forces
//- Force driving the motion in inertial coord. sys.
dimensionedVector force_;
//- Moment driving the motion in inertial coord. sys.
dimensionedVector moment_;
//- Force driving the motion in relative coord. sys.
dimensionedVector forceRelative_;
//- Moment driving the motion in relative coord. sys.
dimensionedVector momentRelative_;
// Private Member Functions
//- Disallow default bitwise copy construct
sixDOFODE(const sixDOFODE&);
//- Disallow default bitwise assignment
void operator=(const sixDOFODE&);
protected:
// Protected Member Functions
//- Set ODE coefficients from position and rotation
virtual void setCoeffs() = 0;
//- Update Aitkens relaxation parameters
void aitkensRelaxation(const scalar min, const scalar max);
public:
//- Run-time type information
TypeName("sixDOFODE");
// Declare run-time constructor selection table
declareRunTimeSelectionTable
(
autoPtr,
sixDOFODE,
dictionary,
(const IOobject& io),
(io)
);
// Constructors
//- Construct from dictionary
sixDOFODE(const IOobject& io);
//- Return a clone, changing the name
virtual autoPtr<sixDOFODE> clone(const word& name) const = 0;
// Selectors
//- Return autoPtr to the selected sixDOFODE
static autoPtr<sixDOFODE> New(const IOobject& io);
// Destructor
virtual ~sixDOFODE();
// Member Functions
// Access to common data
//- Return mass
inline const dimensionedScalar& mass() const;
//- Return access to mass
inline dimensionedScalar& mass();
//- Return moment of inertia
inline const dimensionedDiagTensor& momentOfInertia() const;
//- Return access to moment of inertia
inline dimensionedDiagTensor& momentOfInertia();
//- Return equilibrium position of origin
inline const dimensionedVector& Xequilibrium() const;
//- Return access to equilibrium position of origin
inline dimensionedVector& Xequilibrium();
// Access to forces and moments
//- Return force in inertial coord. sys.
inline const dimensionedVector& force() const;
//- Return access to force in inertial coord. sys.
inline dimensionedVector& force();
//- Return moment in inertial coord. sys.
inline const dimensionedVector& moment() const;
//- Return access to moment in inertial coord. sys.
inline dimensionedVector& moment();
//- Return force in relative coord. sys.
inline const dimensionedVector& forceRelative() const;
//- Return access to force in relative coord. sys.
inline dimensionedVector& forceRelative();
//- Return moment in relative coord. sys.
inline const dimensionedVector& momentRelative() const;
//- Return access to moment in relative coord. sys.
inline dimensionedVector& momentRelative();
//- Return total force in inertial coord. sys.
inline dimensionedVector forceTotal() const;
//- Return total moment in inertial coord. sys.
inline dimensionedVector momentTotal() const;
//- Relax the force (acceleration) using Aitkens or fixed relaxation
void relaxAcceleration
(
const scalar minRelFactor,
const scalar maxRelFactor
);
// Virtual interface for 6DOF motion state
// Variables in relative coordinate system
//- Return displacement in translated coordinate system
// relative to spring equilibrium
virtual const dimensionedVector& Xrel() const = 0;
//- Return rotational velocity in relative coordinate system
virtual const dimensionedVector& omega() const = 0;
// Displacement and rotation in the absolute coordinate system
//- Return position of origin in absolute coordinate system
virtual dimensionedVector X() const = 0;
//- Return velocity of origin
virtual const dimensionedVector& U() const = 0;
//- Return average velocity of origin for the previous time-step
virtual const dimensionedVector& Uaverage() const = 0;
//- Return finite rotation
virtual const finiteRotation& rotation() const = 0;
//- Return rotational vector of body
virtual vector rotVector() const = 0;
//- Return rotation angle of body
virtual dimensionedScalar rotAngle() const = 0;
// Non-access control for setting state variables
//- Set position of origin
virtual void setXrel(const vector& x) = 0;
//- Set velocity of origin
virtual void setU(const vector& u) = 0;
//- Set rotational angle in relative coordinate system
virtual void setRotation(const HamiltonRodriguezRot& rot) = 0;
//- Set rotational velocity in relative coordinate system
virtual void setOmega(const vector& omega) = 0;
//- Set referent coordinate system to apply constraints
virtual void setReferentRotation
(
const HamiltonRodriguezRot& rot
) = 0;
// Average motion per time-step
//- Return average rotational vector of body
virtual vector rotVectorAverage() const = 0;
//- Return average rotational velocity in relative coordinate
// system
virtual const dimensionedVector& omegaAverage() const = 0;
//- Return average rotational velocity in absolute coordinate
// system
virtual const dimensionedVector&
omegaAverageAbsolute() const = 0;
// Rotations
//- Return rotation tensor to relative coordinate system
virtual tensor toRelative() const = 0;
//- Return rotation tensor to absolute coordinate system
virtual tensor toAbsolute() const = 0;
//- Return transformation tensor between new and previous
// rotation
virtual const tensor& rotIncrementTensor() const = 0;
//- Set ODE parameters from another ODE
virtual void setState(const sixDOFODE&) = 0;
// ODE parameters
//- Return number of equations
virtual label nEqns() const = 0;
//- Return access to coefficients
virtual scalarField& coeffs() = 0;
//- Return reference to coefficients
virtual const scalarField& coeffs() const = 0;
//- Evaluate derivatives
virtual void derivatives
(
const scalar x,
const scalarField& y,
scalarField& dydx
) const = 0;
//- Evaluate Jacobian
virtual void jacobian
(
const scalar x,
const scalarField& y,
scalarField& dfdx,
scalarSquareMatrix& dfdy
) const = 0;
//- Update ODE after the solution, advancing by delta
virtual void update(const scalar delta) = 0;
// Write control
//- writeData member function required by regIOobject
virtual bool writeData(Ostream&) const;
// Ostream operator
friend Ostream& operator<<(Ostream&, const sixDOFODE&);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "sixDOFODEI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //

View file

@ -21,74 +21,104 @@ License
You should have received a copy of the GNU General Public License
along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
translationODE
Description
Ordinary differential equation for three degrees of freedom
solid body translation
Author
Hrvoje Jasak
Dubravko Matijasevic
\*---------------------------------------------------------------------------*/
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
const Foam::dimensionedScalar& Foam::translationODE::mass() const
const Foam::dimensionedScalar& Foam::sixDOFODE::mass() const
{
return mass_;
}
const Foam::dimensionedVector& Foam::translationODE::Xrel() const
Foam::dimensionedScalar& Foam::sixDOFODE::mass()
{
return Xrel_;
return mass_;
}
Foam::dimensionedVector Foam::translationODE::X() const
const Foam::dimensionedDiagTensor& Foam::sixDOFODE::momentOfInertia() const
{
return xEquilibrium_ + Xrel();
return momentOfInertia_;
}
const Foam::dimensionedVector& Foam::translationODE::U() const
Foam::dimensionedDiagTensor& Foam::sixDOFODE::momentOfInertia()
{
return U_;
return momentOfInertia_;
}
const Foam::dimensionedVector& Foam::translationODE::Uold() const
const Foam::dimensionedVector& Foam::sixDOFODE::Xequilibrium() const
{
return Uold_;
return Xequilibrium_;
}
Foam::dimensionedVector Foam::translationODE::A() const
Foam::dimensionedVector& Foam::sixDOFODE::Xequilibrium()
{
return A(Xrel(), U());
return Xequilibrium_;
}
Foam::dimensionedVector Foam::translationODE::Uaverage() const
{
return 0.5*(U_ + Uold_);
}
const Foam::dimensionedVector& Foam::translationODE::force() const
const Foam::dimensionedVector& Foam::sixDOFODE::force() const
{
return force_;
}
Foam::dimensionedVector& Foam::translationODE::force()
Foam::dimensionedVector& Foam::sixDOFODE::force()
{
return force_;
}
const Foam::dimensionedVector& Foam::sixDOFODE::moment() const
{
return moment_;
}
Foam::dimensionedVector& Foam::sixDOFODE::moment()
{
return moment_;
}
const Foam::dimensionedVector& Foam::sixDOFODE::forceRelative() const
{
return forceRelative_;
}
Foam::dimensionedVector& Foam::sixDOFODE::forceRelative()
{
return forceRelative_;
}
const Foam::dimensionedVector& Foam::sixDOFODE::momentRelative() const
{
return momentRelative_;
}
Foam::dimensionedVector& Foam::sixDOFODE::momentRelative()
{
return momentRelative_;
}
Foam::dimensionedVector Foam::sixDOFODE::forceTotal() const
{
return force_ + (this->toAbsolute() & forceRelative_);
}
Foam::dimensionedVector Foam::sixDOFODE::momentTotal() const
{
return moment_ + (this->toAbsolute() & momentRelative_);
}
// ************************************************************************* //

View file

@ -1,202 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation, either version 3 of the License, or (at your
option) any later version.
foam-extend is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
General Public License for more details.
You should have received a copy of the GNU General Public License
along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
translationODE
Description
Ordinary differential equation for three degrees of freedom
solid body motion
Author
Hrvoje Jasak
Dubravko Matijasevic
\*---------------------------------------------------------------------------*/
#include "translationODE.H"
#include "objectRegistry.H"
#include "foamTime.H"
// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
// namespace Foam
// {
// defineTypeNameAndDebug(translationODE, 0);
// }
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
void Foam::translationODE::setCoeffs()
{
// Set ODE coefficients from position and rotation
// Linear displacement in relative coordinate system
{
const vector& xVal = Xrel_.value();
coeffs_[0] = xVal.x();
coeffs_[1] = xVal.y();
coeffs_[2] = xVal.z();
}
// Linear velocity in relative coordinate system
{
const vector& uVal = U_.value();
coeffs_[3] = uVal.x();
coeffs_[4] = uVal.y();
coeffs_[5] = uVal.z();
}
}
Foam::dimensionedVector Foam::translationODE::A
(
const dimensionedVector& xR,
const dimensionedVector& uR
) const
{
return
(
- (linSpringCoeffs_ & xR) // spring
- (linDampingCoeffs_ & uR) // damping
+ force()
)/mass_; // gravity
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
// Construct from components
Foam::translationODE::translationODE
(
const IOobject& io
)
:
IOdictionary(io),
mass_(lookup("mass")),
xEquilibrium_(lookup("equilibriumPosition")),
linSpringCoeffs_(lookup("linearSpring")),
linDampingCoeffs_(lookup("linearDamping")),
Xrel_(lookup("Xrel")),
U_(lookup("U")),
Uold_(lookup("Uold")),
force_(lookup("force")),
coeffs_(6, 0.0)
{
setCoeffs();
}
// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
Foam::translationODE::~translationODE()
{}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
void Foam::translationODE::derivatives
(
const scalar x,
const scalarField& y,
scalarField& dydx
) const
{
// Set the derivatives for displacement
dydx[0] = y[3];
dydx[1] = y[4];
dydx[2] = y[5];
dimensionedVector curX("curX", dimLength, vector(y[0], y[1], y[2]));
dimensionedVector curU("curU", dimLength/dimTime, vector(y[3], y[4], y[5]));
const vector& accel = A(curX, curU).value();
dydx[3] = accel.x();
dydx[4] = accel.y();
dydx[5] = accel.z();
}
void Foam::translationODE::jacobian
(
const scalar x,
const scalarField& y,
scalarField& dfdx,
scalarSquareMatrix& dfdy
) const
{
notImplemented("translationODE::jacobian(...) const");
}
void Foam::translationODE::update(const scalar delta)
{
// Update position
vector& Xval = Xrel_.value();
Xval.x() = coeffs_[0];
Xval.y() = coeffs_[1];
Xval.z() = coeffs_[2];
// Update velocity
Uold_ = U_;
vector& Uval = U_.value();
Uval.x() = coeffs_[3];
Uval.y() = coeffs_[4];
Uval.z() = coeffs_[5];
}
bool Foam::translationODE::writeData(Ostream& os) const
{
os << *this;
return os.good();
}
// * * * * * * * * * * * * * * * IOstream Operators * * * * * * * * * * * * //
Foam::Ostream& Foam::operator<<(Ostream& os, const translationODE& sds)
{
os.writeKeyword("mass") << sds.mass_ << token::END_STATEMENT << endl;
os.writeKeyword("equilibriumPosition") << sds.xEquilibrium_
<< token::END_STATEMENT << endl;
os.writeKeyword("linearSpring") << sds.linSpringCoeffs_
<< token::END_STATEMENT << endl;
os.writeKeyword("linearDamping") << sds.linDampingCoeffs_
<< token::END_STATEMENT << endl;
os.writeKeyword("Xrel") << sds.Xrel() << token::END_STATEMENT << endl;
os.writeKeyword("U") << sds.U() << token::END_STATEMENT << endl;
os.writeKeyword("Uold") << tab << sds.Uold() << token::END_STATEMENT << nl;
os.writeKeyword("force") << sds.force() << token::END_STATEMENT << endl;
return os;
}
// ************************************************************************* //

View file

@ -1,251 +0,0 @@
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | foam-extend: Open Source CFD
\\ / O peration | Version: 4.0
\\ / A nd | Web: http://www.foam-extend.org
\\/ M anipulation | For copyright notice see file Copyright
-------------------------------------------------------------------------------
License
This file is part of foam-extend.
foam-extend is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by the
Free Software Foundation, either version 3 of the License, or (at your
option) any later version.
foam-extend is distributed in the hope that it will be useful, but
WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
General Public License for more details.
You should have received a copy of the GNU General Public License
along with foam-extend. If not, see <http://www.gnu.org/licenses/>.
Class
translationODE
Description
Ordinary differential equation for three degrees of freedom solid
body translation
Author
Hrvoje Jasak
Dubravko Matijasevic
SourceFiles
translationODEI.H
translationODE.C
\*---------------------------------------------------------------------------*/
#ifndef translationODE_H
#define translationODE_H
#include "ODE.H"
#include "IOdictionary.H"
#include "dimensionedTypes.H"
#include "dimensionedDiagTensor.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
/*---------------------------------------------------------------------------*\
Class translationODE Declaration
\*---------------------------------------------------------------------------*/
class translationODE
:
public IOdictionary,
public ODE
{
// Private data
// Body data
//- Mass
dimensionedScalar mass_;
// Platform variables
//- Spring equilibrium position for translation
const dimensionedVector xEquilibrium_;
//- Linear spring coeffs
const dimensionedDiagTensor linSpringCoeffs_;
//- Linear damping coeffs
const dimensionedDiagTensor linDampingCoeffs_;
// Body position and rotation variables
//- Displacement relative to spring equilibrium
dimensionedVector Xrel_;
//- Velocity of mass centroid
dimensionedVector U_;
//- Velocity of mass centroid at previous time-step
dimensionedVector Uold_;
// External forces
//- Force driving the motion
dimensionedVector force_;
//- ODE coefficients
scalarField coeffs_;
// Private Member Functions
//- Disallow default bitwise copy construct
translationODE(const translationODE&);
//- Disallow default bitwise assignment
void operator=(const translationODE&);
//- Set ODE coefficients from position and rotation
inline void setCoeffs();
// Variables in relative coordinate system (solved for)
//- Return acceleration in relative coordinate system
// given current values of relative displacement and velocity
dimensionedVector A
(
const dimensionedVector& xR,
const dimensionedVector& uR
) const;
public:
// //- Runtime type information
// TypeName("translationODE");
// Constructors
//- Construct from dictionary
translationODE(const IOobject& io);
// Destructor
virtual ~translationODE();
// Member Functions
//- Return mass
inline const dimensionedScalar& mass() const;
// Variables in relative coordinate system (solved for)
//- Return displacement in relative coordinate system
inline const dimensionedVector& Xrel() const;
// Displacement and rotation in the absolute coordinate system
//- Return position of origin in absolute coordinate system
inline dimensionedVector X() const;
//- Return velocity of origin
inline const dimensionedVector& U() const;
//- Return velocity of origin for the previous time-step
inline const dimensionedVector& Uold() const;
//- Return acceleration of origin
inline dimensionedVector A() const;
// Average motion per time-step
//- Return average velocity of origin
inline dimensionedVector Uaverage() const;
// Force
//- Return force
inline const dimensionedVector& force() const;
//- Return access to force
inline dimensionedVector& force();
// ODE parameters
//- Return number of equations
virtual label nEqns() const
{
return 6;
}
//- Return reference to interpolation coefficients
virtual scalarField& coeffs()
{
return coeffs_;
}
//- Return reference to interpolation coefficients
virtual const scalarField& coeffs() const
{
return coeffs_;
}
//- Return derivatives
virtual void derivatives
(
const scalar x,
const scalarField& y,
scalarField& dydx
) const;
//- Return Jacobian
virtual void jacobian
(
const scalar x,
const scalarField& y,
scalarField& dfdx,
scalarSquareMatrix& dfdy
) const;
//- Update ODE after the solution, advancing by delta
virtual void update(const scalar delta);
//- WriteData member function required by regIOobject
bool writeData(Ostream&) const;
// Ostream operator
friend Ostream& operator<<(Ostream&, const translationODE&);
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#include "translationODEI.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //