diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C index c001027d46..2f52c1b744 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C @@ -203,7 +203,7 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() // calculate the forces on the motion object from this data, then // update the positions - motion_.updatePosition(t.deltaTValue()); + motion_.updatePosition(t.deltaTValue(), t.deltaT0Value()); dictionary forcesDict; diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C index 1ab6cada87..1debd0b6a9 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C @@ -174,6 +174,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion() initialQ_(I), momentOfInertia_(diagTensor::one*VSMALL), mass_(VSMALL), + cDamp_(0.0), + aLim_(VGREAT), report_(false) {} @@ -190,6 +192,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion const point& initialCentreOfMass, const tensor& initialQ, const diagTensor& momentOfInertia, + scalar cDamp, + scalar aLim, bool report ) : @@ -211,6 +215,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion initialQ_(initialQ), momentOfInertia_(momentOfInertia), mass_(mass), + cDamp_(cDamp), + aLim_(aLim), report_(report) {} @@ -233,6 +239,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion(const dictionary& dict) ), momentOfInertia_(dict.lookup("momentOfInertia")), mass_(readScalar(dict.lookup("mass"))), + cDamp_(dict.lookupOrDefault("accelerationDampingCoeff", 0.0)), + aLim_(dict.lookupOrDefault("accelerationLimit", VGREAT)), report_(dict.lookupOrDefault("report", false)) { addRestraints(dict); @@ -246,17 +254,19 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion const sixDoFRigidBodyMotion& sDoFRBM ) : - motionState_(sDoFRBM.motionState()), - restraints_(sDoFRBM.restraints()), - restraintNames_(sDoFRBM.restraintNames()), - constraints_(sDoFRBM.constraints()), - constraintNames_(sDoFRBM.constraintNames()), - maxConstraintIterations_(sDoFRBM.maxConstraintIterations()), - initialCentreOfMass_(sDoFRBM.initialCentreOfMass()), - initialQ_(sDoFRBM.initialQ()), - momentOfInertia_(sDoFRBM.momentOfInertia()), - mass_(sDoFRBM.mass()), - report_(sDoFRBM.report()) + motionState_(sDoFRBM.motionState_), + restraints_(sDoFRBM.restraints_), + restraintNames_(sDoFRBM.restraintNames_), + constraints_(sDoFRBM.constraints_), + constraintNames_(sDoFRBM.constraintNames_), + maxConstraintIterations_(sDoFRBM.maxConstraintIterations_), + initialCentreOfMass_(sDoFRBM.initialCentreOfMass_), + initialQ_(sDoFRBM.initialQ_), + momentOfInertia_(sDoFRBM.momentOfInertia_), + mass_(sDoFRBM.mass_), + cDamp_(sDoFRBM.cDamp_), + aLim_(sDoFRBM.aLim_), + report_(sDoFRBM.report_) {} @@ -358,7 +368,8 @@ void Foam::sixDoFRigidBodyMotion::addConstraints void Foam::sixDoFRigidBodyMotion::updatePosition ( - scalar deltaT + scalar deltaT, + scalar deltaT0 ) { // First leapfrog velocity adjust and motion part, required before @@ -366,9 +377,32 @@ void Foam::sixDoFRigidBodyMotion::updatePosition if (Pstream::master()) { - v() += 0.5*deltaT*a(); + vector aClip = a(); + scalar aMag = mag(aClip); - pi() += 0.5*deltaT*tau(); + if (aMag > SMALL) + { + aClip /= aMag; + } + + if (aMag > aLim_) + { + WarningIn + ( + "void Foam::sixDoFRigidBodyMotion::updatePosition" + "(" + "scalar deltaT, " + "scalar deltaT0" + ")" + ) + << "Limited acceleration " << a() + << " to " << aClip*aLim_ + << endl; + } + + v() += 0.5*(1 - cDamp_)*deltaT0*aClip*min(aMag, aLim_); + + pi() += 0.5*(1 - cDamp_)*deltaT0*tau(); // Leapfrog move part centreOfMass() += deltaT*v(); @@ -401,9 +435,33 @@ void Foam::sixDoFRigidBodyMotion::updateForce applyConstraints(deltaT); - v() += 0.5*deltaT*a(); + vector aClip = a(); + scalar aMag = mag(aClip); - pi() += 0.5*deltaT*tau(); + if (aMag > SMALL) + { + aClip /= aMag; + } + + if (aMag > aLim_) + { + WarningIn + ( + "void Foam::sixDoFRigidBodyMotion::updateForce" + "(" + "const vector& fGlobal, " + "const vector& tauGlobal, " + "scalar deltaT" + ")" + ) + << "Limited acceleration " << a() + << " to " << aClip*aLim_ + << endl; + } + + v() += 0.5*(1 - cDamp_)*deltaT*aClip*min(aMag, aLim_); + + pi() += 0.5*(1 - cDamp_)*deltaT*tau(); if(report_) { diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H index a26c2378d1..bd0943fb1d 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H @@ -116,6 +116,15 @@ class sixDoFRigidBodyMotion //- Mass of the body scalar mass_; + //- Acceleration damping coefficient. Modify applied acceleration: + // v1 = v0 + a*dt - cDamp*a*dt + // = v0 + dt*f*(1 - cDamp)/m + // Increases effective mass by 1/(1 - cDamp). + scalar cDamp_; + + //- Acceleration magnitude limit - clips large accelerations + scalar aLim_; + //- Switch to turn reporting of motion data on and off Switch report_; @@ -235,6 +244,8 @@ public: const point& initialCentreOfMass, const tensor& initialQ, const diagTensor& momentOfInertia, + scalar cDamp = 0.0, + scalar aLim = VGREAT, bool report = false ); @@ -260,10 +271,12 @@ public: void addConstraints(const dictionary& dict); //- First leapfrog velocity adjust and motion part, required - // before force calculation + // before force calculation. Takes old timestep for variable + // timestep cases. void updatePosition ( - scalar deltaT + scalar deltaT, + scalar deltaT0 ); //- Second leapfrog velocity adjust part, required after motion and diff --git a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionIO.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionIO.C index 0312f396fe..bb3d009701 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionIO.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionIO.C @@ -40,6 +40,10 @@ void Foam::sixDoFRigidBodyMotion::write(Ostream& os) const << momentOfInertia_ << token::END_STATEMENT << nl; os.writeKeyword("mass") << mass_ << token::END_STATEMENT << nl; + os.writeKeyword("accelerationDampingCoeff") + << cDamp_ << token::END_STATEMENT << nl; + os.writeKeyword("accelerationLimit") + << aLim_ << token::END_STATEMENT << nl; os.writeKeyword("report") << report_ << token::END_STATEMENT << nl; 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 7a045c04d4..fb168c8813 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 @@ -83,15 +83,20 @@ Foam::sixDoFRigidBodyMotionRestraints::linearAxialAngularSpring::restrain ) const { vector refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 1, 0); + vector oldDir = refQ_ & refDir; + vector newDir = motion.orientation() & refDir; if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95) { // Directions getting close to the axis, change reference + refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 0, 1); - oldDir = refQ_ & refDir; - newDir = motion.orientation() & refDir; + + vector oldDir = refQ_ & 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/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionRestraint/tabulatedAxialAngularSpring/tabulatedAxialAngularSpring.C index 8581b3dd50..50cdb1460b 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 @@ -87,6 +87,7 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain vector refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 1, 0); vector oldDir = refQ_ & refDir; + vector newDir = motion.orientation() & refDir; if (mag(oldDir & axis_) > 0.95 || mag(newDir & axis_) > 0.95) @@ -95,8 +96,9 @@ Foam::sixDoFRigidBodyMotionRestraints::tabulatedAxialAngularSpring::restrain refDir = rotationTensor(vector(1, 0 ,0), axis_) & vector(0, 0, 1); - oldDir = refQ_ & refDir; - newDir = motion.orientation() & refDir; + vector oldDir = refQ_ & 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/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C index 75771c7d84..8927b8f4c5 100644 --- a/src/postProcessing/functionObjects/forces/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C +++ b/src/postProcessing/functionObjects/forces/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C @@ -147,7 +147,7 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() const polyMesh& mesh = this->dimensionedInternalField().mesh()(); const Time& t = mesh.time(); - motion_.updatePosition(t.deltaTValue()); + motion_.updatePosition(t.deltaTValue(), t.deltaT0Value()); vector gravity = vector::zero;