From cbc4cee9351dfd2a8366ea9ba5205166c736dd7c Mon Sep 17 00:00:00 2001 From: graham Date: Mon, 7 Mar 2011 10:52:08 +0000 Subject: [PATCH] ENH: Using proper expression for accelerating frame terms. No "axis" required - when considering omegaDot the concept becomes meaningless. Centre of rotation can be specified, to avoid it needing to be the origin, and to allow it to move. --- .../createNonInertialFrameFields.H | 50 ++--- .../NonInertialFrame/NonInertialFrameForce.C | 173 +++++++----------- .../NonInertialFrame/NonInertialFrameForce.H | 40 ++-- .../NonInertialFrame/NonInertialFrameForceI.H | 21 +-- .../Kinematic/ParticleForces/SRF/SRFForce.C | 3 +- 5 files changed, 97 insertions(+), 190 deletions(-) diff --git a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createNonInertialFrameFields.H b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createNonInertialFrameFields.H index 42304ddcb8..abaede9385 100644 --- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createNonInertialFrameFields.H +++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createNonInertialFrameFields.H @@ -22,24 +22,24 @@ } - IOobject angularVelocityMagHeader + IOobject angularVelocityHeader ( - "angularVelocityMag", + "angularVelocity", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE ); - autoPtr angularVelocityMagPtr; + autoPtr angularVelocityPtr; - if (angularVelocityMagHeader.headerOk()) + if (angularVelocityHeader.headerOk()) { - Info<< " Reading " << angularVelocityMagHeader.name() << endl; + Info<< " Reading " << angularVelocityHeader.name() << endl; - angularVelocityMagPtr.reset + angularVelocityPtr.reset ( - new uniformDimensionedScalarField(angularVelocityMagHeader) + new uniformDimensionedVectorField(angularVelocityHeader) ); } @@ -66,45 +66,23 @@ } - IOobject axisHeader + IOobject centreOfRotationHeader ( - "axis", + "centreOfRotation", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE ); - autoPtr axisPtr; + autoPtr centreOfRotationPtr; - if (axisHeader.headerOk()) + if (centreOfRotationHeader.headerOk()) { - Info<< " Reading " << axisHeader.name() << endl; + Info<< " Reading " << centreOfRotationHeader.name() << endl; - axisPtr.reset + centreOfRotationPtr.reset ( - new uniformDimensionedVectorField(axisHeader) - ); - } - - - IOobject axisRefPointHeader - ( - "axisRefPoint", - runTime.timeName(), - mesh, - IOobject::MUST_READ, - IOobject::AUTO_WRITE - ); - - autoPtr axisRefPointPtr; - - if (axisRefPointHeader.headerOk()) - { - Info<< " Reading " << axisRefPointHeader.name() << endl; - - axisRefPointPtr.reset - ( - new uniformDimensionedVectorField(axisRefPointHeader) + new uniformDimensionedVectorField(centreOfRotationHeader) ); } diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.C b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.C index 807a07aa16..77d8c3e037 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.C +++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.C @@ -47,15 +47,15 @@ Foam::NonInertialFrameForce::NonInertialFrameForce ) ), W_(vector::zero), - omegaMagName_ + omegaName_ ( this->coeffs().template lookupOrDefault ( - "angularVelocityMagName", - "angularVelocityMag" + "angularVelocityName", + "angularVelocity" ) ), - omegaMag_(0.0), + omega_(vector::zero), omegaDotName_ ( this->coeffs().template lookupOrDefault @@ -65,25 +65,15 @@ Foam::NonInertialFrameForce::NonInertialFrameForce ) ), omegaDot_(vector::zero), - axisName_ + centreOfRotationName_ ( this->coeffs().template lookupOrDefault ( - "axisName", - "axis" + "centreOfRotationName", + "centreOfRotation" ) ), - axis_(vector::zero), - hasAxis_(false), - axisRefPointName_ - ( - this->coeffs().template lookupOrDefault - ( - "axisRefPointName", - "axisRefPoint" - ) - ), - axisRefPoint_(vector::zero) + centreOfRotation_(vector::zero) {} @@ -96,15 +86,12 @@ Foam::NonInertialFrameForce::NonInertialFrameForce ParticleForce(niff), WName_(niff.WName_), W_(niff.W_), - omegaMagName_(niff.omegaMagName_), - omegaMag_(niff.omegaMag_), + omegaName_(niff.omegaName_), + omega_(niff.omega_), omegaDotName_(niff.omegaDotName_), omegaDot_(niff.omegaDot_), - axisName_(niff.axisName_), - axis_(niff.axis_), - hasAxis_(niff.hasAxis_), - axisRefPointName_(niff.axisRefPointName_), - axisRefPoint_(niff.axisRefPoint_) + centreOfRotationName_(niff.centreOfRotationName_), + centreOfRotation_(niff.centreOfRotation_) {} @@ -121,71 +108,18 @@ template void Foam::NonInertialFrameForce::cacheFields(const bool store) { W_ = vector::zero; - omegaMag_ = 0.0; + omega_ = vector::zero; omegaDot_ = vector::zero; - axis_ = vector::zero; - hasAxis_ = false; - axisRefPoint_ = vector::zero; + centreOfRotation_ = vector::zero; if (store) { if ( - this->mesh().template - foundObject(omegaMagName_) - ) - { - uniformDimensionedScalarField omegaMag = this->mesh().template - lookupObject(omegaMagName_); - - omegaMag_ = omegaMag.value(); - - // If omegaMag is found, require that the axis and axisRefPoint is - // found. - uniformDimensionedVectorField a = this->mesh().template - lookupObject(axisName_); - - axis_ = a.value(); - - hasAxis_ = true; - - scalar axisMag = mag(axis_); - - if (mag(axis_) < SMALL) - { - FatalErrorIn - ( - "void Foam::NonInertialFrameForce::" - "cacheFields(const bool store)" - ) << axisName_ << " " << axis_ << " too small." - << abort(FatalError); - } - - axis_ /= axisMag; - - uniformDimensionedVectorField axisRefPoint = this->mesh().template - lookupObject(axisRefPointName_); - - axisRefPoint_ = axisRefPoint.value(); - - // Only look for omegaDot is omegaMag is found, optional. - if + this->mesh().template foundObject ( - this->mesh().template - foundObject(omegaDotName_) + WName_ ) - { - uniformDimensionedVectorField omegaDot = this->mesh().template - lookupObject(omegaDotName_); - - omegaDot_ = omegaDot.value(); - } - } - - if - ( - this->mesh().template - foundObject(WName_) ) { uniformDimensionedVectorField W = this->mesh().template @@ -193,13 +127,50 @@ void Foam::NonInertialFrameForce::cacheFields(const bool store) W_ = W.value(); } - else if (!hasAxis_) - { - WarningIn + + if + ( + this->mesh().template foundObject ( - "void Foam::NonInertialFrameForce::" - "cacheFields(const bool store)" - ) << "No " << typeName << " variables found." << endl; + omegaName_ + ) + ) + { + uniformDimensionedVectorField omega = this->mesh().template + lookupObject(omegaName_); + + omega_ = omega.value(); + } + + if + ( + this->mesh().template foundObject + ( + omegaDotName_ + ) + ) + { + uniformDimensionedVectorField omegaDot = this->mesh().template + lookupObject(omegaDotName_); + + omegaDot_ = omegaDot.value(); + } + + if + ( + this->mesh().template foundObject + ( + centreOfRotationName_ + ) + ) + { + uniformDimensionedVectorField omegaDot = this->mesh().template + lookupObject + ( + centreOfRotationName_ + ); + + centreOfRotation_ = omegaDot.value(); } } } @@ -217,24 +188,16 @@ Foam::forceSuSp Foam::NonInertialFrameForce::calcNonCoupled { forceSuSp value(vector::zero, 0.0); - value.Su() += -mass*W_; + const vector& r = p.position() - centreOfRotation_; - if (hasAxis_) - { - const vector pRel = p.position() - axisRefPoint_; - - const vector r = pRel - axis_*(axis_ & pRel); - - vector omega = axis_*omegaMag_; - - value.Su() += - mass - *( - (r ^ omegaDot_) - + 2.0*(p.U() ^ omega) - + (omega ^ (r ^ omega)) - ); - } + value.Su() = + mass + *( + -W_ + + (r ^ omegaDot_) + + 2.0*(p.U() ^ omega_) + + (omega_ ^ (r ^ omega_)) + ); return value; } diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.H b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.H index 093be4520d..aaf8f1065a 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.H +++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForce.H @@ -62,36 +62,23 @@ class NonInertialFrameForce //- The linear acceleration of the reference frame vector W_; - //- Name of the angular velocity magnitude field - word omegaMagName_; + //- Name of the angular velocity field + word omegaName_; - //- The magnitude of the angular velocity of the reference frame, - // combines with axis to give omega - scalar omegaMag_; + //- The angular velocity of the reference frame + vector omega_; //- Name of the angular acceleration field word omegaDotName_; - //- Pointer to the angular acceleration of the reference frame + //- The angular acceleration of the reference frame vector omegaDot_; - //- Name of the axis field - word axisName_; + //- Name of the centre of rotation field + word centreOfRotationName_; - //- Pointer to the axis of rotation - assumed to be a unit vector. - // duplication of omega to allow situations where omega = 0 and - // omegaDot != 0. - vector axis_; - - // Boolean flag for whether rotational motion is active - bool hasAxis_; - - //- Name of the axisRefPoint field - word axisRefPointName_; - - //- Pointer to a point in non-inertial space on the axis of rotation - // (omega), used to calculate r. - vector axisRefPoint_; + //- The centre of rotation of the reference frame + vector centreOfRotation_; public: @@ -141,14 +128,9 @@ public: //- Return the angular acceleration of the reference frame inline const vector& omegaDot() const; - //- Return the axis of rotation - inline const vector& axis() const; + //- Return the centre of rotation of the reference frame + inline const vector& centreOfRotation() const; - //- Return bool stating if the frame is rotating - inline bool hasAxis() const; - - //- Return the point in non-inertial space on the axis of rotation - inline const vector& axisRefPoint() const; // Evaluation diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForceI.H b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForceI.H index ae4ffc13a7..fb59297222 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForceI.H +++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/NonInertialFrame/NonInertialFrameForceI.H @@ -35,7 +35,7 @@ inline const Foam::vector& Foam::NonInertialFrameForce::W() const template inline const Foam::vector& Foam::NonInertialFrameForce::omega() const { - return omegaMag_*axis_; + return omega_; } @@ -49,24 +49,9 @@ Foam::NonInertialFrameForce::omegaDot() const template inline const Foam::vector& -Foam::NonInertialFrameForce::axis() const +Foam::NonInertialFrameForce::centreOfRotation() const { - return axis_; -} - - -template -inline bool Foam::NonInertialFrameForce::hasAxis() const -{ - return hasAxis_; -} - - -template -inline const Foam::vector& -Foam::NonInertialFrameForce::axisRefPoint() const -{ - return axisRefPoint_; + return centreOfRotation_; } diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/SRF/SRFForce.C b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/SRF/SRFForce.C index e06789f3fb..e001d29f19 100644 --- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/SRF/SRFForce.C +++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/SRF/SRFForce.C @@ -92,9 +92,8 @@ Foam::forceSuSp Foam::SRFForce::calcNonCoupled const typename SRF::SRFModel& srf = *srfPtr_; const vector& omega = srf.omega().value(); - const vector& axis = srf.axis(); - const vector r = p.position() - axis*(axis & p.position()); + const vector& r = p.position(); // Coriolis and centrifugal acceleration terms value.Su() = mass*(2.0*(p.U() ^ omega) + (omega ^ (r ^ omega)));