diff --git a/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C b/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C
index aa1488badee267f646100ac76d34dd0514e94d1b..d2f25bb46f5dc4e39532aeac4e08a37a2aa37b59 100644
--- a/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C
+++ b/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C
@@ -133,7 +133,11 @@ int main(int argc, char *argv[])
{
Info<< "Time = " << runTime.timeName() << endl;
- mesh.update();
+ for (int i = 0; i<2; i++)
+ {
+ mesh.update();
+ }
+
mesh.checkMesh(true);
if (checkAMI)
diff --git a/src/sixDoFRigidBodyMotion/Make/files b/src/sixDoFRigidBodyMotion/Make/files
index 8e825f0a010760696e12c38e24e7e3f5f5bf36dc..1184c6ad72d33ab7e051c816f76e94b9d53e0fc0 100644
--- a/src/sixDoFRigidBodyMotion/Make/files
+++ b/src/sixDoFRigidBodyMotion/Make/files
@@ -31,4 +31,10 @@ sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C
sixDoFRigidBodyMotionSolver/externalPointEdgePoint.C
sixDoFRigidBodyMotionSolver/pointPatchDist.C
+sixDoFSolvers/sixDoFSolver/sixDoFSolver.C
+sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C
+sixDoFSolvers/symplectic/symplectic.C
+sixDoFSolvers/CrankNicolson/CrankNicolson.C
+sixDoFSolvers/Newmark/Newmark.C
+
LIB = $(FOAM_LIBBIN)/libsixDoFRigidBodyMotion
diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C
index a22e4e586551af465a9fc12833c7112ca161d613..7f51b9551d45ea72eb1fd0affbe6d7035bea5655 100644
--- a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C
+++ b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C
@@ -224,11 +224,6 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
f.calcForcesMoment();
- // Patch force data is valid for the current positions, so
- // calculate the forces on the motion object from this data, then
- // update the positions
- motion_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value());
-
// Get the forces on the patch faces at the current positions
if (lookupGravity_ == 1)
@@ -242,7 +237,7 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
// scalar ramp = min(max((t.value() - 5)/10, 0), 1);
scalar ramp = 1.0;
- motion_.updateAcceleration
+ motion_.update
(
firstIter,
ramp*(f.forceEff() + motion_.mass()*g_),
diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C
index 03bdfcd77935d87eb114edc2354c58536cd54ff0..d567e3e04aa9249b364552a5ee5576e35bc8b6a3 100644
--- a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C
+++ b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C
@@ -160,8 +160,6 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
firstIter = true;
}
- motion_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value());
-
vector gravity = vector::zero;
if (db().foundObject<uniformDimensionedVectorField>("g"))
@@ -173,7 +171,7 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs()
}
// Do not modify the accelerations, except with gravity, where available
- motion_.updateAcceleration
+ motion_.update
(
firstIter,
gravity*motion_.mass(),
diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C
index 9bb33ecbf9f677df7662840903c012a779be9b2f..e7a8b92d339880ba8820abd4946cd0bd8abf96f6 100644
--- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C
+++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C
@@ -24,6 +24,7 @@ License
\*---------------------------------------------------------------------------*/
#include "sixDoFRigidBodyMotion.H"
+#include "sixDoFSolver.H"
#include "septernion.H"
// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
@@ -84,7 +85,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion()
momentOfInertia_(diagTensor::one*VSMALL),
aRelax_(1.0),
aDamp_(1.0),
- report_(false)
+ report_(false),
+ solver_(NULL)
{}
@@ -121,7 +123,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion
momentOfInertia_(dict.lookup("momentOfInertia")),
aRelax_(dict.lookupOrDefault<scalar>("accelerationRelaxation", 1.0)),
aDamp_(dict.lookupOrDefault<scalar>("accelerationDamping", 1.0)),
- report_(dict.lookupOrDefault<Switch>("report", false))
+ report_(dict.lookupOrDefault<Switch>("report", false)),
+ solver_(sixDoFSolver::New(dict.subDict("solver"), *this))
{
addRestraints(dict);
@@ -171,6 +174,12 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion
{}
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::sixDoFRigidBodyMotion::~sixDoFRigidBodyMotion()
+{}
+
+
// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
void Foam::sixDoFRigidBodyMotion::addRestraints
@@ -257,50 +266,35 @@ void Foam::sixDoFRigidBodyMotion::addConstraints
}
-void Foam::sixDoFRigidBodyMotion::updatePosition
+void Foam::sixDoFRigidBodyMotion::updateAcceleration
(
- bool firstIter,
- scalar deltaT,
- scalar deltaT0
+ const vector& fGlobal,
+ const vector& tauGlobal
)
{
- if (Pstream::master())
- {
- if (firstIter)
- {
- // First simplectic step:
- // Half-step for linear and angular velocities
- // Update position and orientation
-
- v() = tConstraints_ & (v0() + aDamp_*0.5*deltaT0*a());
- pi() = rConstraints_ & (pi0() + aDamp_*0.5*deltaT0*tau());
-
- centreOfRotation() = centreOfRotation0() + deltaT*v();
- }
- else
- {
- // For subsequent iterations use Crank-Nicolson
+ static bool first = false;
- v() = tConstraints_
- & (v0() + aDamp_*0.5*deltaT*(a() + motionState0_.a()));
- pi() = rConstraints_
- & (pi0() + aDamp_*0.5*deltaT*(tau() + motionState0_.tau()));
+ // Save the previous iteration accelerations for relaxation
+ vector aPrevIter = a();
+ vector tauPrevIter = tau();
- centreOfRotation() =
- centreOfRotation0() + 0.5*deltaT*(v() + motionState0_.v());
- }
+ // Calculate new accelerations
+ a() = fGlobal/mass_;
+ tau() = (Q().T() & tauGlobal);
+ applyRestraints();
- // Correct orientation
- Tuple2<tensor, vector> Qpi = rotate(Q0(), pi(), deltaT);
- Q() = Qpi.first();
- pi() = rConstraints_ & Qpi.second();
+ // Relax accelerations on all but first iteration
+ if (!first)
+ {
+ a() = aRelax_*a() + (1 - aRelax_)*aPrevIter;
+ tau() = aRelax_*tau() + (1 - aRelax_)*tauPrevIter;
}
- Pstream::scatter(motionState_);
+ first = false;
}
-void Foam::sixDoFRigidBodyMotion::updateAcceleration
+void Foam::sixDoFRigidBodyMotion::update
(
bool firstIter,
const vector& fGlobal,
@@ -309,44 +303,9 @@ void Foam::sixDoFRigidBodyMotion::updateAcceleration
scalar deltaT0
)
{
- static bool first = false;
-
if (Pstream::master())
{
- // Save the previous iteration accelerations for relaxation
- vector aPrevIter = a();
- vector tauPrevIter = tau();
-
- // Calculate new accelerations
- a() = fGlobal/mass_;
- tau() = (Q().T() & tauGlobal);
- applyRestraints();
-
- // Relax accelerations on all but first iteration
- if (!first)
- {
- a() = aRelax_*a() + (1 - aRelax_)*aPrevIter;
- tau() = aRelax_*tau() + (1 - aRelax_)*tauPrevIter;
- }
- first = false;
-
- if (firstIter)
- {
- // Second simplectic step:
- // Complete update of linear and angular velocities
-
- v() += tConstraints_ & aDamp_*0.5*deltaT*a();
- pi() += rConstraints_ & aDamp_*0.5*deltaT*tau();
- }
- else
- {
- // For subsequent iterations use Crank-Nicolson
-
- v() = tConstraints_
- & (v0() + aDamp_*0.5*deltaT*(a() + motionState0_.a()));
- pi() = rConstraints_
- & (pi0() + aDamp_*0.5*deltaT*(tau() + motionState0_.tau()));
- }
+ solver_->solve(firstIter, fGlobal, tauGlobal, deltaT, deltaT0);
if (report_)
{
@@ -358,7 +317,6 @@ void Foam::sixDoFRigidBodyMotion::updateAcceleration
}
-
void Foam::sixDoFRigidBodyMotion::status() const
{
Info<< "6-DoF rigid body motion" << nl
diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H
index 70f2e7674f60b37e885eb6125f2878a40857918b..a91f17e23dcba0650fae3c8bcde5830b6db32714 100644
--- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H
+++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H
@@ -25,28 +25,16 @@ Class
Foam::sixDoFRigidBodyMotion
Description
- Six degree of freedom motion for a rigid body. Angular momentum
- stored in body fixed reference frame. Reference orientation of
- the body (where Q = I) must align with the cartesian axes such
- that the Inertia tensor is in principle component form.
+ Six degree of freedom motion for a rigid body.
- Symplectic motion as per:
+ Angular momentum stored in body fixed reference frame. Reference
+ orientation of the body (where Q = I) must align with the cartesian axes
+ such that the Inertia tensor is in principle component form. Can add
+ restraints (e.g. a spring) and constraints (e.g. motion may only be on a
+ plane).
- title = {Symplectic splitting methods for rigid body molecular dynamics},
- publisher = {AIP},
- year = {1997},
- journal = {The Journal of Chemical Physics},
- volume = {107},
- number = {15},
- pages = {5840-5851},
- url = {http://link.aip.org/link/?JCP/107/5840/1},
- doi = {10.1063/1.474310}
-
- Changes to Crank-Nicolson integration if subsequent iterations to handle
- implicit coupling between pressure and motion.
-
- Can add restraints (e.g. a spring)
- and constraints (e.g. motion may only be on a plane).
+ The time-integrator for the motion is run-time selectable with options for
+ symplectic (explicit), Crank-Nicolson and Newmark schemes.
SourceFiles
sixDoFRigidBodyMotionI.H
@@ -69,12 +57,17 @@ SourceFiles
namespace Foam
{
+// Forward declarations
+class sixDoFSolver;
+
/*---------------------------------------------------------------------------*\
Class sixDoFRigidBodyMotion Declaration
\*---------------------------------------------------------------------------*/
class sixDoFRigidBodyMotion
{
+ friend class sixDoFSolver;
+
// Private data
//- Motion state data object
@@ -120,6 +113,9 @@ class sixDoFRigidBodyMotion
//- Switch to turn reporting of motion data on and off
Switch report_;
+ //- Motion solver
+ autoPtr<sixDoFSolver> solver_;
+
// Private Member Functions
@@ -147,6 +143,9 @@ class sixDoFRigidBodyMotion
//- Apply the restraints to the object
void applyRestraints();
+ //- Update and relax accelerations from the force and torque
+ void updateAcceleration(const vector& fGlobal, const vector& tauGlobal);
+
// Access functions retained as private because of the risk of
// confusion over what is a body local frame vector and what is global
@@ -179,21 +178,6 @@ class sixDoFRigidBodyMotion
//- Return the current torque
inline const vector& tau() const;
- //- Return the orientation at previous time-step
- inline const tensor& Q0() const;
-
- //- Return the velocity at previous time-step
- inline const vector& v0() const;
-
- //- Return the acceleration at previous time-step
- inline const vector& a0() const;
-
- //- Return the angular momentum at previous time-step
- inline const vector& pi0() const;
-
- //- Return the torque at previous time-step
- inline const vector& tau0() const;
-
// Edit
@@ -237,6 +221,10 @@ public:
sixDoFRigidBodyMotion(const sixDoFRigidBodyMotion&);
+ //- Destructor
+ ~sixDoFRigidBodyMotion();
+
+
// Member Functions
// Access
@@ -253,9 +241,6 @@ public:
//- Return the current centre of rotation
inline const point& centreOfRotation() const;
- //- Return the centre of rotation at previous time-step
- inline const point& centreOfRotation0() const;
-
//- Return the initial centre of mass
inline const point& initialCentreOfMass() const;
@@ -301,16 +286,9 @@ public:
// Update state
- //- First symplectic velocity adjust and motion part, required
- // before force calculation. Takes old timestep for variable
- // timestep cases.
- // Changes to Crank-Nicolson integration for subsequent iterations.
- void updatePosition(bool firstIter, scalar deltaT, scalar deltaT0);
-
- //- Second symplectic velocity adjust part
- // required after motion and force calculation.
+ //- Symplectic integration of velocities, orientation and position.
// Changes to Crank-Nicolson integration for subsequent iterations.
- void updateAcceleration
+ void update
(
bool firstIter,
const vector& fGlobal,
diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H
index d09e579e8cca4a3f5e51ed724bf8cd819ff8681c..ffe89819fd11e5028dcf1d7bd78c83b113ce15a8 100644
--- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H
+++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H
@@ -2,7 +2,7 @@
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
- \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation
+ \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
@@ -161,36 +161,6 @@ inline const Foam::vector& Foam::sixDoFRigidBodyMotion::tau() const
}
-inline const Foam::tensor& Foam::sixDoFRigidBodyMotion::Q0() const
-{
- return motionState0_.Q();
-}
-
-
-inline const Foam::vector& Foam::sixDoFRigidBodyMotion::v0() const
-{
- return motionState0_.v();
-}
-
-
-inline const Foam::vector& Foam::sixDoFRigidBodyMotion::a0() const
-{
- return motionState0_.a();
-}
-
-
-inline const Foam::vector& Foam::sixDoFRigidBodyMotion::pi0() const
-{
- return motionState0_.pi();
-}
-
-
-inline const Foam::vector& Foam::sixDoFRigidBodyMotion::tau0() const
-{
- return motionState0_.tau();
-}
-
-
inline Foam::point& Foam::sixDoFRigidBodyMotion::initialCentreOfRotation()
{
return initialCentreOfRotation_;
@@ -261,12 +231,6 @@ inline const Foam::point& Foam::sixDoFRigidBodyMotion::centreOfRotation() const
}
-inline const Foam::point& Foam::sixDoFRigidBodyMotion::centreOfRotation0() const
-{
- return motionState0_.centreOfRotation();
-}
-
-
inline const Foam::point&
Foam::sixDoFRigidBodyMotion::initialCentreOfMass() const
{
diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C
index 39453bb9d7353d95211e663119552b128b05a7ab..1d3d01856d0a0322c6badfa8bdda2a3e910d71c7 100644
--- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C
+++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C
@@ -205,9 +205,7 @@ void Foam::sixDoFRigidBodyMotionSolver::solve()
if (test_)
{
- motion_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value());
-
- motion_.updateAcceleration
+ motion_.update
(
firstIter,
ramp*(motion_.mass()*g.value()),
@@ -230,12 +228,7 @@ void Foam::sixDoFRigidBodyMotionSolver::solve()
f.calcForcesMoment();
- // Patch force data is valid for the current positions, so
- // calculate the forces on the motion object from this data, then
- // update the positions
- motion_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value());
-
- motion_.updateAcceleration
+ motion_.update
(
firstIter,
ramp*(f.forceEff() + motion_.mass()*g.value()),
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C
new file mode 100644
index 0000000000000000000000000000000000000000..366266e014b866dd63d9d4fae0a110e85bbb35f0
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C
@@ -0,0 +1,94 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CrankNicolson.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+ defineTypeNameAndDebug(CrankNicolson, 0);
+ addToRunTimeSelectionTable(sixDoFSolver, CrankNicolson, dictionary);
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::CrankNicolson::CrankNicolson
+(
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+)
+:
+ sixDoFSolver(body),
+ aoc_(dict.lookupOrDefault<scalar>("aoc", 0.5)),
+ voc_(dict.lookupOrDefault<scalar>("voc", 0.5))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::CrankNicolson::~CrankNicolson()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::sixDoFSolvers::CrankNicolson::solve
+(
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+)
+{
+ // Update the linear acceleration and torque
+ updateAcceleration(fGlobal, tauGlobal);
+
+ // Correct linear velocity
+ v() = tConstraints()
+ & (v0() + aDamp()*deltaT*(aoc_*a() + (1 - aoc_)*a0()));
+
+ // Correct angular momentum
+ pi() = rConstraints()
+ & (pi0() + aDamp()*deltaT*(aoc_*tau() + (1 - aoc_)*tau0()));
+
+ // Correct position
+ centreOfRotation() =
+ centreOfRotation0() + deltaT*(voc_*v() + (1 - voc_)*v0());
+
+ // Correct orientation
+ Tuple2<tensor, vector> Qpi =
+ rotate(Q0(), (voc_*pi() + (1 - voc_)*pi0()), deltaT);
+ Q() = Qpi.first();
+}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H
new file mode 100644
index 0000000000000000000000000000000000000000..a8dd2ccbff3d5c11a7553b1f7475de82fd3bada8
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H
@@ -0,0 +1,126 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::sixDoFSolvers::CrankNicolson
+
+Description
+ Crank-Nicolson 2nd-order time-integrator for 6DoF solid-body motion.
+
+ The off-centering coefficients for acceleration (velocity integration) and
+ velocity (position/orientation integration) may be specified but default
+ values of 0.5 for each are used if they are not specified. With the default
+ off-centering this scheme is equivalent to the Newmark scheme with default
+ coefficients.
+
+ Example specification in dynamicMeshDict:
+ \verbatim
+ solver
+ {
+ type CrankNicolson;
+ aoc 0.5; // Acceleration off-centering coefficient
+ voc 0.5; // Velocity off-centering coefficient
+ }
+ \endverbatim
+
+SeeAlso
+ Foam::sixDoFSolvers::Newmark
+
+SourceFiles
+ CrankNicolson.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef CrankNicolson_H
+#define CrankNicolson_H
+
+#include "sixDoFSolver.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+
+/*---------------------------------------------------------------------------*\
+ Class CrankNicolson Declaration
+\*---------------------------------------------------------------------------*/
+
+class CrankNicolson
+:
+ public sixDoFSolver
+{
+ // Private data
+
+ //- Acceleration off-centering coefficient (default: 0.5)
+ const scalar aoc_;
+
+ //- Velocity off-centering coefficient (default: 0.5)
+ const scalar voc_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("CrankNicolson");
+
+
+ // Constructors
+
+ //- Construct from a dictionary and the body
+ CrankNicolson
+ (
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+ );
+
+
+ //- Destructor
+ virtual ~CrankNicolson();
+
+
+ // Member Functions
+
+ //- Drag coefficient
+ virtual void solve
+ (
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sixDoFSolvers
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C
new file mode 100644
index 0000000000000000000000000000000000000000..4f39b3266c92c8a9ef0fe8c670dc74437d6d07fc
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C
@@ -0,0 +1,115 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "Newmark.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+ defineTypeNameAndDebug(Newmark, 0);
+ addToRunTimeSelectionTable(sixDoFSolver, Newmark, dictionary);
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::Newmark::Newmark
+(
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+)
+:
+ sixDoFSolver(body),
+ gamma_(dict.lookupOrDefault<scalar>("gamma", 0.5)),
+ beta_
+ (
+ max
+ (
+ 0.25*sqr(gamma_ + 0.5),
+ dict.lookupOrDefault<scalar>("beta", 0.25)
+ )
+ )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::Newmark::~Newmark()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::sixDoFSolvers::Newmark::solve
+(
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+)
+{
+ // Update the linear acceleration and torque
+ updateAcceleration(fGlobal, tauGlobal);
+
+ // Correct linear velocity
+ v() =
+ tConstraints()
+ & (v0() + aDamp()*deltaT*(gamma_*a() + (1 - gamma_)*a0()));
+
+ // Correct angular momentum
+ pi() =
+ rConstraints()
+ & (pi0() + aDamp()*deltaT*(gamma_*tau() + (1 - gamma_)*tau0()));
+
+ // Correct position
+ centreOfRotation() =
+ centreOfRotation0()
+ + (
+ tConstraints()
+ & (
+ deltaT*v0()
+ + sqr(deltaT)*beta_*a() + sqr(deltaT)*(0.5 - beta_)*a0()
+ )
+ );
+
+ // Correct orientation
+ vector piDeltaT =
+ rConstraints()
+ & (
+ deltaT*pi0()
+ + sqr(deltaT)*beta_*tau() + sqr(deltaT)*(0.5 - beta_)*tau0()
+ );
+ Tuple2<tensor, vector> Qpi = rotate(Q0(), piDeltaT, 1);
+ Q() = Qpi.first();
+}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H
new file mode 100644
index 0000000000000000000000000000000000000000..e605fcc294a3103afb4a4dafa633636ca622b577
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H
@@ -0,0 +1,124 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::sixDoFSolvers::Newmark
+
+Description
+ Newmark 2nd-order time-integrator for 6DoF solid-body motion.
+
+ Reference:
+ \verbatim
+ Newmark, N. M. (1959).
+ A method of computation for structural dynamics.
+ Journal of the Engineering Mechanics Division, 85(3), 67-94.
+ \endverbatim
+
+ Example specification in dynamicMeshDict:
+ \verbatim
+ solver
+ {
+ type Newmark;
+ gamma 0.5; // Velocity integration coefficient
+ beta 0.25; // Position integration coefficient
+ }
+ \endverbatim
+
+SourceFiles
+ Newmark.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Newmark_H
+#define Newmark_H
+
+#include "sixDoFSolver.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+
+/*---------------------------------------------------------------------------*\
+ Class Newmark Declaration
+\*---------------------------------------------------------------------------*/
+
+class Newmark
+:
+ public sixDoFSolver
+{
+ // Private data
+
+ //- Coefficient for velocity integration (default: 0.5)
+ const scalar gamma_;
+
+ //- Coefficient for position and orientation integration (default: 0.25)
+ const scalar beta_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("Newmark");
+
+
+ // Constructors
+
+ //- Construct from a dictionary and the body
+ Newmark
+ (
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+ );
+
+
+ //- Destructor
+ virtual ~Newmark();
+
+
+ // Member Functions
+
+ //- Drag coefficient
+ virtual void solve
+ (
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sixDoFSolvers
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C
new file mode 100644
index 0000000000000000000000000000000000000000..329a9aaea185bca039d59b8ea579e4da598f963d
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C
@@ -0,0 +1,57 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sixDoFSolver.H"
+
+// * * * * * * * * * * * * * * * * Selector * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::sixDoFSolver> Foam::sixDoFSolver::New
+(
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+)
+{
+ word sixDoFSolverType(dict.lookup("type"));
+
+ Info<< "Selecting sixDoFSolver " << sixDoFSolverType << endl;
+
+ dictionaryConstructorTable::iterator cstrIter =
+ dictionaryConstructorTablePtr_->find(sixDoFSolverType);
+
+ if (cstrIter == dictionaryConstructorTablePtr_->end())
+ {
+ FatalErrorIn("sixDoFSolver::New")
+ << "Unknown sixDoFSolverType type "
+ << sixDoFSolverType << endl << endl
+ << "Valid sixDoFSolver types are : " << endl
+ << dictionaryConstructorTablePtr_->sortedToc()
+ << exit(FatalError);
+ }
+
+ return cstrIter()(dict, body);
+}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C
new file mode 100644
index 0000000000000000000000000000000000000000..48fbba16f688070a2f40eb3f29dbbffd9a30d00f
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C
@@ -0,0 +1,51 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "sixDoFSolver.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(sixDoFSolver, 0);
+ defineRunTimeSelectionTable(sixDoFSolver, dictionary);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolver::sixDoFSolver(sixDoFRigidBodyMotion& body)
+:
+ body_(body)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolver::~sixDoFSolver()
+{}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H
new file mode 100644
index 0000000000000000000000000000000000000000..cb00660314a1be4042070224c4cc0aed1f75cfc1
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H
@@ -0,0 +1,190 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::sixDoFSolver
+
+Description
+
+SourceFiles
+ sixDoFSolver.C
+ newSixDoFSolver.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef sixDoFSolver_H
+#define sixDoFSolver_H
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "sixDoFRigidBodyMotion.H"
+#include "runTimeSelectionTables.H"
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class sixDoFSolver Declaration
+\*---------------------------------------------------------------------------*/
+
+class sixDoFSolver
+{
+protected:
+
+ // Protected data
+
+ //- The rigid body
+ sixDoFRigidBodyMotion& body_;
+
+
+ // Protected member functions
+
+ //- Return the current centre of rotation
+ inline point& centreOfRotation();
+
+ //- Return the orientation
+ inline tensor& Q();
+
+ //- Return non-const access to vector
+ inline vector& v();
+
+ //- Return non-const access to acceleration
+ inline vector& a();
+
+ //- Return non-const access to angular momentum
+ inline vector& pi();
+
+ //- Return non-const access to torque
+ inline vector& tau();
+
+ //- Return the centre of rotation at previous time-step
+ inline const point& centreOfRotation0() const;
+
+ //- Return the orientation at previous time-step
+ inline const tensor& Q0() const;
+
+ //- Return the velocity at previous time-step
+ inline const vector& v0() const;
+
+ //- Return the acceleration at previous time-step
+ inline const vector& a0() const;
+
+ //- Return the angular momentum at previous time-step
+ inline const vector& pi0() const;
+
+ //- Return the torque at previous time-step
+ inline const vector& tau0() const;
+
+ //- Acceleration damping coefficient (for steady-state simulations)
+ inline scalar aDamp() const;
+
+ //- Translational constraint tensor
+ inline tensor tConstraints() const;
+
+ //- Rotational constraint tensor
+ inline tensor rConstraints() const;
+
+ //- Apply rotation tensors to Q0 for the given torque (pi) and deltaT
+ // and return the rotated Q and pi as a tuple
+ inline Tuple2<tensor, vector> rotate
+ (
+ const tensor& Q0,
+ const vector& pi,
+ const scalar deltaT
+ ) const;
+
+ //- Update and relax accelerations from the force and torque
+ inline void updateAcceleration
+ (
+ const vector& fGlobal,
+ const vector& tauGlobal
+ );
+
+
+public:
+
+ //- Runtime type information
+ TypeName("sixDoFSolver");
+
+
+ // Declare runtime construction
+
+ declareRunTimeSelectionTable
+ (
+ autoPtr,
+ sixDoFSolver,
+ dictionary,
+ (
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+ ),
+ (dict, body)
+ );
+
+
+ // Constructors
+
+ // Construct for given body
+ sixDoFSolver(sixDoFRigidBodyMotion& body);
+
+
+ //- Destructor
+ virtual ~sixDoFSolver();
+
+
+ // Selectors
+
+ static autoPtr<sixDoFSolver> New
+ (
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+ );
+
+
+ // Member Functions
+
+ //- Drag coefficient
+ virtual void solve
+ (
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+ ) = 0;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "sixDoFSolverI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H
new file mode 100644
index 0000000000000000000000000000000000000000..0301b2b3b2bbb1a66df188122ffe7e4f583c2bdb
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H
@@ -0,0 +1,131 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
+
+inline Foam::point& Foam::sixDoFSolver::centreOfRotation()
+{
+ return body_.motionState_.centreOfRotation();
+}
+
+inline Foam::tensor& Foam::sixDoFSolver::Q()
+{
+ return body_.motionState_.Q();
+}
+
+inline Foam::vector& Foam::sixDoFSolver::v()
+{
+ return body_.motionState_.v();
+}
+
+inline Foam::vector& Foam::sixDoFSolver::a()
+{
+ return body_.motionState_.a();
+}
+
+inline Foam::vector& Foam::sixDoFSolver::pi()
+{
+ return body_.motionState_.pi();
+}
+
+inline Foam::vector& Foam::sixDoFSolver::tau()
+{
+ return body_.motionState_.tau();
+}
+
+
+inline const Foam::point& Foam::sixDoFSolver::centreOfRotation0() const
+{
+ return body_.motionState0_.centreOfRotation();
+}
+
+inline const Foam::tensor& Foam::sixDoFSolver::Q0() const
+{
+ return body_.motionState0_.Q();
+}
+
+
+inline const Foam::vector& Foam::sixDoFSolver::v0() const
+{
+ return body_.motionState0_.v();
+}
+
+
+inline const Foam::vector& Foam::sixDoFSolver::a0() const
+{
+ return body_.motionState0_.a();
+}
+
+
+inline const Foam::vector& Foam::sixDoFSolver::pi0() const
+{
+ return body_.motionState0_.pi();
+}
+
+
+inline const Foam::vector& Foam::sixDoFSolver::tau0() const
+{
+ return body_.motionState0_.tau();
+}
+
+inline Foam::scalar Foam::sixDoFSolver::aDamp() const
+{
+ return body_.aDamp_;
+}
+
+inline Foam::tensor Foam::sixDoFSolver::tConstraints() const
+{
+ return body_.tConstraints_;
+}
+
+inline Foam::tensor Foam::sixDoFSolver::rConstraints() const
+{
+ return body_.rConstraints_;
+}
+
+//- Apply rotation tensors to Q0 for the given torque (pi) and deltaT
+// and return the rotated Q and pi as a tuple
+inline Foam::Tuple2<Foam::tensor, Foam::vector> Foam::sixDoFSolver::rotate
+(
+ const tensor& Q0,
+ const vector& pi,
+ const scalar deltaT
+) const
+{
+ return body_.rotate(Q0, pi, deltaT);
+}
+
+//- Update and relax accelerations from the force and torque
+inline void Foam::sixDoFSolver::updateAcceleration
+(
+ const vector& fGlobal,
+ const vector& tauGlobal
+)
+{
+ body_.updateAcceleration(fGlobal, tauGlobal);
+}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C
new file mode 100644
index 0000000000000000000000000000000000000000..61d94abc7a976d3af278b91eea2e051c6ece6ceb
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C
@@ -0,0 +1,102 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "symplectic.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+ defineTypeNameAndDebug(symplectic, 0);
+ addToRunTimeSelectionTable(sixDoFSolver, symplectic, dictionary);
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::symplectic::symplectic
+(
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+)
+:
+ sixDoFSolver(body)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::sixDoFSolvers::symplectic::~symplectic()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::sixDoFSolvers::symplectic::solve
+(
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+)
+{
+ if (!firstIter)
+ {
+ FatalErrorIn("sixDoFSolvers::symplectic::solve")
+ << "The symplectic integrator is explicit "
+ "and can only be solved once per time-step"
+ << exit(FatalError);
+ }
+
+ // First simplectic step:
+ // Half-step for linear and angular velocities
+ // Update position and orientation
+
+ v() = tConstraints() & (v0() + aDamp()*0.5*deltaT0*a0());
+ pi() = rConstraints() & (pi0() + aDamp()*0.5*deltaT0*tau0());
+
+ centreOfRotation() = centreOfRotation0() + deltaT*v();
+
+ Tuple2<tensor, vector> Qpi = rotate(Q0(), pi(), deltaT);
+ Q() = Qpi.first();
+ pi() = rConstraints() & Qpi.second();
+
+ // Update the linear acceleration and torque
+ updateAcceleration(fGlobal, tauGlobal);
+
+ // Second simplectic step:
+ // Complete update of linear and angular velocities
+
+ v() += tConstraints() & aDamp()*0.5*deltaT*a();
+ pi() += rConstraints() & aDamp()*0.5*deltaT*tau();
+}
+
+
+// ************************************************************************* //
diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H
new file mode 100644
index 0000000000000000000000000000000000000000..d03b77bbd59593ea9e914fa828caba5dc2a40a30
--- /dev/null
+++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H
@@ -0,0 +1,123 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
+ \\/ 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 3 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, see <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::sixDoFSolvers::symplectic
+
+Description
+ Symplectic 2nd-order explicit time-integrator for 6DoF solid-body motion.
+
+ Reference:
+ \verbatim
+ Dullweber, A., Leimkuhler, B., & McLachlan, R. (1997).
+ Symplectic splitting methods for rigid body molecular dynamics.
+ The Journal of chemical physics, 107(15), 5840-5851.
+ \endverbatim
+
+ Can only be used for explicit integration of the motion of the body,
+ i.e. may only be called once per time-step, no outer-correctors may be
+ applied. For implicit integration with outer-correctors choose either
+ CrankNicolson or Newmark schemes.
+
+ Example specification in dynamicMeshDict:
+ \verbatim
+ solver
+ {
+ type symplectic;
+ }
+ \endverbatim
+
+SeeAlso
+ Foam::sixDoFSolvers::CrankNicolson
+ Foam::sixDoFSolvers::Newmark
+
+SourceFiles
+ symplectic.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef symplectic_H
+#define symplectic_H
+
+#include "sixDoFSolver.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace sixDoFSolvers
+{
+
+/*---------------------------------------------------------------------------*\
+ Class symplectic Declaration
+\*---------------------------------------------------------------------------*/
+
+class symplectic
+:
+ public sixDoFSolver
+{
+
+public:
+
+ //- Runtime type information
+ TypeName("symplectic");
+
+
+ // Constructors
+
+ //- Construct from a dictionary and the body
+ symplectic
+ (
+ const dictionary& dict,
+ sixDoFRigidBodyMotion& body
+ );
+
+
+ //- Destructor
+ virtual ~symplectic();
+
+
+ // Member Functions
+
+ //- Drag coefficient
+ virtual void solve
+ (
+ bool firstIter,
+ const vector& fGlobal,
+ const vector& tauGlobal,
+ scalar deltaT,
+ scalar deltaT0
+ );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace sixDoFSolvers
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict b/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict
index c2720340129d77645e823c7bf9773cadb2202671..13d0937968571b27eab91dacefac70fe0672c97a 100644
--- a/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict
+++ b/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict
@@ -41,6 +41,11 @@ sixDoFRigidBodyMotionCoeffs
rhoInf 1;
report on;
+ solver
+ {
+ type symplectic;
+ }
+
constraints
{
yLine
diff --git a/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict b/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict
index 7e7f4969af5564d68e522fad2b52d023ff5db4f9..02088e1075be36c5154c92b08985f908b1f97488 100644
--- a/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict
+++ b/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict
@@ -31,9 +31,16 @@ sixDoFRigidBodyMotionCoeffs
momentOfInertia (40 921 921);
rhoInf 1;
report on;
- accelerationRelaxation 0.3;
+
value uniform (0 0 0);
+ accelerationRelaxation 0.4;
+
+ solver
+ {
+ type Newmark;
+ }
+
constraints
{
zAxis
diff --git a/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution
index afd4c5c21931650b8f612585f08af661207917b6..4aaf52f6d20284b6028d368341a38367bb852767 100644
--- a/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution
+++ b/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution
@@ -60,7 +60,7 @@ solvers
cacheAgglomeration true;
tolerance 5e-8;
- relTol 0.001;
+ relTol 0;
};
p_rghFinal
@@ -86,8 +86,8 @@ PIMPLE
{
momentumPredictor no;
- nOuterCorrectors 1;
- nCorrectors 3;
+ nOuterCorrectors 3;
+ nCorrectors 1;
nNonOrthogonalCorrectors 0;
correctPhi yes;
diff --git a/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict b/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict
index 68f291d7762b8f7e65b310b8999f4530cfdad3fa..090aa63d16be1aaffe93f0b4c86e12f9a5a857da 100644
--- a/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict
+++ b/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict
@@ -59,7 +59,12 @@ sixDoFRigidBodyMotionCoeffs
};
report on;
- accelerationRelaxation 0.3;
+ accelerationRelaxation 0.7;
+
+ solver
+ {
+ type Newmark;
+ }
constraints
{
diff --git a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
index bc278a91fa4a94e324454178c182a4a28e7bc0db..ae6f2be175d2ac5f2f149077cd96c80ac550582d 100644
--- a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
+++ b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
@@ -117,7 +117,7 @@ solvers
PIMPLE
{
momentumPredictor no;
- nOuterCorrectors 5;
+ nOuterCorrectors 3;
nCorrectors 1;
nNonOrthogonalCorrectors 0;
correctPhi yes;
diff --git a/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict b/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict
deleted file mode 100644
index 7e7f4969af5564d68e522fad2b52d023ff5db4f9..0000000000000000000000000000000000000000
--- a/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict
+++ /dev/null
@@ -1,67 +0,0 @@
-/*--------------------------------*- C++ -*----------------------------------*\
-| ========= | |
-| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
-| \\ / O peration | Version: dev |
-| \\ / A nd | Web: www.OpenFOAM.org |
-| \\/ M anipulation | |
-\*---------------------------------------------------------------------------*/
-FoamFile
-{
- version 2.0;
- format ascii;
- class dictionary;
- object dynamicMeshDict;
-}
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-dynamicFvMesh dynamicMotionSolverFvMesh;
-
-motionSolverLibs ("libsixDoFRigidBodyMotion.so");
-
-solver sixDoFRigidBodyMotion;
-
-sixDoFRigidBodyMotionCoeffs
-{
- patches (hull);
- innerDistance 0.3;
- outerDistance 1;
-
- centreOfMass (2.929541 0 0.2);
- mass 412.73;
- momentOfInertia (40 921 921);
- rhoInf 1;
- report on;
- accelerationRelaxation 0.3;
- value uniform (0 0 0);
-
- constraints
- {
- zAxis
- {
- sixDoFRigidBodyMotionConstraint line;
- direction (0 0 1);
- }
- yPlane
- {
- sixDoFRigidBodyMotionConstraint axis;
- axis (0 1 0);
- }
- }
-
- restraints
- {
- translationDamper
- {
- sixDoFRigidBodyMotionRestraint linearDamper;
- coeff 8596;
- }
- rotationDamper
- {
- sixDoFRigidBodyMotionRestraint sphericalAngularDamper;
- coeff 11586;
- }
- }
-}
-
-
-// ************************************************************************* //