diff --git a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C index f14acc6e60..94794dfe50 100644 --- a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C +++ b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.C @@ -178,6 +178,28 @@ void Foam::RBD::rigidBodyMotion::status(const label bodyID) const } +Foam::tmp Foam::RBD::rigidBodyMotion::transformPoints +( + const label bodyID, + const pointField& initialPoints +) const +{ + // Calculate the transform from the initial state in the global frame + // to the current state in the global frame + spatialTransform X(X0(bodyID).inv() & X00(bodyID)); + + tmp tpoints(new pointField(initialPoints.size())); + pointField& points = tpoints.ref(); + + forAll(points, i) + { + points[i] = X.transformPoint(initialPoints[i]); + } + + return tpoints; +} + + Foam::tmp Foam::RBD::rigidBodyMotion::transformPoints ( const label bodyID, diff --git a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H index c386169064..d065d9797c 100644 --- a/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H +++ b/src/rigidBodyDynamics/rigidBodyMotion/rigidBodyMotion.H @@ -183,6 +183,14 @@ public: // Transformations + //- Transform the given initial pointField of the specified body + // to correspond to the current motion state + tmp transformPoints + ( + const label bodyID, + const pointField& initialPoints + ) const; + //- Transform the given initial pointField of the specified body // to correspond to the current motion state scaled using // 'slerp' interpolation diff --git a/src/rigidBodyMeshMotion/Make/files b/src/rigidBodyMeshMotion/Make/files index f7fcbc1361..57443ab01d 100644 --- a/src/rigidBodyMeshMotion/Make/files +++ b/src/rigidBodyMeshMotion/Make/files @@ -1,3 +1,4 @@ -rigidBodyMeshMotion.C +rigidBodyMeshMotion/rigidBodyMeshMotion.C +rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.C LIB = $(FOAM_LIBBIN)/librigidBodyMeshMotion diff --git a/src/rigidBodyMeshMotion/rigidBodyMeshMotion.C b/src/rigidBodyMeshMotion/rigidBodyMeshMotion/rigidBodyMeshMotion.C similarity index 100% rename from src/rigidBodyMeshMotion/rigidBodyMeshMotion.C rename to src/rigidBodyMeshMotion/rigidBodyMeshMotion/rigidBodyMeshMotion.C diff --git a/src/rigidBodyMeshMotion/rigidBodyMeshMotion.H b/src/rigidBodyMeshMotion/rigidBodyMeshMotion/rigidBodyMeshMotion.H similarity index 100% rename from src/rigidBodyMeshMotion/rigidBodyMeshMotion.H rename to src/rigidBodyMeshMotion/rigidBodyMeshMotion/rigidBodyMeshMotion.H diff --git a/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.C b/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.C new file mode 100644 index 0000000000..f6476e46a0 --- /dev/null +++ b/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.C @@ -0,0 +1,334 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +\*---------------------------------------------------------------------------*/ + +#include "rigidBodyMeshMotionSolver.H" +#include "addToRunTimeSelectionTable.H" +#include "polyMesh.H" +#include "pointPatchDist.H" +#include "pointConstraints.H" +#include "uniformDimensionedFields.H" +#include "forces.H" +#include "mathematicalConstants.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(rigidBodyMeshMotionSolver, 0); + + addToRunTimeSelectionTable + ( + motionSolver, + rigidBodyMeshMotionSolver, + dictionary + ); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::rigidBodyMeshMotionSolver::bodyMesh::bodyMesh +( + const polyMesh& mesh, + const word& name, + const label bodyID, + const dictionary& dict +) +: + name_(name), + bodyID_(bodyID), + patches_(wordReList(dict.lookup("patches"))), + patchSet_(mesh.boundaryMesh().patchSet(patches_)) +{} + + +Foam::rigidBodyMeshMotionSolver::rigidBodyMeshMotionSolver +( + const polyMesh& mesh, + const IOdictionary& dict +) +: + motionSolver(mesh, dict, typeName), + model_ + ( + coeffDict(), + IOobject + ( + "rigidBodyMotionState", + mesh.time().timeName(), + "uniform", + mesh + ).headerOk() + ? IOdictionary + ( + IOobject + ( + "rigidBodyMotionState", + mesh.time().timeName(), + "uniform", + mesh, + IOobject::READ_IF_PRESENT, + IOobject::NO_WRITE, + false + ) + ) + : coeffDict() + ), + test_(coeffDict().lookupOrDefault("test", false)), + rhoInf_(1.0), + rhoName_(coeffDict().lookupOrDefault("rho", "rho")), + curTimeIndex_(-1), + meshSolverPtr_ + ( + motionSolver::New + ( + mesh, + IOdictionary + ( + IOobject + ( + "rigidBodyMotionSolver:meshSolver", + mesh.time().constant(), + mesh + ), + coeffDict().subDict("meshSolver") + ) + ) + ), + meshSolver_(refCast(meshSolverPtr_())) +{ + if (rhoName_ == "rhoInf") + { + rhoInf_ = readScalar(coeffDict().lookup("rhoInf")); + } + + const dictionary& bodiesDict = coeffDict().subDict("bodies"); + + forAllConstIter(IDLList, bodiesDict, iter) + { + const dictionary& bodyDict = iter().dict(); + + if (bodyDict.found("patches")) + { + const label bodyID = model_.bodyID(iter().keyword()); + + if (bodyID == -1) + { + FatalErrorInFunction + << "Body " << iter().keyword() + << " has been merged with another body" + " and cannot be assigned a set of patches" + << exit(FatalError); + } + + bodyMeshes_.append + ( + new bodyMesh + ( + mesh, + iter().keyword(), + bodyID, + bodyDict + ) + ); + } + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::rigidBodyMeshMotionSolver::~rigidBodyMeshMotionSolver() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +Foam::tmp +Foam::rigidBodyMeshMotionSolver::curPoints() const +{ + return meshSolverPtr_->curPoints(); +} + + +void Foam::rigidBodyMeshMotionSolver::solve() +{ + const Time& t = mesh().time(); + + if (mesh().nPoints() != meshSolver_.points0().size()) + { + FatalErrorInFunction + << "The number of points in the mesh seems to have changed." << endl + << "In constant/polyMesh there are " << meshSolver_.points0().size() + << " points; in the current mesh there are " << mesh().nPoints() + << " points." << exit(FatalError); + } + + // Store the motion state at the beginning of the time-step + if (curTimeIndex_ != this->db().time().timeIndex()) + { + model_.newTime(); + curTimeIndex_ = this->db().time().timeIndex(); + } + + if (db().foundObject("g")) + { + model_.g() = + db().lookupObject("g").value(); + } + + if (test_) + { + label nIter(readLabel(coeffDict().lookup("nIter"))); + + for (label i=0; i(model_.nBodies(), Zero) + ); + } + } + else + { + Field fx(model_.nBodies(), Zero); + + forAll(bodyMeshes_, bi) + { + const label bodyID = bodyMeshes_[bi].bodyID_; + + dictionary forcesDict; + forcesDict.add("type", functionObjects::forces::typeName); + forcesDict.add("patches", bodyMeshes_[bi].patches_); + forcesDict.add("rhoInf", rhoInf_); + forcesDict.add("rho", rhoName_); + forcesDict.add("CofR", vector::zero); + + functionObjects::forces f("forces", db(), forcesDict); + f.calcForcesMoment(); + + fx[bodyID] = spatialVector(f.momentEff(), f.forceEff()); + } + + model_.solve + ( + t.deltaTValue(), + scalarField(model_.nDoF(), Zero), + fx + ); + } + + if (Pstream::master() && model_.report()) + { + forAll(bodyMeshes_, bi) + { + model_.status(bodyMeshes_[bi].bodyID_); + } + } + + // Update the displacements + forAll(bodyMeshes_, bi) + { + forAllConstIter(labelHashSet, bodyMeshes_[bi].patchSet_, iter) + { + label patchi = iter.key(); + + pointField patchPoints0 + ( + meshSolver_.pointDisplacement().boundaryField()[patchi] + .patchInternalField(meshSolver_.points0()) + ); + + meshSolver_.pointDisplacement().boundaryFieldRef()[patchi] == + ( + model_.transformPoints + ( + bodyMeshes_[bi].bodyID_, + patchPoints0 + ) - patchPoints0 + )(); + } + } + + meshSolverPtr_->solve(); +} + + +bool Foam::rigidBodyMeshMotionSolver::writeObject +( + IOstream::streamFormat fmt, + IOstream::versionNumber ver, + IOstream::compressionType cmp +) const +{ + IOdictionary dict + ( + IOobject + ( + "rigidBodyMotionState", + mesh().time().timeName(), + "uniform", + mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ) + ); + + model_.state().write(dict); + return dict.regIOobject::write(); +} + + +bool Foam::rigidBodyMeshMotionSolver::read() +{ + if (motionSolver::read()) + { + model_.read(coeffDict()); + + return true; + } + else + { + return false; + } +} + + +void Foam::rigidBodyMeshMotionSolver::movePoints(const pointField& points) +{ + meshSolverPtr_->movePoints(points); +} + + +void Foam::rigidBodyMeshMotionSolver::updateMesh(const mapPolyMesh& mpm) +{ + meshSolverPtr_->updateMesh(mpm); +} + + +// ************************************************************************* // diff --git a/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.H b/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.H new file mode 100644 index 0000000000..26dfb21899 --- /dev/null +++ b/src/rigidBodyMeshMotion/rigidBodyMeshMotionSolver/rigidBodyMeshMotionSolver.H @@ -0,0 +1,184 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 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 . + +Class + Foam::rigidBodyMeshMotionSolver + +Description + Rigid-body mesh motion solver for fvMesh. + + Applies septernion interpolation of movement as function of distance to the + object surface. + +SourceFiles + rigidBodyMeshMotionSolver.C + +\*---------------------------------------------------------------------------*/ + +#ifndef rigidBodyMeshMotionSolver_H +#define rigidBodyMeshMotionSolver_H + +#include "displacementMotionSolver.H" +#include "rigidBodyMotion.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class rigidBodyMeshMotionSolver Declaration +\*---------------------------------------------------------------------------*/ + +class rigidBodyMeshMotionSolver +: + public motionSolver +{ + //- Class containing the patches and point motion weighting for each body + class bodyMesh + { + //- Name of the body + const word name_; + + //- ID of the body in the RBD::rigidBodyMotion + const label bodyID_; + + //- List of mesh patches associated with this body + const wordReList patches_; + + //- Patches to integrate forces + const labelHashSet patchSet_; + + + public: + + friend class rigidBodyMeshMotionSolver; + + bodyMesh + ( + const polyMesh& mesh, + const word& name, + const label bodyID, + const dictionary& dict + ); + }; + + + // Private data + + //- Rigid-body model + RBD::rigidBodyMotion model_; + + //- List of the bodyMeshes containing the patches and point motion + // weighting for each body + PtrList bodyMeshes_; + + //- Switch for test-mode in which only the + // gravitational body-force is applied + Switch test_; + + //- Reference density required by the forces object for + // incompressible calculations, required if rho == rhoInf + scalar rhoInf_; + + //- Name of density field, optional unless used for an + // incompressible simulation, when this needs to be specified + // as rhoInf + word rhoName_; + + //- Current time index (used for updating) + label curTimeIndex_; + + autoPtr meshSolverPtr_; + + displacementMotionSolver& meshSolver_; + + + // Private Member Functions + + //- Disallow default bitwise copy construct + rigidBodyMeshMotionSolver + ( + const rigidBodyMeshMotionSolver& + ); + + //- Disallow default bitwise assignment + void operator=(const rigidBodyMeshMotionSolver&); + + +public: + + //- Runtime type information + TypeName("rigidBodyMotionSolver"); + + + // Constructors + + //- Construct from polyMesh and IOdictionary + rigidBodyMeshMotionSolver + ( + const polyMesh&, + const IOdictionary& dict + ); + + + //- Destructor + ~rigidBodyMeshMotionSolver(); + + + // Member Functions + + //- Return point location obtained from the current motion field + virtual tmp curPoints() const; + + //- Solve for motion + virtual void solve(); + + //- Write state using given format, version and compression + virtual bool writeObject + ( + IOstream::streamFormat fmt, + IOstream::versionNumber ver, + IOstream::compressionType cmp + ) const; + + //- Read dynamicMeshDict dictionary + virtual bool read(); + + //- Update local data for geometry changes + virtual void movePoints(const pointField&); + + //- Update local data for topology changes + virtual void updateMesh(const mapPolyMesh&); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* //