From a950efe49ea2fe43158407289f5b0065c5f8452f Mon Sep 17 00:00:00 2001 From: graham Date: Tue, 9 Feb 2010 15:02:27 +0000 Subject: [PATCH] BUG: sixDoFRigidBodyMotion. Fixing problems that only occur when orientation != I at start. Being more careful with when using directions relative to the absolute axes or the direction of something relative to its initial direction. Renaming p0, v0 to pInitial, vInitial to be clearer. fixedOrientation constraint does not need a reference orientation, as it is not to move at all, relative to any orientation, so the Cartesian axes are fine. --- .../sixDoFRigidBodyMotion.C | 18 +++++++------ .../sixDoFRigidBodyMotion.H | 15 ++++++----- .../fixedAxis/fixedAxis.C | 2 +- .../fixedOrientation/fixedOrientation.C | 26 +------------------ .../fixedOrientation/fixedOrientation.H | 6 ----- .../sixDoFRigidBodyMotionI.H | 20 +++++++++----- .../linearAxialAngularSpring.C | 4 +-- .../sphericalAngularSpring.C | 2 +- .../tabulatedAxialAngularSpring.C | 4 +-- 9 files changed, 39 insertions(+), 58 deletions(-) diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C index d8d2baaddd..c6a738678c 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C @@ -139,11 +139,13 @@ void Foam::sixDoFRigidBodyMotion::applyConstraints(scalar deltaT) << ") exceeded." << nl << exit(FatalError); } - else if (report_) + + Info<< "sixDoFRigidBodyMotion constraints converged in " + << iteration << " iterations" << endl; + + if (report_) { - Info<< "sixDoFRigidBodyMotion constraints converged in " - << iteration << " iterations" - << nl << "Constraint force: " << cFA << nl + Info<< "Constraint force: " << cFA << nl << "Constraint moment: " << cMA << endl; } @@ -444,7 +446,7 @@ void Foam::sixDoFRigidBodyMotion::updateForce Foam::point Foam::sixDoFRigidBodyMotion::predictedPosition ( - const point& p0, + const point& pInitial, const vector& deltaForce, const vector& deltaMoment, scalar deltaT @@ -463,14 +465,14 @@ Foam::point Foam::sixDoFRigidBodyMotion::predictedPosition return ( centreOfMassTemp - + (QTemp & initialQ_.T() & (p0 - initialCentreOfMass_)) + + (QTemp & initialQ_.T() & (pInitial - initialCentreOfMass_)) ); } Foam::vector Foam::sixDoFRigidBodyMotion::predictedOrientation ( - const vector& v0, + const vector& vInitial, const vector& deltaMoment, scalar deltaT ) const @@ -481,7 +483,7 @@ Foam::vector Foam::sixDoFRigidBodyMotion::predictedOrientation rotate(QTemp, piTemp, deltaT); - return (QTemp & initialQ_.T() & v0); + return (QTemp & initialQ_.T() & vInitial); } diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H index aee340a81a..20cbdae929 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H @@ -287,27 +287,30 @@ public: //- Transform the given initial state pointField by the current // motion state - inline tmp currentPosition(const pointField& p0) const; + inline tmp currentPosition + ( + const pointField& pInitial + ) const; //- Transform the given initial state point by the current motion // state - inline point currentPosition(const point& p0) const; + inline point currentPosition(const point& pInitial) const; //- Transform the given initial state direction by the current // motion state - inline vector currentOrientation(const vector& dir0) const; + inline vector currentOrientation(const vector& vInitial) const; //- Access the orientation tensor, Q. // globalVector = Q & bodyLocalVector // bodyLocalVector = Q.T() & globalVector - inline const tensor& currentOrientation() const; + inline const tensor& orientation() const; //- Predict the position of the supplied initial state point // after deltaT given the current motion state and the // additional supplied force and moment point predictedPosition ( - const point& p0, + const point& pInitial, const vector& deltaForce, const vector& deltaMoment, scalar deltaT @@ -318,7 +321,7 @@ public: // additional supplied moment vector predictedOrientation ( - const vector& v0, + const vector& vInitial, const vector& deltaMoment, scalar deltaT ) const; diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedAxis/fixedAxis.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedAxis/fixedAxis.C index 699849929a..5d9d4800d9 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedAxis/fixedAxis.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedAxis/fixedAxis.C @@ -97,7 +97,7 @@ bool Foam::sixDoFRigidBodyMotionConstraints::fixedAxis::constrain { rotationAxis /= magRotationAxis; - const tensor& Q = motion.currentOrientation(); + const tensor& Q = motion.orientation(); // Transform rotationAxis to body local system rotationAxis = Q.T() & rotationAxis; diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.C index eb634772c4..14f8fdadaa 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.C @@ -52,8 +52,7 @@ Foam::sixDoFRigidBodyMotionConstraints::fixedOrientation::fixedOrientation const dictionary& sDoFRBMCDict ) : - sixDoFRigidBodyMotionConstraint(sDoFRBMCDict), - fixedOrientation_() + sixDoFRigidBodyMotionConstraint(sDoFRBMCDict) { read(sDoFRBMCDict); } @@ -99,10 +98,6 @@ bool Foam::sixDoFRigidBodyMotionConstraints::fixedOrientation::constrain deltaT ); - axis = (fixedOrientation_ & axis); - - refDir = (fixedOrientation_ & refDir); - // Removing any axis component from predictedDir predictedDir -= (axis & predictedDir)*axis; @@ -177,23 +172,6 @@ bool Foam::sixDoFRigidBodyMotionConstraints::fixedOrientation::read { sixDoFRigidBodyMotionConstraint::read(sDoFRBMCDict); - fixedOrientation_ = - sDoFRBMCCoeffs_.lookupOrDefault("fixedOrientation", I); - - if (mag(mag(fixedOrientation_) - sqrt(3.0)) > 1e-9) - { - FatalErrorIn - ( - "Foam::sixDoFRigidBodyMotionConstraints::fixedOrientation::read" - "(" - "const dictionary& sDoFRBMCDict" - ")" - ) - << "fixedOrientation " << fixedOrientation_ - << " is not a rotation tensor." - << exit(FatalError); - } - return true; } @@ -203,8 +181,6 @@ void Foam::sixDoFRigidBodyMotionConstraints::fixedOrientation::write Ostream& os ) const { - os.writeKeyword("fixedOrientation") - << fixedOrientation_ << token::END_STATEMENT << nl; } // ************************************************************************* // diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.H b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.H index 8975f60862..6aaa50aa6f 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.H +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionConstraint/fixedOrientation/fixedOrientation.H @@ -59,12 +59,6 @@ class fixedOrientation public sixDoFRigidBodyMotionConstraint { - // Private data - - //- Reference orientation where there is no moment - tensor fixedOrientation_; - - public: //- Runtime type information diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H index abf8a5d31c..2335dcf8eb 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H @@ -222,34 +222,40 @@ inline Foam::vector& Foam::sixDoFRigidBodyMotion::tau() // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // inline Foam::tmp -Foam::sixDoFRigidBodyMotion::currentPosition(const pointField& p0) const +Foam::sixDoFRigidBodyMotion::currentPosition(const pointField& pInitial) const { return - (centreOfMass() + (Q() & initialQ_.T() & (p0 - initialCentreOfMass_))); + ( + centreOfMass() + + (Q() & initialQ_.T() & (pInitial - initialCentreOfMass_))) + ; } inline Foam::point Foam::sixDoFRigidBodyMotion::currentPosition ( - const point& p0 + const point& pInitial ) const { return - (centreOfMass() + (Q() & initialQ_.T() & (p0 - initialCentreOfMass_))); + ( + centreOfMass() + + (Q() & initialQ_.T() & (pInitial - initialCentreOfMass_)) + ); } inline Foam::vector Foam::sixDoFRigidBodyMotion::currentOrientation ( - const vector& dir0 + const vector& vInitial ) const { - return (Q() & initialQ_.T() & dir0); + return (Q() & initialQ_.T() & vInitial); } inline const Foam::tensor& -Foam::sixDoFRigidBodyMotion::currentOrientation() const +Foam::sixDoFRigidBodyMotion::orientation() const { return Q(); } diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/linearAxialAngularSpring/linearAxialAngularSpring.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/linearAxialAngularSpring/linearAxialAngularSpring.C index 16c72280d6..67f7206e59 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/linearAxialAngularSpring/linearAxialAngularSpring.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/linearAxialAngularSpring/linearAxialAngularSpring.C @@ -86,7 +86,7 @@ Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::restrain vector oldDir = refQ_ & refDir; - vector newDir = motion.currentOrientation(refDir); + vector newDir = motion.orientation() & refDir; if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95) { @@ -96,7 +96,7 @@ Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::restrain vector oldDir = refQ_ & refDir; - vector newDir = motion.currentOrientation(refDir); + vector newDir = motion.orientation() & refDir; } // Removing any axis component from oldDir and newDir and normalising diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/sphericalAngularSpring/sphericalAngularSpring.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/sphericalAngularSpring/sphericalAngularSpring.C index c3f4292e7b..43a5e795b6 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/sphericalAngularSpring/sphericalAngularSpring.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/sphericalAngularSpring/sphericalAngularSpring.C @@ -92,7 +92,7 @@ Foam::sixDoFRigidBodyMotionRestraints::sphericalAngularSpring::restrain refDir[(cmpt + 1) % 3] = 1; - vector newDir = motion.currentOrientation(refDir); + vector newDir = motion.orientation() & refDir; axis = (refQ_ & axis); diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C index d5c4a50d66..1a173c2445 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C @@ -88,7 +88,7 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain vector oldDir = refQ_ & refDir; - vector newDir = motion.currentOrientation(refDir); + vector newDir = motion.orientation() & refDir; if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95) { @@ -98,7 +98,7 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain vector oldDir = refQ_ & refDir; - vector newDir = motion.currentOrientation(refDir); + vector newDir = motion.orientation() & refDir; } // Removing any axis component from oldDir and newDir and normalising