From 488ccffca96ef9f142adb47727c3d2161ff0bc30 Mon Sep 17 00:00:00 2001 From: Vuko Vukcevic Date: Tue, 7 Mar 2017 10:28:43 +0100 Subject: [PATCH] Implementation of constraints, part 2 --- src/ODE/Make/files | 5 + .../constantAngularAcceleration.C | 29 ++- .../constantAngularAcceleration.H | 8 +- ...ionConstraint.C => rotationalConstraint.C} | 2 +- ...ionConstraint.H => rotationalConstraint.H} | 8 +- .../constantTranslationalAcceleration.C | 84 ++++++++ .../constantTranslationalAcceleration.H | 126 +++++++++++ .../translationalConstraint.C | 99 +++++++++ .../translationalConstraint.H | 200 ++++++++++++++++++ .../sixDOF/geometricSixDOF/geometricSixDOF.C | 50 +++-- .../sixDOF/geometricSixDOF/geometricSixDOF.H | 14 +- 11 files changed, 588 insertions(+), 37 deletions(-) rename src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/{rotationConstraint.C => rotationalConstraint.C} (98%) rename src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/{rotationConstraint.H => rotationalConstraint.H} (95%) create mode 100644 src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.C create mode 100644 src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.H create mode 100644 src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.C create mode 100644 src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.H diff --git a/src/ODE/Make/files b/src/ODE/Make/files index dfcfcc898..f6d20418f 100644 --- a/src/ODE/Make/files +++ b/src/ODE/Make/files @@ -16,6 +16,11 @@ sixDOF = sixDOF $(sixDOF)/finiteRotation/finiteRotation.C $(sixDOF)/sixDOFqODE/sixDOFqODE.C +$(sixDOF)/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.C +$(sixDOF)/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.C + +$(sixDOF)/constraints/translationalConstraints/translationalConstraint/translationalConstraint.C + $(sixDOF)/sixDOFODE/sixDOFODE.C $(sixDOF)/sixDOFODE/newSixDOFODE.C $(sixDOF)/quaternionSixDOF/quaternionSixDOF.C diff --git a/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.C b/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.C index d431d18e0..6eccb0e35 100644 --- a/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.C +++ b/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.C @@ -63,15 +63,36 @@ Foam::constantAngularAcceleration::~constantAngularAcceleration() // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // -const Foam::vector -Foam::constantAngularAcceleration::matrixContribution() const +Foam::vector Foam::constantAngularAcceleration::matrixContribution +( + const scalar, + const tensor& toRelative +) const { + vector mc; + + if (inGlobal_) + { + // Constraint is given in global (inertial) coordinate system, transform + // it to local + mc = toRelative & dir_; + } + else + { + // Constraint already in local (body) coordinate system + mc = dir_; + } + + return mc; } -const Foam::vector -Foam::constantAngularAcceleration::sourceContribution() const +Foam::scalar Foam::constantAngularAcceleration::sourceContribution +( + const scalar +) const { + return alpha_; } diff --git a/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.H b/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.H index b666226b6..6f1a7494b 100644 --- a/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.H +++ b/src/ODE/sixDOF/constraints/rotationalConstraints/constantAngularAcceleration/constantAngularAcceleration.H @@ -105,10 +105,14 @@ public: // Constraint specific functions //- Return matrix contribution defined by constraint, f(t) - const vector matrixContribution() const = 0; + virtual vector matrixContribution + ( + const scalar, + const tensor& toRelative + ) const; //- Return source contribution defined by constraint, a(t) - const vector sourceContribution() const = 0; + virtual scalar sourceContribution(const scalar) const; }; diff --git a/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.C b/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.C similarity index 98% rename from src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.C rename to src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.C index 8967a0a02..1bbf07157 100644 --- a/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.C +++ b/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.C @@ -65,7 +65,7 @@ Foam::autoPtr Foam::rotationalConstraint::New const word constraintType(dict.lookup("type")); wordConstructorTable::iterator cstrIter = - wordConstructorTablePtr_->find(constraintTypeType); + wordConstructorTablePtr_->find(constraintType); if (cstrIter == wordConstructorTablePtr_->end()) { diff --git a/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.H b/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.H similarity index 95% rename from src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.H rename to src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.H index 9e85b8b3f..f897cbac6 100644 --- a/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationConstraint.H +++ b/src/ODE/sixDOF/constraints/rotationalConstraints/rotationalConstraint/rotationalConstraint.H @@ -178,10 +178,14 @@ public: // Constraint specific functions //- Return matrix contribution defined by constraint, f(t) - const vector matrixContribution() const = 0; + virtual vector matrixContribution + ( + const scalar t, + const tensor& toRelative + ) const = 0; //- Return source contribution defined by constraint, a(t) - const vector sourceContribution() const = 0; + virtual scalar sourceContribution(const scalar t) const = 0; }; diff --git a/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.C b/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.C new file mode 100644 index 000000000..9289f6c1a --- /dev/null +++ b/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.C @@ -0,0 +1,84 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +\*---------------------------------------------------------------------------*/ + +#include "constantTranslationalAcceleration.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(constantTranslationalAcceleration, 0); + addToRunTimeSelectionTable + ( + translationalConstraint, + constantTranslationalAcceleration, + word + ); +} + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::constantTranslationalAcceleration::constantTranslationalAcceleration +( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict +) +: + translationalConstraint(name, sixDOFODE, dict), + dir_(dict.lookup("constraintDirection")), + a_(readScalar(dict.lookup("translationalAcceleration"))) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::constantTranslationalAcceleration::~constantTranslationalAcceleration() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +Foam::vector Foam::constantTranslationalAcceleration::matrixContribution +( + const scalar, + const tensor& +) const +{ + return dir_; +} + + +Foam::scalar Foam::constantTranslationalAcceleration::sourceContribution +( + const scalar +) const +{ + return a_; +} + + +// ************************************************************************* // diff --git a/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.H b/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.H new file mode 100644 index 000000000..2670eda97 --- /dev/null +++ b/src/ODE/sixDOF/constraints/translationalConstraints/constantTranslationalAcceleration/constantTranslationalAcceleration.H @@ -0,0 +1,126 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +Class + Foam::constantTranslationalAcceleration + +Description + Translational constraint defined by constant translational acceleration + + g(vDot) = vDot - a = 0. + +Author + Viktor Pandza, FSB Zagreb. All rights reserved. + Vuko Vukcevic, FSB Zagreb. All rights reserved. + +SourceFiles + constantTranslationalAcceleration.C + +\*---------------------------------------------------------------------------*/ + +#ifndef constantTranslationalAcceleration_H +#define constantTranslationalAcceleration_H + +#include "translationalConstraint.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class constantTranslationalAcceleration Declaration +\*---------------------------------------------------------------------------*/ + +class constantTranslationalAcceleration +: + public translationalConstraint +{ + // Private data + + //- Direction of the constraint (unit vector) + const vector dir_; + + //- Constant value of the translational acceleration + const scalar a_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + constantTranslationalAcceleration + ( + const constantTranslationalAcceleration& + ); + + //- Disallow default bitwise assignment + void operator=(const constantTranslationalAcceleration&); + + +public: + + //- Runtime type information + TypeName("constantTranslationalAcceleration"); + + + // Constructors + + //- Construct from dictionary + constantTranslationalAcceleration + ( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict + ); + + + // Destructor + + virtual ~constantTranslationalAcceleration(); + + + // Member Functions + + // Constraint specific functions + + //- Return matrix contribution defined by constraint, f(t) + virtual vector matrixContribution + ( + const scalar, + const tensor& + ) const; + + //- Return source contribution defined by constraint, a(t) + virtual scalar sourceContribution(const scalar) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.C b/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.C new file mode 100644 index 000000000..e2b2d805a --- /dev/null +++ b/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.C @@ -0,0 +1,99 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +\*---------------------------------------------------------------------------*/ + +#include "translationalConstraint.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(translationalConstraint, 0); + defineRunTimeSelectionTable(translationalConstraint, word); +} + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::translationalConstraint::translationalConstraint +( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict +) +: + name_(name), + sixDOFODE_(sixDOFODE) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::translationalConstraint::~translationalConstraint() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +Foam::autoPtr Foam::translationalConstraint::New +( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict +) +{ + const word constraintType(dict.lookup("type")); + + wordConstructorTable::iterator cstrIter = + wordConstructorTablePtr_->find(constraintType); + + if (cstrIter == wordConstructorTablePtr_->end()) + { + FatalErrorIn + ( + "translationalConstraint::New" + "\n(" + "\n const word& name," + "\n const sixDOFODE& sixDOFODE," + "\n const dictionary& dict," + "\n)" + ) << "Unknown translation constraint type: " << constraintType + << endl << endl + << "Valid translation constraint types are: " << endl + << wordConstructorTablePtr_->sortedToc() + << exit(FatalError); + } + + return autoPtr + ( + cstrIter() + ( + name, + sixDOFODE, + dict + ) + ); +} + + +// ************************************************************************* // diff --git a/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.H b/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.H new file mode 100644 index 000000000..32d509ec6 --- /dev/null +++ b/src/ODE/sixDOF/constraints/translationalConstraints/translationalConstraint/translationalConstraint.H @@ -0,0 +1,200 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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 . + +Class + Foam::translationalConstraint + +Description + Abstract base class containing interface for translational constraints used + within sixDOFODE classes. + + The constraint is implicitly defined as: + + g(vDot, t) = f(t)*vDot + a(t) = 0, + + where vDot is translational acceleration. + + Interface provides all the necessary data for inserting the constraint into + the resulting linear system via Lagrangian multipliers: + 1. matrixContribution() - corresponding to f(t) (prefactor multiplying + vDot), to be inserted into the matrix. + 2. sourceContribution() - corresponding to a(t), to be inserted into right + hand side vector. + +Author + Viktor Pandza, FSB Zagreb. All rights reserved. + Vuko Vukcevic, FSB Zagreb. All rights reserved. + +SourceFiles + translationalConstraint.C + +\*---------------------------------------------------------------------------*/ + +#ifndef translationalConstraint_H +#define translationalConstraint_H + +#include "sixDOFODE.H" +#include "typeInfo.H" +#include "runTimeSelectionTables.H" +#include "autoPtr.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class translationalConstraint Declaration +\*---------------------------------------------------------------------------*/ + +class translationalConstraint +{ + // Private data + + //- Name of the constraint + const word name_; + + //- Reference to underlying sixDOFODE + const sixDOFODE& sixDOFODE_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + translationalConstraint(const translationalConstraint&); + + //- Disallow default bitwise assignment + void operator=(const translationalConstraint&); + + +public: + + //- Runtime type information + TypeName("translationalConstraint"); + + + // Declare run-time constructor selection table + + declareRunTimeSelectionTable + ( + autoPtr, + translationalConstraint, + word, + ( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict + ), + (name, sixDOFODE, dict) + ); + + + //- Class used for the read-construction of + // PtrLists of translationalConstraint + class iNew + { + const sixDOFODE& sixDOFODE_; + + public: + + iNew(const sixDOFODE& sixDOFODE) + : + sixDOFODE_(sixDOFODE) + {} + + autoPtr operator()(Istream& is) const + { + word name(is); + dictionary dict(is); + return translationalConstraint::New(name, sixDOFODE_, dict); + } + }; + + + // Constructors + + //- Construct from dictionary + translationalConstraint + ( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict + ); + + + // Selectors + + //- Return a reference to the selected translationalConstraint + static autoPtr New + ( + const word& name, + const sixDOFODE& sixDOFODE, + const dictionary& dict + ); + + + // Destructor + + virtual ~translationalConstraint(); + + + // Member Functions + + // Access functions + + //- Return const reference to name of the constraint + const word& name() const + { + return name_; + } + + //- Return const reference to underlying sixDOFODE object + const sixDOFODE& sixDOF() const + { + return sixDOFODE_; + } + + + // Constraint specific functions + + //- Return matrix contribution defined by constraint, f(t) + virtual vector matrixContribution + ( + const scalar t, + const tensor& toRelative + ) const = 0; + + //- Return source contribution defined by constraint, a(t) + virtual scalar sourceContribution(const scalar t) const = 0; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.C b/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.C index c99d2c45e..fe27835e2 100644 --- a/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.C +++ b/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.C @@ -72,7 +72,8 @@ Foam::dimensionedVector Foam::geometricSixDOF::A ( const dimensionedVector& xR, const dimensionedVector& uR, - const tensor& R + const tensor& R, + const scalar t ) const { // Create a scalar square matrix representing Newton equations with @@ -87,7 +88,7 @@ Foam::dimensionedVector Foam::geometricSixDOF::A ( force() // External force - (linSpringCoeffs() & xR) // Spring force - - (linDampingCoeffs() & uR); // Damping force + - (linDampingCoeffs() & uR) // Damping force ); const vector& efVal = explicitForcing.value(); const scalar& m = mass().value(); @@ -103,23 +104,23 @@ Foam::dimensionedVector Foam::geometricSixDOF::A forAll(translationalConstraints_, tcI) { // Get reference to current constraint - const translationalConstraint& curTc = translationalConstraint_[tcI]; + const translationalConstraint& curTc = translationalConstraints_[tcI]; // Get matrix contribution from constraint - const vector mc = curTc.matrixContribution(); + const vector mc = curTc.matrixContribution(t, R); // Get matrix index const label index = tcI + 3; // Insert contributions into the matrix - forAll(lower, dirI) + forAll(mc, dirI) { M[dirI][index] = mc[dirI]; M[index][dirI] = mc[dirI]; } // Insert source contribution (remainder of the constraint function) - rhs[index] = curTc.sourceContribution(); + rhs[index] = curTc.sourceContribution(t); } // Solve the matrix using LU decomposition. Note: solution is in the rhs and @@ -132,7 +133,7 @@ Foam::dimensionedVector Foam::geometricSixDOF::A ( "A", force().dimensions()/mass().dimensions(), - vector(rhs.x(), rhs.y(), rhs.z()) + vector(rhs[0], rhs[1], rhs[2]) ); } @@ -140,7 +141,8 @@ Foam::dimensionedVector Foam::geometricSixDOF::A Foam::dimensionedVector Foam::geometricSixDOF::OmegaDot ( const tensor& R, - const dimensionedVector& omega + const dimensionedVector& omega, + const scalar t ) const { // Create a scalar square matrix representing Euler equations with @@ -149,11 +151,14 @@ Foam::dimensionedVector Foam::geometricSixDOF::OmegaDot scalarField rhs(rotationalConstraints_.size() + 3, 0.0); scalarSquareMatrix J(rhs.size(), 0.0); + // Get current inertial-to-local transformation + const dimensionedTensor RT("RT", dimless, R.T()); + // Insert moment of inertia and explicit forcing into the system const dimensionedVector explicitForcing ( E(omega) // Euler part - + (dimensionedTensor("R", dimless, R) & moment()) // External torque + + (RT & moment()) // External torque ); const vector& efVal = explicitForcing.value(); const diagTensor& I = momentOfInertia().value(); @@ -172,20 +177,20 @@ Foam::dimensionedVector Foam::geometricSixDOF::OmegaDot const rotationalConstraint& curRc = rotationalConstraints_[rcI]; // Get matrix contribution from the constraint - const vector mc = curRc.matrixContribution(); + const vector mc = curRc.matrixContribution(t, RT.value()); // Get matrix index const label index = rcI + 3; // Insert contributions into the matrix - forAll(upper, dirI) + forAll(mc, dirI) { J[dirI][index] = mc[dirI]; J[index][dirI] = mc[dirI]; } // Insert source contribution (remainder of the constraint function) - rhs[index] = curRc.sourceContribution(); + rhs[index] = curRc.sourceContribution(t); } // Solve the matrix using LU decomposition. Note: solution is in the rhs and @@ -198,7 +203,7 @@ Foam::dimensionedVector Foam::geometricSixDOF::OmegaDot ( "OmegaDot", moment().dimensions()/momentOfInertia().dimensions(), - vector(rhs.x(), rhs.y(), rhs.z()) + vector(rhs[0], rhs[1], rhs[2]) ); } @@ -307,8 +312,8 @@ void Foam::geometricSixDOF::setState(const sixDOFODE& sd) // HJ, 23/Mar/2015 coeffs_ = gsd.coeffs_; - translationalConstraints = gsd.translationalConstraints_; - rotationalConstraints = gsd.rotationalConstraints_; + translationalConstraints_ = gsd.translationalConstraints_; + rotationalConstraints_ = gsd.rotationalConstraints_; } @@ -365,10 +370,10 @@ Foam::geometricSixDOF::geometricSixDOF(const IOobject& io) // Read translation constraints if they are present if (dict().found("translationalConstraints")) { - PtrList tcList + PtrList tcList ( dict().lookup("translationalConstraints"), - translationalConstraints::iNew(*this) + translationalConstraint::iNew(*this) ); translationalConstraints_.transfer(tcList); } @@ -376,10 +381,10 @@ Foam::geometricSixDOF::geometricSixDOF(const IOobject& io) // Read rotation constraints if they are present if (dict().found("rotationalConstraints")) { - PtrList tcList + PtrList tcList ( dict().lookup("rotationalConstraints"), - rotationalConstraints::iNew(*this) + rotationalConstraint::iNew(*this) ); rotationalConstraints_.transfer(tcList); } @@ -404,7 +409,7 @@ Foam::geometricSixDOF::geometricSixDOF coeffs_(gsd.coeffs_), translationalConstraints_(gsd.translationalConstraints_), - rotationalConstraints_(gsd.rotationalConstraints_), + rotationalConstraints_(gsd.rotationalConstraints_) {} @@ -511,7 +516,7 @@ void Foam::geometricSixDOF::derivatives const tensor curRot = (rotation_ & expMap(rotIncrementVector)); // Calculate translational acceleration using current rotation - const vector accel = A(curX, curU, curRot).value(); + const vector accel = A(curX, curU, curRot, x).value(); // Set the derivatives for velocity dydx[3] = accel.x(); @@ -528,7 +533,7 @@ void Foam::geometricSixDOF::derivatives ); // Calculate rotational acceleration using current rotation - const vector omegaDot = OmegaDot(curRot, curOmega).value(); + const vector omegaDot = OmegaDot(curRot, curOmega, x).value(); dydx[6] = omegaDot.x(); dydx[7] = omegaDot.y(); @@ -569,7 +574,6 @@ void Foam::geometricSixDOF::update(const scalar delta) Uval.z() = coeffs_[5]; // Constrain velocity and re-set coefficients - constrainTranslation(Uval); coeffs_[3] = Uval.x(); coeffs_[4] = Uval.y(); coeffs_[5] = Uval.z(); diff --git a/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.H b/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.H index 7d7d8b7f3..494652802 100644 --- a/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.H +++ b/src/ODE/sixDOF/geometricSixDOF/geometricSixDOF.H @@ -116,7 +116,7 @@ class geometricSixDOF PtrList translationalConstraints_; //- List of rotational constraints - Ptrlist rotationalConstraints_; + PtrList rotationalConstraints_; // Private Member Functions @@ -131,21 +131,25 @@ class geometricSixDOF // Variables in relative coordinate system (solved for) //- Return acceleration in relative coordinate system - // given current values of relative displacement and velocity + // given current values of relative displacement, velocity, + // orientation (rotation tensor) and time dimensionedVector A ( const dimensionedVector& xR, const dimensionedVector& uR, - const tensor& R + const tensor& R, + const scalar t ) const; //- Return rotational acceleration in relative coordinate system - // given current values for relative rotational velocity + // given current values for relative rotational velocity, + // orientation (rotation tensor) and time dimensionedVector OmegaDot ( const tensor& R, - const dimensionedVector& omega + const dimensionedVector& omega, + const scalar t ) const; //- Return the Euler part of moment equation