From e60f9f457f4e42a3a499e86644610a04658ca460 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Mon, 18 Apr 2016 15:40:23 +0100 Subject: [PATCH] src/rigidBodyDynamics/rigidBodyMotion: Added support for acceleration relaxation --- src/rigidBodyDynamics/Make/files | 1 + .../rigidBodyMotion/rigidBodyMotion.C | 13 +++++++++++++ .../rigidBodyMotion/rigidBodyMotion.H | 10 ++++++++++ .../rigidBodySolvers/CrankNicolson/CrankNicolson.C | 2 +- .../rigidBodySolvers/Newmark/Newmark.C | 4 ++-- .../rigidBodySolver/rigidBodySolver.H | 4 ---- .../rigidBodySolver/rigidBodySolverI.H | 6 ------ .../rigidBodySolvers/symplectic/symplectic.C | 4 ++-- src/rigidBodyMeshMotion/rigidBodyMeshMotion.C | 7 ++----- .../sixDoFSolvers/Newmark/Newmark.C | 4 ++-- 10 files changed, 33 insertions(+), 22 deletions(-) diff --git a/src/rigidBodyDynamics/Make/files b/src/rigidBodyDynamics/Make/files index 6352eabd93..07612c1038 100644 --- a/src/rigidBodyDynamics/Make/files +++ b/src/rigidBodyDynamics/Make/files @@ -4,6 +4,7 @@ bodies/jointBody/jointBody.C bodies/compositeBody/compositeBody.C bodies/subBody/subBody.C bodies/sphere/sphere.C +bodies/cuboid/cuboid.C joints/joint/joint.C joints/null/nullJoint.C diff --git a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C index 30ad316f6c..23146ef869 100644 --- a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C +++ b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C @@ -125,6 +125,19 @@ Foam::spatialTransform Foam::RBD::rigidBodyMotion::X00 } +void Foam::RBD::rigidBodyMotion::forwardDynamics +( + rigidBodyModelState& state, + const scalarField& tau, + const Field& fx +) const +{ + scalarField qDdotPrev = state.qDdot(); + rigidBodyModel::forwardDynamics(state, tau, fx); + state.qDdot() = aDamp_*(aRelax_*state.qDdot() + (1 - aRelax_)*qDdotPrev); +} + + void Foam::RBD::rigidBodyMotion::solve ( scalar deltaT, diff --git a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H index 4e3203de59..a88763e15d 100644 --- a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H +++ b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H @@ -158,6 +158,16 @@ public: // Update state + //- Calculate and optionally relax the joint acceleration qDdot from + // the joint state q, velocity qDot, internal force tau (in the + // joint frame) and external force fx (in the global frame) + void forwardDynamics + ( + rigidBodyModelState& state, + const scalarField& tau, + const Field& fx + ) const; + //- Integrate velocities, orientation and position // for the given time-step void solve diff --git a/src/rigidBodyDynamics/rigidBodySolvers/CrankNicolson/CrankNicolson.C b/src/rigidBodyDynamics/rigidBodySolvers/CrankNicolson/CrankNicolson.C index 6d48971be6..fc15b3e6f7 100644 --- a/src/rigidBodyDynamics/rigidBodySolvers/CrankNicolson/CrankNicolson.C +++ b/src/rigidBodyDynamics/rigidBodySolvers/CrankNicolson/CrankNicolson.C @@ -77,7 +77,7 @@ void Foam::RBD::rigidBodySolvers::CrankNicolson::solve model_.forwardDynamics(state(), tau, rfx); // Correct velocity - qDot() = qDot0() + aDamp()*deltaT()*(aoc_*qDdot() + (1 - aoc_)*qDdot0()); + qDot() = qDot0() + deltaT()*(aoc_*qDdot() + (1 - aoc_)*qDdot0()); // Correct position q() = q0() + deltaT()*(voc_*qDot() + (1 - voc_)*qDot0()); diff --git a/src/rigidBodyDynamics/rigidBodySolvers/Newmark/Newmark.C b/src/rigidBodyDynamics/rigidBodySolvers/Newmark/Newmark.C index 2d5e5d90ea..11f80a3732 100644 --- a/src/rigidBodyDynamics/rigidBodySolvers/Newmark/Newmark.C +++ b/src/rigidBodyDynamics/rigidBodySolvers/Newmark/Newmark.C @@ -85,12 +85,12 @@ void Foam::RBD::rigidBodySolvers::Newmark::solve // Correct velocity qDot() = qDot0() - + aDamp()*deltaT()*(gamma_*qDdot() + (1 - gamma_)*qDdot0()); + + deltaT()*(gamma_*qDdot() + (1 - gamma_)*qDdot0()); // Correct position q() = q0() + deltaT()*qDot0() - + sqr(deltaT())*beta_*qDdot() + sqr(deltaT())*(0.5 - beta_)*qDdot0(); + + sqr(deltaT())*(beta_*qDdot() + (0.5 - beta_)*qDdot0()); correctQuaternionJoints(); diff --git a/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolver.H b/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolver.H index 1ab119bcee..6c2c3c678b 100644 --- a/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolver.H +++ b/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolver.H @@ -94,10 +94,6 @@ protected: //- Return the previous time-step inline scalar deltaT0() const; - - //- Acceleration damping coefficient (for steady-state simulations) - inline scalar aDamp() const; - //- Correct the quaternion joints based on the current change in q void correctQuaternionJoints(); diff --git a/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolverI.H b/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolverI.H index 53736a02ab..afac337637 100644 --- a/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolverI.H +++ b/src/rigidBodyDynamics/rigidBodySolvers/rigidBodySolver/rigidBodySolverI.H @@ -85,10 +85,4 @@ inline Foam::scalar Foam::RBD::rigidBodySolver::deltaT0() const } -inline Foam::scalar Foam::RBD::rigidBodySolver::aDamp() const -{ - return model_.aDamp_; -} - - // ************************************************************************* // diff --git a/src/rigidBodyDynamics/rigidBodySolvers/symplectic/symplectic.C b/src/rigidBodyDynamics/rigidBodySolvers/symplectic/symplectic.C index 18b78ad89f..f0dbb3e0df 100644 --- a/src/rigidBodyDynamics/rigidBodySolvers/symplectic/symplectic.C +++ b/src/rigidBodyDynamics/rigidBodySolvers/symplectic/symplectic.C @@ -70,7 +70,7 @@ void Foam::RBD::rigidBodySolvers::symplectic::solve // First simplectic step: // Half-step for linear and angular velocities // Update position and orientation - qDot() = qDot0() + aDamp()*0.5*deltaT0()*qDdot(); + qDot() = qDot0() + 0.5*deltaT0()*qDdot(); q() = q0() + deltaT()*qDot(); correctQuaternionJoints(); @@ -88,7 +88,7 @@ void Foam::RBD::rigidBodySolvers::symplectic::solve // Second simplectic step: // Complete update of linear and angular velocities - qDot() += aDamp()*0.5*deltaT()*qDdot(); + qDot() += 0.5*deltaT()*qDdot(); } diff --git a/src/rigidBodyMeshMotion/rigidBodyMeshMotion.C b/src/rigidBodyMeshMotion/rigidBodyMeshMotion.C index f114cfc6b6..6397e637ed 100644 --- a/src/rigidBodyMeshMotion/rigidBodyMeshMotion.C +++ b/src/rigidBodyMeshMotion/rigidBodyMeshMotion.C @@ -250,19 +250,16 @@ void Foam::rigidBodyMeshMotion::solve() const label bodyID = bodyMeshes_[bi].bodyID_; dictionary forcesDict; - forcesDict.add("type", forces::typeName); forcesDict.add("patches", bodyMeshes_[bi].patches_); forcesDict.add("rhoInf", rhoInf_); forcesDict.add("rhoName", rhoName_); - forcesDict.add("CofR", model_.X0(bodyID).r()); + forcesDict.add("CofR", vector::zero); forces f("forces", db(), forcesDict); - f.calcForcesMoment(); - fx[bodyID].l() = f.forceEff(); - fx[bodyID].w() = f.momentEff(); + fx[bodyID] = spatialVector(f.momentEff(), f.forceEff()); } model_.solve diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C index 4f39b3266c..b1685e1528 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C @@ -96,7 +96,7 @@ void Foam::sixDoFSolvers::Newmark::solve tConstraints() & ( deltaT*v0() - + sqr(deltaT)*beta_*a() + sqr(deltaT)*(0.5 - beta_)*a0() + + aDamp()*sqr(deltaT)*(beta_*a() + (0.5 - beta_)*a0()) ) ); @@ -105,7 +105,7 @@ void Foam::sixDoFSolvers::Newmark::solve rConstraints() & ( deltaT*pi0() - + sqr(deltaT)*beta_*tau() + sqr(deltaT)*(0.5 - beta_)*tau0() + + aDamp()*sqr(deltaT)*(beta_*tau() + (0.5 - beta_)*tau0()) ); Tuple2 Qpi = rotate(Q0(), piDeltaT, 1); Q() = Qpi.first();