diff --git a/applications/solvers/incompressible/boundaryFoam/boundaryFoam.C b/applications/solvers/incompressible/boundaryFoam/boundaryFoam.C index 39d214a6e1319e5375b34823cb72dd28d07a7a5d..4e35c5258f70a019fcfb3cec9fc6db8d013073bc 100644 --- a/applications/solvers/incompressible/boundaryFoam/boundaryFoam.C +++ b/applications/solvers/incompressible/boundaryFoam/boundaryFoam.C @@ -45,6 +45,8 @@ Description int main(int argc, char *argv[]) { + argList::noParallel(); + #include "setRootCase.H" #include "createTime.H" diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/reactingTwoPhaseEulerFoam.C b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/reactingTwoPhaseEulerFoam.C index d723a7ce5e44e82ff3346ba8079d624d938bbb3e..d18171a35d3549a38ab290efd70189c429069e85 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/reactingTwoPhaseEulerFoam.C +++ b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/reactingTwoPhaseEulerFoam.C @@ -128,7 +128,7 @@ int main(int argc, char *argv[]) } } - #include "write.H" + runTime.write(); Info<< "ExecutionTime = " << runTime.elapsedCpuTime() diff --git a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/write.H b/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/write.H deleted file mode 100644 index 303661beb6492a7e608f6d5e2385c824bf8dc55f..0000000000000000000000000000000000000000 --- a/applications/solvers/multiphase/reactingEulerFoam/reactingTwoPhaseEulerFoam/write.H +++ /dev/null @@ -1,17 +0,0 @@ - if (runTime.outputTime()) - { - volVectorField Ur - ( - IOobject - ( - "Ur", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - U1 - U2 - ); - - runTime.write(); - } diff --git a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C index 987e26d2e3dc7cb365d6ac35b74b679d827ac715..68c316f24d297383e906f16888c831eddd6a2b45 100644 --- a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C +++ b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C @@ -493,12 +493,15 @@ Foam::label Foam::checkGeometry(const polyMesh& mesh, const bool allGeometry) // Min length scalar minDistSqr = magSqr(1e-6 * globalBb.span()); - // Non-empty directions + // Geometric directions const Vector<label> validDirs = (mesh.geometricD() + Vector<label>::one)/2; - Info<< " Mesh (non-empty, non-wedge) directions " << validDirs << endl; + Info<< " Mesh has " << mesh.nGeometricD() + << " geometric (non-empty/wedge) directions " << validDirs << endl; + // Solution directions const Vector<label> solDirs = (mesh.solutionD() + Vector<label>::one)/2; - Info<< " Mesh (non-empty) directions " << solDirs << endl; + Info<< " Mesh has " << mesh.nSolutionD() + << " solution (non-empty) directions " << solDirs << endl; if (mesh.nGeometricD() < 3) { diff --git a/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C b/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C index aa1488badee267f646100ac76d34dd0514e94d1b..d2f25bb46f5dc4e39532aeac4e08a37a2aa37b59 100644 --- a/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C +++ b/applications/utilities/mesh/manipulation/moveDynamicMesh/moveDynamicMesh.C @@ -133,7 +133,11 @@ int main(int argc, char *argv[]) { Info<< "Time = " << runTime.timeName() << endl; - mesh.update(); + for (int i = 0; i<2; i++) + { + mesh.update(); + } + mesh.checkMesh(true); if (checkAMI) diff --git a/etc/templates/axisymmetricJet/constant/transportProperties b/etc/templates/axisymmetricJet/constant/transportProperties index 455ea9b0bfb25b077cb8aee2342180ee2dddfe85..2a62bb7f6cfa71517180dff74bbc87b7820e7707 100644 --- a/etc/templates/axisymmetricJet/constant/transportProperties +++ b/etc/templates/axisymmetricJet/constant/transportProperties @@ -22,7 +22,7 @@ BirdCarreauCoeffs { nu0 [0 2 -1 0 0 0 0] 1e-03; nuInf [0 2 -1 0 0 0 0] 1e-05; - m [0 0 1 0 0 0 0] 1; + k [0 0 1 0 0 0 0] 1; n [0 0 0 0 0 0 0] 0.5; } diff --git a/etc/templates/closedVolumeRotating/constant/transportProperties b/etc/templates/closedVolumeRotating/constant/transportProperties index 455ea9b0bfb25b077cb8aee2342180ee2dddfe85..2a62bb7f6cfa71517180dff74bbc87b7820e7707 100644 --- a/etc/templates/closedVolumeRotating/constant/transportProperties +++ b/etc/templates/closedVolumeRotating/constant/transportProperties @@ -22,7 +22,7 @@ BirdCarreauCoeffs { nu0 [0 2 -1 0 0 0 0] 1e-03; nuInf [0 2 -1 0 0 0 0] 1e-05; - m [0 0 1 0 0 0 0] 1; + k [0 0 1 0 0 0 0] 1; n [0 0 0 0 0 0 0] 0.5; } diff --git a/etc/templates/inflowOutflow/constant/transportProperties b/etc/templates/inflowOutflow/constant/transportProperties index 455ea9b0bfb25b077cb8aee2342180ee2dddfe85..2a62bb7f6cfa71517180dff74bbc87b7820e7707 100644 --- a/etc/templates/inflowOutflow/constant/transportProperties +++ b/etc/templates/inflowOutflow/constant/transportProperties @@ -22,7 +22,7 @@ BirdCarreauCoeffs { nu0 [0 2 -1 0 0 0 0] 1e-03; nuInf [0 2 -1 0 0 0 0] 1e-05; - m [0 0 1 0 0 0 0] 1; + k [0 0 1 0 0 0 0] 1; n [0 0 0 0 0 0 0] 0.5; } diff --git a/etc/templates/inflowOutflowRotating/constant/transportProperties b/etc/templates/inflowOutflowRotating/constant/transportProperties index 455ea9b0bfb25b077cb8aee2342180ee2dddfe85..2a62bb7f6cfa71517180dff74bbc87b7820e7707 100644 --- a/etc/templates/inflowOutflowRotating/constant/transportProperties +++ b/etc/templates/inflowOutflowRotating/constant/transportProperties @@ -22,7 +22,7 @@ BirdCarreauCoeffs { nu0 [0 2 -1 0 0 0 0] 1e-03; nuInf [0 2 -1 0 0 0 0] 1e-05; - m [0 0 1 0 0 0 0] 1; + k [0 0 1 0 0 0 0] 1; n [0 0 0 0 0 0 0] 0.5; } diff --git a/src/OpenFOAM/meshes/polyMesh/polyMesh.C b/src/OpenFOAM/meshes/polyMesh/polyMesh.C index 1c9a6203f80315b1c3446f9d54594bc0503d0374..4b8012b71d843b075e12ef16ec358098bf0d1f3a 100644 --- a/src/OpenFOAM/meshes/polyMesh/polyMesh.C +++ b/src/OpenFOAM/meshes/polyMesh/polyMesh.C @@ -1483,11 +1483,7 @@ Foam::label Foam::polyMesh::findCell return -1; } - if - ( - Pstream::parRun() - && (decompMode == FACE_DIAG_TRIS || decompMode == CELL_TETS) - ) + if (Pstream::parRun() && decompMode == FACE_DIAG_TRIS) { // Force construction of face-diagonal decomposition before testing // for zero cells. diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index c7100518ab6efaf1cd1578410e3ff6c54dd91476..1b9c242765d4fa990d84b093cef7d78a213e8824 100644 --- a/src/finiteVolume/Make/files +++ b/src/finiteVolume/Make/files @@ -199,6 +199,7 @@ $(derivedFvPatchFields)/waveTransmissive/waveTransmissiveFvPatchFields.C $(derivedFvPatchFields)/waveSurfacePressure/waveSurfacePressureFvPatchScalarField.C $(derivedFvPatchFields)/interstitialInletVelocity/interstitialInletVelocityFvPatchVectorField.C $(derivedFvPatchFields)/prghPressure/prghPressureFvPatchScalarField.C +$(derivedFvPatchFields)/prghTotalPressure/prghTotalPressureFvPatchScalarField.C fvsPatchFields = fields/fvsPatchFields $(fvsPatchFields)/fvsPatchField/fvsPatchFields.C diff --git a/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.C b/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.C new file mode 100644 index 0000000000000000000000000000000000000000..7b8789df44832621abc39b1cf6bc583a55117b79 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.C @@ -0,0 +1,213 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "prghTotalPressureFvPatchScalarField.H" +#include "addToRunTimeSelectionTable.H" +#include "fvPatchFieldMapper.H" +#include "volFields.H" +#include "surfaceFields.H" +#include "uniformDimensionedFields.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF +) +: + fixedValueFvPatchScalarField(p, iF), + UName_("U"), + phiName_("phi"), + rhoName_("rho"), + p0_(p.size(), 0.0) +{} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF, + const dictionary& dict +) +: + fixedValueFvPatchScalarField(p, iF), + UName_(dict.lookupOrDefault<word>("U", "U")), + phiName_(dict.lookupOrDefault<word>("phi", "phi")), + rhoName_(dict.lookupOrDefault<word>("rhoName", "rho")), + p0_("p0", dict, p.size()) +{ + if (dict.found("value")) + { + fvPatchScalarField::operator= + ( + scalarField("value", dict, p.size()) + ); + } + else + { + fvPatchField<scalar>::operator=(p0_); + } +} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const prghTotalPressureFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField<scalar, volMesh>& iF, + const fvPatchFieldMapper& mapper +) +: + fixedValueFvPatchScalarField(ptf, p, iF, mapper), + UName_(ptf.UName_), + phiName_(ptf.phiName_), + rhoName_(ptf.rhoName_), + p0_(ptf.p0_, mapper) +{} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const prghTotalPressureFvPatchScalarField& ptf +) +: + fixedValueFvPatchScalarField(ptf), + UName_(ptf.UName_), + phiName_(ptf.phiName_), + rhoName_(ptf.rhoName_), + p0_(ptf.p0_) +{} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const prghTotalPressureFvPatchScalarField& ptf, + const DimensionedField<scalar, volMesh>& iF +) +: + fixedValueFvPatchScalarField(ptf, iF), + UName_(ptf.UName_), + phiName_(ptf.phiName_), + rhoName_(ptf.rhoName_), + p0_(ptf.p0_) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::prghTotalPressureFvPatchScalarField::autoMap +( + const fvPatchFieldMapper& m +) +{ + fixedValueFvPatchScalarField::autoMap(m); + p0_.autoMap(m); +} + + +void Foam::prghTotalPressureFvPatchScalarField::rmap +( + const fvPatchScalarField& ptf, + const labelList& addr +) +{ + fixedValueFvPatchScalarField::rmap(ptf, addr); + + const prghTotalPressureFvPatchScalarField& tiptf = + refCast<const prghTotalPressureFvPatchScalarField>(ptf); + + p0_.rmap(tiptf.p0_, addr); +} + + +void Foam::prghTotalPressureFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + const scalarField& rhop = + patch().lookupPatchField<volScalarField, scalar>(rhoName_); + + const scalarField& phip = + patch().lookupPatchField<surfaceScalarField, scalar>(phiName_); + + const vectorField& Up = + patch().lookupPatchField<volVectorField, vector>(UName_); + + const uniformDimensionedVectorField& g = + db().lookupObject<uniformDimensionedVectorField>("g"); + + const uniformDimensionedScalarField& hRef = + db().lookupObject<uniformDimensionedScalarField>("hRef"); + + dimensionedScalar ghRef + ( + mag(g.value()) > SMALL + ? g & (cmptMag(g.value())/mag(g.value()))*hRef + : dimensionedScalar("ghRef", g.dimensions()*dimLength, 0) + ); + + operator== + ( + p0_ + - 0.5*rhop*(1.0 - pos(phip))*magSqr(Up) + - rhop*((g.value() & patch().Cf()) - ghRef.value()) + ); + + fixedValueFvPatchScalarField::updateCoeffs(); +} + + +void Foam::prghTotalPressureFvPatchScalarField::write(Ostream& os) const +{ + fvPatchScalarField::write(os); + writeEntryIfDifferent<word>(os, "U", "U", UName_); + writeEntryIfDifferent<word>(os, "phi", "phi", phiName_); + writeEntryIfDifferent<word>(os, "rho", "rho", rhoName_); + p0_.writeEntry("p0", os); + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + makePatchTypeField + ( + fvPatchScalarField, + prghTotalPressureFvPatchScalarField + ); +} + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.H b/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.H new file mode 100644 index 0000000000000000000000000000000000000000..f809029261bebe7383be54b74fd5f804a1455d55 --- /dev/null +++ b/src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.H @@ -0,0 +1,243 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +Class + Foam::prghTotalPressureFvPatchScalarField + +Group + grpGenericBoundaryConditions + +Description + This boundary condition provides static pressure condition for p_rgh, + calculated as: + + \f[ + p_rgh = p - \rho g (h - hRef) + \f] + + where + \vartable + p_rgh | Pseudo hydrostatic pressure [Pa] + p | Static pressure [Pa] + h | Height in the opposite direction to gravity + hRef | Reference height in the opposite direction to gravity + \rho | density + g | acceleration due to gravity [m/s^2] + \endtable + + \heading Patch usage + + \table + Property | Description | Required | Default value + rhoName | rho field name | no | rho + p | static pressure | yes | + \endtable + + Example of the boundary condition specification: + \verbatim + myPatch + { + type prghTotalPressure; + rhoName rho; + p uniform 0; + value uniform 0; // optional initial value + } + \endverbatim + +SeeAlso + Foam::fixedValueFvPatchScalarField + +SourceFiles + prghTotalPressureFvPatchScalarField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef prghTotalPressureFvPatchScalarField_H +#define prghTotalPressureFvPatchScalarField_H + +#include "fixedValueFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class prghTotalPressureFvPatchScalarField Declaration +\*---------------------------------------------------------------------------*/ + +class prghTotalPressureFvPatchScalarField +: + public fixedValueFvPatchScalarField +{ + +protected: + + // Protected data + + //- Name of the velocity field + word UName_; + + //- Name of the flux transporting the field + word phiName_; + + //- Name of phase-fraction field + word rhoName_; + + //- Total pressure + scalarField p0_; + + +public: + + //- Runtime type information + TypeName("prghTotalPressure"); + + + // Constructors + + //- Construct from patch and internal field + prghTotalPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField<scalar, volMesh>& + ); + + //- Construct from patch, internal field and dictionary + prghTotalPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField<scalar, volMesh>&, + const dictionary& + ); + + //- Construct by mapping given + // prghTotalPressureFvPatchScalarField onto a new patch + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField&, + const fvPatch&, + const DimensionedField<scalar, volMesh>&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp<fvPatchScalarField> clone() const + { + return tmp<fvPatchScalarField > + ( + new prghTotalPressureFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField&, + const DimensionedField<scalar, volMesh>& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp<fvPatchScalarField> clone + ( + const DimensionedField<scalar, volMesh>& iF + ) const + { + return tmp<fvPatchScalarField> + ( + new prghTotalPressureFvPatchScalarField(*this, iF) + ); + } + + + // Member functions + + // Access + + //- Return the rhoName + const word& rhoName() const + { + return rhoName_; + } + + //- Return reference to the rhoName to allow adjustment + word& rhoName() + { + return rhoName_; + } + + //- Return the total pressure + const scalarField& p0() const + { + return p0_; + } + + //- Return reference to the total pressure to allow adjustment + scalarField& p0() + { + return p0_; + } + + + // Mapping functions + + //- Map (and resize as needed) from self given a mapping object + virtual void autoMap + ( + const fvPatchFieldMapper& + ); + + //- Reverse map the given fvPatchField onto this fvPatchField + virtual void rmap + ( + const fvPatchScalarField&, + const labelList& + ); + + + // Evaluation functions + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C index 11c8a886858010e448b941b751ae5f224e8eeab4..e5870388908325c1df50c5665659fc8557066bd4 100644 --- a/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C +++ b/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C @@ -81,15 +81,8 @@ Foam::uniformFixedValueFvPatchField<Type>::uniformFixedValueFvPatchField fixedValueFvPatchField<Type>(p, iF), uniformValue_(DataEntry<Type>::New("uniformValue", dict)) { - if (dict.found("value")) - { - fvPatchField<Type>::operator==(Field<Type>("value", dict, p.size())); - } - else - { - const scalar t = this->db().time().timeOutputValue(); - fvPatchField<Type>::operator==(uniformValue_->value(t)); - } + const scalar t = this->db().time().timeOutputValue(); + fvPatchField<Type>::operator==(uniformValue_->value(t)); } diff --git a/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C b/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C index a3351652f9f0204a622ae5297d65a7bf19ba4de9..2841fb85e347f5d3a4a7d54292d7f125b8820a08 100644 --- a/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C +++ b/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C @@ -501,7 +501,7 @@ void Foam::fvSchemes::setFluxRequired(const word& name) const Info<< "Setting fluxRequired for " << name << endl; } - fluxRequired_.add(name, true); + fluxRequired_.add(name, true, true); } diff --git a/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.C b/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.C index cb1fb77a364a67b8d51aacb08420b554cb83508b..6a16f35c1379cffe6a322bccd50466ba2d0a2762 100644 --- a/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.C +++ b/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.C @@ -35,27 +35,31 @@ namespace Foam defineTypeNameAndDebug(turbulenceFields, 0); template<> - const char* NamedEnum<turbulenceFields::compressibleField, 6>::names[] = + const char* NamedEnum<turbulenceFields::compressibleField, 8>::names[] = { - "R", - "devRhoReff", + "k", + "epsilon", "mut", "muEff", "alphat", - "alphaEff" + "alphaEff", + "R", + "devRhoReff" }; - const NamedEnum<turbulenceFields::compressibleField, 6> + const NamedEnum<turbulenceFields::compressibleField, 8> turbulenceFields::compressibleFieldNames_; template<> - const char* NamedEnum<turbulenceFields::incompressibleField, 4>::names[] = + const char* NamedEnum<turbulenceFields::incompressibleField, 6>::names[] = { - "R", - "devReff", + "k", + "epsilon", "nut", - "nuEff" + "nuEff", + "R", + "devReff" }; - const NamedEnum<turbulenceFields::incompressibleField, 4> + const NamedEnum<turbulenceFields::incompressibleField, 6> turbulenceFields::incompressibleFieldNames_; const word turbulenceFields::modelName = turbulenceModel::propertiesName; @@ -174,14 +178,14 @@ void Foam::turbulenceFields::execute() const word& f = iter.key(); switch (compressibleFieldNames_[f]) { - case cfR: + case cfK: { - processField<symmTensor>(f, model.R()); + processField<scalar>(f, model.k()); break; } - case cfDevRhoReff: + case cfEpsilon: { - processField<symmTensor>(f, model.devRhoReff()); + processField<scalar>(f, model.epsilon()); break; } case cfMut: @@ -204,6 +208,16 @@ void Foam::turbulenceFields::execute() processField<scalar>(f, model.alphaEff()); break; } + case cfR: + { + processField<symmTensor>(f, model.R()); + break; + } + case cfDevRhoReff: + { + processField<symmTensor>(f, model.devRhoReff()); + break; + } default: { FatalErrorIn("void Foam::turbulenceFields::execute()") @@ -222,14 +236,14 @@ void Foam::turbulenceFields::execute() const word& f = iter.key(); switch (incompressibleFieldNames_[f]) { - case ifR: + case ifK: { - processField<symmTensor>(f, model.R()); + processField<scalar>(f, model.k()); break; } - case ifDevReff: + case ifEpsilon: { - processField<symmTensor>(f, model.devReff()); + processField<scalar>(f, model.epsilon()); break; } case ifNut: @@ -242,6 +256,16 @@ void Foam::turbulenceFields::execute() processField<scalar>(f, model.nuEff()); break; } + case ifR: + { + processField<symmTensor>(f, model.R()); + break; + } + case ifDevReff: + { + processField<symmTensor>(f, model.devReff()); + break; + } default: { FatalErrorIn("void Foam::turbulenceFields::execute()") @@ -263,15 +287,11 @@ void Foam::turbulenceFields::end() void Foam::turbulenceFields::timeSet() -{ - // Do nothing -} +{} void Foam::turbulenceFields::write() -{ - // Do nothing -} +{} // ************************************************************************* // diff --git a/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.H b/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.H index 45e581cec18e60acf603b6c73c9f4ab607148800..422af29812cc79fce3625a8f0349933c38a02c1e 100644 --- a/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.H +++ b/src/postProcessing/functionObjects/utilities/turbulenceFields/turbulenceFields.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -53,22 +53,24 @@ Description \heading Function object usage \table - Property | Description | Required | Default value - type | type name: processorField | yes | - fields | fields to store (see below) | yes | + Property | Description | Required | Default value + type | type name: processorField | yes | + fields | fields to store (see below) | yes | \endtable Where \c fields can include: \plaintable - R | Stress tensor - devRhoReff | + k | turbulence kinetic energy + epsilon | turbulence kinetic energy dissipation rate + nut | turbulence viscosity (incompressible) + nuEff | effective turbulence viscosity (incompressible) mut | turbulence viscosity (compressible) muEff | effective turbulence viscosity (compressible) alphat | turbulence thermal diffusivity (compressible) alphaEff | effective turbulence thermal diffusivity (compressible) - devReff | - nut | turbulence viscosity (incompressible) - nuEff | effective turbulence viscosity (incompressible) + R | Reynolds stress tensor + devReff | Deviatoric part of the effective Reynolds stress + devRhoReff | Divergence of the Reynolds stress \endplaintable SeeAlso @@ -109,23 +111,27 @@ public: enum compressibleField { - cfR, - cfDevRhoReff, + cfK, + cfEpsilon, cfMut, cfMuEff, cfAlphat, - cfAlphaEff + cfAlphaEff, + cfR, + cfDevRhoReff }; - static const NamedEnum<compressibleField, 6> compressibleFieldNames_; + static const NamedEnum<compressibleField, 8> compressibleFieldNames_; enum incompressibleField { - ifR, - ifDevReff, + ifK, + ifEpsilon, ifNut, - ifNuEff + ifNuEff, + ifR, + ifDevReff }; - static const NamedEnum<incompressibleField, 4> incompressibleFieldNames_; + static const NamedEnum<incompressibleField, 6> incompressibleFieldNames_; static const word modelName; diff --git a/src/sixDoFRigidBodyMotion/Make/files b/src/sixDoFRigidBodyMotion/Make/files index 8e825f0a010760696e12c38e24e7e3f5f5bf36dc..1184c6ad72d33ab7e051c816f76e94b9d53e0fc0 100644 --- a/src/sixDoFRigidBodyMotion/Make/files +++ b/src/sixDoFRigidBodyMotion/Make/files @@ -31,4 +31,10 @@ sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C sixDoFRigidBodyMotionSolver/externalPointEdgePoint.C sixDoFRigidBodyMotionSolver/pointPatchDist.C +sixDoFSolvers/sixDoFSolver/sixDoFSolver.C +sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C +sixDoFSolvers/symplectic/symplectic.C +sixDoFSolvers/CrankNicolson/CrankNicolson.C +sixDoFSolvers/Newmark/Newmark.C + LIB = $(FOAM_LIBBIN)/libsixDoFRigidBodyMotion diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C index dc731f060a39cef715c1c43e1299302974eac315..7f51b9551d45ea72eb1fd0affbe6d7035bea5655 100644 --- a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C +++ b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -204,18 +204,14 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() const pointPatch& ptPatch = this->patch(); // Store the motion state at the beginning of the time-step + bool firstIter = false; if (curTimeIndex_ != t.timeIndex()) { motion_.newTime(); curTimeIndex_ = t.timeIndex(); + firstIter = true; } - // Patch force data is valid for the current positions, so - // calculate the forces on the motion object from this data, then - // update the positions - - motion_.updatePosition(t.deltaTValue(), t.deltaT0Value()); - dictionary forcesDict; forcesDict.add("type", forces::typeName); @@ -241,11 +237,13 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() // scalar ramp = min(max((t.value() - 5)/10, 0), 1); scalar ramp = 1.0; - motion_.updateAcceleration + motion_.update ( + firstIter, ramp*(f.forceEff() + motion_.mass()*g_), ramp*(f.momentEff() + motion_.mass()*(motion_.momentArm() ^ g_)), - t.deltaTValue() + t.deltaTValue(), + t.deltaT0Value() ); Field<vector>::operator= diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C index fb832ca6c80aed91d41ecfd95f26b59dba4987f5..d567e3e04aa9249b364552a5ee5576e35bc8b6a3 100644 --- a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C +++ b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -152,14 +152,14 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() const Time& t = mesh.time(); // Store the motion state at the beginning of the time-step + bool firstIter = false; if (curTimeIndex_ != t.timeIndex()) { motion_.newTime(); curTimeIndex_ = t.timeIndex(); + firstIter = true; } - motion_.updatePosition(t.deltaTValue(), t.deltaT0Value()); - vector gravity = vector::zero; if (db().foundObject<uniformDimensionedVectorField>("g")) @@ -171,11 +171,13 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() } // Do not modify the accelerations, except with gravity, where available - motion_.updateAcceleration + motion_.update ( + firstIter, gravity*motion_.mass(), vector::zero, - t.deltaTValue() + t.deltaTValue(), + t.deltaT0Value() ); Field<vector>::operator= diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C index 8155e218dd74e223dcb8ed591003e258f9cff2ea..e7a8b92d339880ba8820abd4946cd0bd8abf96f6 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C @@ -24,6 +24,7 @@ License \*---------------------------------------------------------------------------*/ #include "sixDoFRigidBodyMotion.H" +#include "sixDoFSolver.H" #include "septernion.H" // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // @@ -84,7 +85,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion() momentOfInertia_(diagTensor::one*VSMALL), aRelax_(1.0), aDamp_(1.0), - report_(false) + report_(false), + solver_(NULL) {} @@ -121,7 +123,8 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion momentOfInertia_(dict.lookup("momentOfInertia")), aRelax_(dict.lookupOrDefault<scalar>("accelerationRelaxation", 1.0)), aDamp_(dict.lookupOrDefault<scalar>("accelerationDamping", 1.0)), - report_(dict.lookupOrDefault<Switch>("report", false)) + report_(dict.lookupOrDefault<Switch>("report", false)), + solver_(sixDoFSolver::New(dict.subDict("solver"), *this)) { addRestraints(dict); @@ -171,6 +174,12 @@ Foam::sixDoFRigidBodyMotion::sixDoFRigidBodyMotion {} +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sixDoFRigidBodyMotion::~sixDoFRigidBodyMotion() +{} + + // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * // void Foam::sixDoFRigidBodyMotion::addRestraints @@ -257,67 +266,46 @@ void Foam::sixDoFRigidBodyMotion::addConstraints } -void Foam::sixDoFRigidBodyMotion::updatePosition +void Foam::sixDoFRigidBodyMotion::updateAcceleration ( - scalar deltaT, - scalar deltaT0 + const vector& fGlobal, + const vector& tauGlobal ) { - // First leapfrog velocity adjust and motion part, required before - // force calculation + static bool first = false; - if (Pstream::master()) - { - v() = tConstraints_ & (v0() + aDamp_*0.5*deltaT0*a()); - pi() = rConstraints_ & (pi0() + aDamp_*0.5*deltaT0*tau()); + // Save the previous iteration accelerations for relaxation + vector aPrevIter = a(); + vector tauPrevIter = tau(); - // Leapfrog move part - centreOfRotation() = centreOfRotation0() + deltaT*v(); + // Calculate new accelerations + a() = fGlobal/mass_; + tau() = (Q().T() & tauGlobal); + applyRestraints(); - // Leapfrog orientation adjustment - Tuple2<tensor, vector> Qpi = rotate(Q0(), pi(), deltaT); - Q() = Qpi.first(); - pi() = rConstraints_ & Qpi.second(); + // Relax accelerations on all but first iteration + if (!first) + { + a() = aRelax_*a() + (1 - aRelax_)*aPrevIter; + tau() = aRelax_*tau() + (1 - aRelax_)*tauPrevIter; } - Pstream::scatter(motionState_); + first = false; } -void Foam::sixDoFRigidBodyMotion::updateAcceleration +void Foam::sixDoFRigidBodyMotion::update ( + bool firstIter, const vector& fGlobal, const vector& tauGlobal, - scalar deltaT + scalar deltaT, + scalar deltaT0 ) { - static bool first = false; - - // Second leapfrog velocity adjust part, required after motion and - // acceleration calculation - if (Pstream::master()) { - // Save the previous iteration accelerations for relaxation - vector aPrevIter = a(); - vector tauPrevIter = tau(); - - // Calculate new accelerations - a() = fGlobal/mass_; - tau() = (Q().T() & tauGlobal); - applyRestraints(); - - // Relax accelerations on all but first iteration - if (!first) - { - a() = aRelax_*a() + (1 - aRelax_)*aPrevIter; - tau() = aRelax_*tau() + (1 - aRelax_)*tauPrevIter; - } - first = false; - - // Correct velocities - v() += tConstraints_ & aDamp_*0.5*deltaT*a(); - pi() += rConstraints_ & aDamp_*0.5*deltaT*tau(); + solver_->solve(firstIter, fGlobal, tauGlobal, deltaT, deltaT0); if (report_) { diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H index b08c7d011ff70dfb1d19546ee835fc3d976cbd60..a91f17e23dcba0650fae3c8bcde5830b6db32714 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -25,25 +25,16 @@ Class Foam::sixDoFRigidBodyMotion Description - Six degree of freedom motion for a rigid body. Angular momentum - stored in body fixed reference frame. Reference orientation of - the body (where Q = I) must align with the cartesian axes such - that the Inertia tensor is in principle component form. - - Symplectic motion as per: - - title = {Symplectic splitting methods for rigid body molecular dynamics}, - publisher = {AIP}, - year = {1997}, - journal = {The Journal of Chemical Physics}, - volume = {107}, - number = {15}, - pages = {5840-5851}, - url = {http://link.aip.org/link/?JCP/107/5840/1}, - doi = {10.1063/1.474310} - - Can add restraints (e.g. a spring) - and constraints (e.g. motion may only be on a plane). + Six degree of freedom motion for a rigid body. + + Angular momentum stored in body fixed reference frame. Reference + orientation of the body (where Q = I) must align with the cartesian axes + such that the Inertia tensor is in principle component form. Can add + restraints (e.g. a spring) and constraints (e.g. motion may only be on a + plane). + + The time-integrator for the motion is run-time selectable with options for + symplectic (explicit), Crank-Nicolson and Newmark schemes. SourceFiles sixDoFRigidBodyMotionI.H @@ -66,12 +57,17 @@ SourceFiles namespace Foam { +// Forward declarations +class sixDoFSolver; + /*---------------------------------------------------------------------------*\ Class sixDoFRigidBodyMotion Declaration \*---------------------------------------------------------------------------*/ class sixDoFRigidBodyMotion { + friend class sixDoFSolver; + // Private data //- Motion state data object @@ -117,6 +113,9 @@ class sixDoFRigidBodyMotion //- Switch to turn reporting of motion data on and off Switch report_; + //- Motion solver + autoPtr<sixDoFSolver> solver_; + // Private Member Functions @@ -144,6 +143,9 @@ class sixDoFRigidBodyMotion //- Apply the restraints to the object void applyRestraints(); + //- Update and relax accelerations from the force and torque + void updateAcceleration(const vector& fGlobal, const vector& tauGlobal); + // Access functions retained as private because of the risk of // confusion over what is a body local frame vector and what is global @@ -176,21 +178,6 @@ class sixDoFRigidBodyMotion //- Return the current torque inline const vector& tau() const; - //- Return the orientation at previous time-step - inline const tensor& Q0() const; - - //- Return the velocity at previous time-step - inline const vector& v0() const; - - //- Return the acceleration at previous time-step - inline const vector& a0() const; - - //- Return the angular momentum at previous time-step - inline const vector& pi0() const; - - //- Return the torque at previous time-step - inline const vector& tau0() const; - // Edit @@ -234,6 +221,10 @@ public: sixDoFRigidBodyMotion(const sixDoFRigidBodyMotion&); + //- Destructor + ~sixDoFRigidBodyMotion(); + + // Member Functions // Access @@ -250,9 +241,6 @@ public: //- Return the current centre of rotation inline const point& centreOfRotation() const; - //- Return the centre of rotation at previous time-step - inline const point& centreOfRotation0() const; - //- Return the initial centre of mass inline const point& initialCentreOfMass() const; @@ -298,18 +286,15 @@ public: // Update state - //- First leapfrog velocity adjust and motion part, required - // before force calculation. Takes old timestep for variable - // timestep cases. - void updatePosition(scalar deltaT, scalar deltaT0); - - //- Second leapfrog velocity adjust part - // required after motion and force calculation - void updateAcceleration + //- Symplectic integration of velocities, orientation and position. + // Changes to Crank-Nicolson integration for subsequent iterations. + void update ( + bool firstIter, const vector& fGlobal, const vector& tauGlobal, - scalar deltaT + scalar deltaT, + scalar deltaT0 ); //- Report the status of the motion diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H index d09e579e8cca4a3f5e51ed724bf8cd819ff8681c..ffe89819fd11e5028dcf1d7bd78c83b113ce15a8 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionI.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2011-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -161,36 +161,6 @@ inline const Foam::vector& Foam::sixDoFRigidBodyMotion::tau() const } -inline const Foam::tensor& Foam::sixDoFRigidBodyMotion::Q0() const -{ - return motionState0_.Q(); -} - - -inline const Foam::vector& Foam::sixDoFRigidBodyMotion::v0() const -{ - return motionState0_.v(); -} - - -inline const Foam::vector& Foam::sixDoFRigidBodyMotion::a0() const -{ - return motionState0_.a(); -} - - -inline const Foam::vector& Foam::sixDoFRigidBodyMotion::pi0() const -{ - return motionState0_.pi(); -} - - -inline const Foam::vector& Foam::sixDoFRigidBodyMotion::tau0() const -{ - return motionState0_.tau(); -} - - inline Foam::point& Foam::sixDoFRigidBodyMotion::initialCentreOfRotation() { return initialCentreOfRotation_; @@ -261,12 +231,6 @@ inline const Foam::point& Foam::sixDoFRigidBodyMotion::centreOfRotation() const } -inline const Foam::point& Foam::sixDoFRigidBodyMotion::centreOfRotation0() const -{ - return motionState0_.centreOfRotation(); -} - - inline const Foam::point& Foam::sixDoFRigidBodyMotion::initialCentreOfMass() const { diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C index 15bcc0df9316f9c586f235700317d4a2c0f71e9f..1d3d01856d0a0322c6badfa8bdda2a3e910d71c7 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013-2014 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -85,6 +85,7 @@ Foam::sixDoFRigidBodyMotionSolver::sixDoFRigidBodyMotionSolver patchSet_(mesh.boundaryMesh().patchSet(patches_)), di_(readScalar(coeffDict().lookup("innerDistance"))), do_(readScalar(coeffDict().lookup("outerDistance"))), + test_(coeffDict().lookupOrDefault<Switch>("test", false)), rhoInf_(1.0), rhoName_(coeffDict().lookupOrDefault<word>("rhoName", "rho")), scale_ @@ -179,31 +180,15 @@ void Foam::sixDoFRigidBodyMotionSolver::solve() << " points." << exit(FatalError); } - // Store the motion state at the beginning of the time-step + // Store the motion state at the beginning of the time-stepbool + bool firstIter = false; if (curTimeIndex_ != this->db().time().timeIndex()) { motion_.newTime(); curTimeIndex_ = this->db().time().timeIndex(); + firstIter = true; } - // Patch force data is valid for the current positions, so - // calculate the forces on the motion object from this data, then - // update the positions - - motion_.updatePosition(t.deltaTValue(), t.deltaT0Value()); - - dictionary forcesDict; - - forcesDict.add("type", forces::typeName); - forcesDict.add("patches", patches_); - forcesDict.add("rhoInf", rhoInf_); - forcesDict.add("rhoName", rhoName_); - forcesDict.add("CofR", motion_.centreOfRotation()); - - forces f("forces", db(), forcesDict); - - f.calcForcesMoment(); - dimensionedVector g("g", dimAcceleration, vector::zero); if (db().foundObject<uniformDimensionedVectorField>("g")) @@ -218,12 +203,44 @@ void Foam::sixDoFRigidBodyMotionSolver::solve() // scalar ramp = min(max((this->db().time().value() - 5)/10, 0), 1); scalar ramp = 1.0; - motion_.updateAcceleration - ( - ramp*(f.forceEff() + motion_.mass()*g.value()), - ramp*(f.momentEff() + motion_.mass()*(motion_.momentArm() ^ g.value())), - t.deltaTValue() - ); + if (test_) + { + motion_.update + ( + firstIter, + ramp*(motion_.mass()*g.value()), + ramp*(motion_.mass()*(motion_.momentArm() ^ g.value())), + t.deltaTValue(), + t.deltaT0Value() + ); + } + else + { + dictionary forcesDict; + + forcesDict.add("type", forces::typeName); + forcesDict.add("patches", patches_); + forcesDict.add("rhoInf", rhoInf_); + forcesDict.add("rhoName", rhoName_); + forcesDict.add("CofR", motion_.centreOfRotation()); + + forces f("forces", db(), forcesDict); + + f.calcForcesMoment(); + + motion_.update + ( + firstIter, + ramp*(f.forceEff() + motion_.mass()*g.value()), + ramp + *( + f.momentEff() + + motion_.mass()*(motion_.momentArm() ^ g.value()) + ), + t.deltaTValue(), + t.deltaT0Value() + ); + } // Update the displacements pointDisplacement_.internalField() = diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H index 4f4e1ca31c5a80e14eaee7ab06f00cbfa8ae8685..2ad4921aaa961183ed9b1040a7da90fdde1c764f 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2013 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2013-2015 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -70,6 +70,10 @@ class sixDoFRigidBodyMotionSolver //- Outer morphing distance (limit of linear interpolation region) const scalar do_; + //- 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 rhoName == rhoInf scalar rhoInf_; diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C new file mode 100644 index 0000000000000000000000000000000000000000..366266e014b866dd63d9d4fae0a110e85bbb35f0 --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.C @@ -0,0 +1,94 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "CrankNicolson.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + defineTypeNameAndDebug(CrankNicolson, 0); + addToRunTimeSelectionTable(sixDoFSolver, CrankNicolson, dictionary); +} +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::CrankNicolson::CrankNicolson +( + const dictionary& dict, + sixDoFRigidBodyMotion& body +) +: + sixDoFSolver(body), + aoc_(dict.lookupOrDefault<scalar>("aoc", 0.5)), + voc_(dict.lookupOrDefault<scalar>("voc", 0.5)) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::CrankNicolson::~CrankNicolson() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::sixDoFSolvers::CrankNicolson::solve +( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 +) +{ + // Update the linear acceleration and torque + updateAcceleration(fGlobal, tauGlobal); + + // Correct linear velocity + v() = tConstraints() + & (v0() + aDamp()*deltaT*(aoc_*a() + (1 - aoc_)*a0())); + + // Correct angular momentum + pi() = rConstraints() + & (pi0() + aDamp()*deltaT*(aoc_*tau() + (1 - aoc_)*tau0())); + + // Correct position + centreOfRotation() = + centreOfRotation0() + deltaT*(voc_*v() + (1 - voc_)*v0()); + + // Correct orientation + Tuple2<tensor, vector> Qpi = + rotate(Q0(), (voc_*pi() + (1 - voc_)*pi0()), deltaT); + Q() = Qpi.first(); +} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H new file mode 100644 index 0000000000000000000000000000000000000000..a8dd2ccbff3d5c11a7553b1f7475de82fd3bada8 --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/CrankNicolson/CrankNicolson.H @@ -0,0 +1,126 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +Class + Foam::sixDoFSolvers::CrankNicolson + +Description + Crank-Nicolson 2nd-order time-integrator for 6DoF solid-body motion. + + The off-centering coefficients for acceleration (velocity integration) and + velocity (position/orientation integration) may be specified but default + values of 0.5 for each are used if they are not specified. With the default + off-centering this scheme is equivalent to the Newmark scheme with default + coefficients. + + Example specification in dynamicMeshDict: + \verbatim + solver + { + type CrankNicolson; + aoc 0.5; // Acceleration off-centering coefficient + voc 0.5; // Velocity off-centering coefficient + } + \endverbatim + +SeeAlso + Foam::sixDoFSolvers::Newmark + +SourceFiles + CrankNicolson.C + +\*---------------------------------------------------------------------------*/ + +#ifndef CrankNicolson_H +#define CrankNicolson_H + +#include "sixDoFSolver.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + +/*---------------------------------------------------------------------------*\ + Class CrankNicolson Declaration +\*---------------------------------------------------------------------------*/ + +class CrankNicolson +: + public sixDoFSolver +{ + // Private data + + //- Acceleration off-centering coefficient (default: 0.5) + const scalar aoc_; + + //- Velocity off-centering coefficient (default: 0.5) + const scalar voc_; + + +public: + + //- Runtime type information + TypeName("CrankNicolson"); + + + // Constructors + + //- Construct from a dictionary and the body + CrankNicolson + ( + const dictionary& dict, + sixDoFRigidBodyMotion& body + ); + + + //- Destructor + virtual ~CrankNicolson(); + + + // Member Functions + + //- Drag coefficient + virtual void solve + ( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 + ); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace sixDoFSolvers +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C new file mode 100644 index 0000000000000000000000000000000000000000..4f39b3266c92c8a9ef0fe8c670dc74437d6d07fc --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.C @@ -0,0 +1,115 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "Newmark.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + defineTypeNameAndDebug(Newmark, 0); + addToRunTimeSelectionTable(sixDoFSolver, Newmark, dictionary); +} +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::Newmark::Newmark +( + const dictionary& dict, + sixDoFRigidBodyMotion& body +) +: + sixDoFSolver(body), + gamma_(dict.lookupOrDefault<scalar>("gamma", 0.5)), + beta_ + ( + max + ( + 0.25*sqr(gamma_ + 0.5), + dict.lookupOrDefault<scalar>("beta", 0.25) + ) + ) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::Newmark::~Newmark() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::sixDoFSolvers::Newmark::solve +( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 +) +{ + // Update the linear acceleration and torque + updateAcceleration(fGlobal, tauGlobal); + + // Correct linear velocity + v() = + tConstraints() + & (v0() + aDamp()*deltaT*(gamma_*a() + (1 - gamma_)*a0())); + + // Correct angular momentum + pi() = + rConstraints() + & (pi0() + aDamp()*deltaT*(gamma_*tau() + (1 - gamma_)*tau0())); + + // Correct position + centreOfRotation() = + centreOfRotation0() + + ( + tConstraints() + & ( + deltaT*v0() + + sqr(deltaT)*beta_*a() + sqr(deltaT)*(0.5 - beta_)*a0() + ) + ); + + // Correct orientation + vector piDeltaT = + rConstraints() + & ( + deltaT*pi0() + + sqr(deltaT)*beta_*tau() + sqr(deltaT)*(0.5 - beta_)*tau0() + ); + Tuple2<tensor, vector> Qpi = rotate(Q0(), piDeltaT, 1); + Q() = Qpi.first(); +} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H new file mode 100644 index 0000000000000000000000000000000000000000..e605fcc294a3103afb4a4dafa633636ca622b577 --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/Newmark/Newmark.H @@ -0,0 +1,124 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +Class + Foam::sixDoFSolvers::Newmark + +Description + Newmark 2nd-order time-integrator for 6DoF solid-body motion. + + Reference: + \verbatim + Newmark, N. M. (1959). + A method of computation for structural dynamics. + Journal of the Engineering Mechanics Division, 85(3), 67-94. + \endverbatim + + Example specification in dynamicMeshDict: + \verbatim + solver + { + type Newmark; + gamma 0.5; // Velocity integration coefficient + beta 0.25; // Position integration coefficient + } + \endverbatim + +SourceFiles + Newmark.C + +\*---------------------------------------------------------------------------*/ + +#ifndef Newmark_H +#define Newmark_H + +#include "sixDoFSolver.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + +/*---------------------------------------------------------------------------*\ + Class Newmark Declaration +\*---------------------------------------------------------------------------*/ + +class Newmark +: + public sixDoFSolver +{ + // Private data + + //- Coefficient for velocity integration (default: 0.5) + const scalar gamma_; + + //- Coefficient for position and orientation integration (default: 0.25) + const scalar beta_; + + +public: + + //- Runtime type information + TypeName("Newmark"); + + + // Constructors + + //- Construct from a dictionary and the body + Newmark + ( + const dictionary& dict, + sixDoFRigidBodyMotion& body + ); + + + //- Destructor + virtual ~Newmark(); + + + // Member Functions + + //- Drag coefficient + virtual void solve + ( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 + ); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace sixDoFSolvers +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C new file mode 100644 index 0000000000000000000000000000000000000000..329a9aaea185bca039d59b8ea579e4da598f963d --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/newSixDoFSolver.C @@ -0,0 +1,57 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "sixDoFSolver.H" + +// * * * * * * * * * * * * * * * * Selector * * * * * * * * * * * * * * * * // + +Foam::autoPtr<Foam::sixDoFSolver> Foam::sixDoFSolver::New +( + const dictionary& dict, + sixDoFRigidBodyMotion& body +) +{ + word sixDoFSolverType(dict.lookup("type")); + + Info<< "Selecting sixDoFSolver " << sixDoFSolverType << endl; + + dictionaryConstructorTable::iterator cstrIter = + dictionaryConstructorTablePtr_->find(sixDoFSolverType); + + if (cstrIter == dictionaryConstructorTablePtr_->end()) + { + FatalErrorIn("sixDoFSolver::New") + << "Unknown sixDoFSolverType type " + << sixDoFSolverType << endl << endl + << "Valid sixDoFSolver types are : " << endl + << dictionaryConstructorTablePtr_->sortedToc() + << exit(FatalError); + } + + return cstrIter()(dict, body); +} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C new file mode 100644 index 0000000000000000000000000000000000000000..48fbba16f688070a2f40eb3f29dbbffd9a30d00f --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.C @@ -0,0 +1,51 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "sixDoFSolver.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ + defineTypeNameAndDebug(sixDoFSolver, 0); + defineRunTimeSelectionTable(sixDoFSolver, dictionary); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sixDoFSolver::sixDoFSolver(sixDoFRigidBodyMotion& body) +: + body_(body) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sixDoFSolver::~sixDoFSolver() +{} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H new file mode 100644 index 0000000000000000000000000000000000000000..cb00660314a1be4042070224c4cc0aed1f75cfc1 --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolver.H @@ -0,0 +1,190 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +Class + Foam::sixDoFSolver + +Description + +SourceFiles + sixDoFSolver.C + newSixDoFSolver.C + +\*---------------------------------------------------------------------------*/ + +#ifndef sixDoFSolver_H +#define sixDoFSolver_H + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "sixDoFRigidBodyMotion.H" +#include "runTimeSelectionTables.H" + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class sixDoFSolver Declaration +\*---------------------------------------------------------------------------*/ + +class sixDoFSolver +{ +protected: + + // Protected data + + //- The rigid body + sixDoFRigidBodyMotion& body_; + + + // Protected member functions + + //- Return the current centre of rotation + inline point& centreOfRotation(); + + //- Return the orientation + inline tensor& Q(); + + //- Return non-const access to vector + inline vector& v(); + + //- Return non-const access to acceleration + inline vector& a(); + + //- Return non-const access to angular momentum + inline vector& pi(); + + //- Return non-const access to torque + inline vector& tau(); + + //- Return the centre of rotation at previous time-step + inline const point& centreOfRotation0() const; + + //- Return the orientation at previous time-step + inline const tensor& Q0() const; + + //- Return the velocity at previous time-step + inline const vector& v0() const; + + //- Return the acceleration at previous time-step + inline const vector& a0() const; + + //- Return the angular momentum at previous time-step + inline const vector& pi0() const; + + //- Return the torque at previous time-step + inline const vector& tau0() const; + + //- Acceleration damping coefficient (for steady-state simulations) + inline scalar aDamp() const; + + //- Translational constraint tensor + inline tensor tConstraints() const; + + //- Rotational constraint tensor + inline tensor rConstraints() const; + + //- Apply rotation tensors to Q0 for the given torque (pi) and deltaT + // and return the rotated Q and pi as a tuple + inline Tuple2<tensor, vector> rotate + ( + const tensor& Q0, + const vector& pi, + const scalar deltaT + ) const; + + //- Update and relax accelerations from the force and torque + inline void updateAcceleration + ( + const vector& fGlobal, + const vector& tauGlobal + ); + + +public: + + //- Runtime type information + TypeName("sixDoFSolver"); + + + // Declare runtime construction + + declareRunTimeSelectionTable + ( + autoPtr, + sixDoFSolver, + dictionary, + ( + const dictionary& dict, + sixDoFRigidBodyMotion& body + ), + (dict, body) + ); + + + // Constructors + + // Construct for given body + sixDoFSolver(sixDoFRigidBodyMotion& body); + + + //- Destructor + virtual ~sixDoFSolver(); + + + // Selectors + + static autoPtr<sixDoFSolver> New + ( + const dictionary& dict, + sixDoFRigidBodyMotion& body + ); + + + // Member Functions + + //- Drag coefficient + virtual void solve + ( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 + ) = 0; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#include "sixDoFSolverI.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H new file mode 100644 index 0000000000000000000000000000000000000000..0301b2b3b2bbb1a66df188122ffe7e4f583c2bdb --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/sixDoFSolver/sixDoFSolverI.H @@ -0,0 +1,131 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +// * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * // + +inline Foam::point& Foam::sixDoFSolver::centreOfRotation() +{ + return body_.motionState_.centreOfRotation(); +} + +inline Foam::tensor& Foam::sixDoFSolver::Q() +{ + return body_.motionState_.Q(); +} + +inline Foam::vector& Foam::sixDoFSolver::v() +{ + return body_.motionState_.v(); +} + +inline Foam::vector& Foam::sixDoFSolver::a() +{ + return body_.motionState_.a(); +} + +inline Foam::vector& Foam::sixDoFSolver::pi() +{ + return body_.motionState_.pi(); +} + +inline Foam::vector& Foam::sixDoFSolver::tau() +{ + return body_.motionState_.tau(); +} + + +inline const Foam::point& Foam::sixDoFSolver::centreOfRotation0() const +{ + return body_.motionState0_.centreOfRotation(); +} + +inline const Foam::tensor& Foam::sixDoFSolver::Q0() const +{ + return body_.motionState0_.Q(); +} + + +inline const Foam::vector& Foam::sixDoFSolver::v0() const +{ + return body_.motionState0_.v(); +} + + +inline const Foam::vector& Foam::sixDoFSolver::a0() const +{ + return body_.motionState0_.a(); +} + + +inline const Foam::vector& Foam::sixDoFSolver::pi0() const +{ + return body_.motionState0_.pi(); +} + + +inline const Foam::vector& Foam::sixDoFSolver::tau0() const +{ + return body_.motionState0_.tau(); +} + +inline Foam::scalar Foam::sixDoFSolver::aDamp() const +{ + return body_.aDamp_; +} + +inline Foam::tensor Foam::sixDoFSolver::tConstraints() const +{ + return body_.tConstraints_; +} + +inline Foam::tensor Foam::sixDoFSolver::rConstraints() const +{ + return body_.rConstraints_; +} + +//- Apply rotation tensors to Q0 for the given torque (pi) and deltaT +// and return the rotated Q and pi as a tuple +inline Foam::Tuple2<Foam::tensor, Foam::vector> Foam::sixDoFSolver::rotate +( + const tensor& Q0, + const vector& pi, + const scalar deltaT +) const +{ + return body_.rotate(Q0, pi, deltaT); +} + +//- Update and relax accelerations from the force and torque +inline void Foam::sixDoFSolver::updateAcceleration +( + const vector& fGlobal, + const vector& tauGlobal +) +{ + body_.updateAcceleration(fGlobal, tauGlobal); +} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C new file mode 100644 index 0000000000000000000000000000000000000000..61d94abc7a976d3af278b91eea2e051c6ece6ceb --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.C @@ -0,0 +1,102 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "symplectic.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + defineTypeNameAndDebug(symplectic, 0); + addToRunTimeSelectionTable(sixDoFSolver, symplectic, dictionary); +} +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::symplectic::symplectic +( + const dictionary& dict, + sixDoFRigidBodyMotion& body +) +: + sixDoFSolver(body) +{} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::sixDoFSolvers::symplectic::~symplectic() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::sixDoFSolvers::symplectic::solve +( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 +) +{ + if (!firstIter) + { + FatalErrorIn("sixDoFSolvers::symplectic::solve") + << "The symplectic integrator is explicit " + "and can only be solved once per time-step" + << exit(FatalError); + } + + // First simplectic step: + // Half-step for linear and angular velocities + // Update position and orientation + + v() = tConstraints() & (v0() + aDamp()*0.5*deltaT0*a0()); + pi() = rConstraints() & (pi0() + aDamp()*0.5*deltaT0*tau0()); + + centreOfRotation() = centreOfRotation0() + deltaT*v(); + + Tuple2<tensor, vector> Qpi = rotate(Q0(), pi(), deltaT); + Q() = Qpi.first(); + pi() = rConstraints() & Qpi.second(); + + // Update the linear acceleration and torque + updateAcceleration(fGlobal, tauGlobal); + + // Second simplectic step: + // Complete update of linear and angular velocities + + v() += tConstraints() & aDamp()*0.5*deltaT*a(); + pi() += rConstraints() & aDamp()*0.5*deltaT*tau(); +} + + +// ************************************************************************* // diff --git a/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H new file mode 100644 index 0000000000000000000000000000000000000000..d03b77bbd59593ea9e914fa828caba5dc2a40a30 --- /dev/null +++ b/src/sixDoFRigidBodyMotion/sixDoFSolvers/symplectic/symplectic.H @@ -0,0 +1,123 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2015 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 <http://www.gnu.org/licenses/>. + +Class + Foam::sixDoFSolvers::symplectic + +Description + Symplectic 2nd-order explicit time-integrator for 6DoF solid-body motion. + + Reference: + \verbatim + Dullweber, A., Leimkuhler, B., & McLachlan, R. (1997). + Symplectic splitting methods for rigid body molecular dynamics. + The Journal of chemical physics, 107(15), 5840-5851. + \endverbatim + + Can only be used for explicit integration of the motion of the body, + i.e. may only be called once per time-step, no outer-correctors may be + applied. For implicit integration with outer-correctors choose either + CrankNicolson or Newmark schemes. + + Example specification in dynamicMeshDict: + \verbatim + solver + { + type symplectic; + } + \endverbatim + +SeeAlso + Foam::sixDoFSolvers::CrankNicolson + Foam::sixDoFSolvers::Newmark + +SourceFiles + symplectic.C + +\*---------------------------------------------------------------------------*/ + +#ifndef symplectic_H +#define symplectic_H + +#include "sixDoFSolver.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace sixDoFSolvers +{ + +/*---------------------------------------------------------------------------*\ + Class symplectic Declaration +\*---------------------------------------------------------------------------*/ + +class symplectic +: + public sixDoFSolver +{ + +public: + + //- Runtime type information + TypeName("symplectic"); + + + // Constructors + + //- Construct from a dictionary and the body + symplectic + ( + const dictionary& dict, + sixDoFRigidBodyMotion& body + ); + + + //- Destructor + virtual ~symplectic(); + + + // Member Functions + + //- Drag coefficient + virtual void solve + ( + bool firstIter, + const vector& fGlobal, + const vector& tauGlobal, + scalar deltaT, + scalar deltaT0 + ); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace sixDoFSolvers +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict b/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict index c2720340129d77645e823c7bf9773cadb2202671..13d0937968571b27eab91dacefac70fe0672c97a 100644 --- a/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict +++ b/tutorials/incompressible/pimpleDyMFoam/wingMotion/wingMotion2D_pimpleDyMFoam/constant/dynamicMeshDict @@ -41,6 +41,11 @@ sixDoFRigidBodyMotionCoeffs rhoInf 1; report on; + solver + { + type symplectic; + } + constraints { yLine diff --git a/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict b/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict index 7e7f4969af5564d68e522fad2b52d023ff5db4f9..02088e1075be36c5154c92b08985f908b1f97488 100644 --- a/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict +++ b/tutorials/multiphase/interDyMFoam/ras/DTCHull/constant/dynamicMeshDict @@ -31,9 +31,16 @@ sixDoFRigidBodyMotionCoeffs momentOfInertia (40 921 921); rhoInf 1; report on; - accelerationRelaxation 0.3; + value uniform (0 0 0); + accelerationRelaxation 0.4; + + solver + { + type Newmark; + } + constraints { zAxis diff --git a/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution index afd4c5c21931650b8f612585f08af661207917b6..4aaf52f6d20284b6028d368341a38367bb852767 100644 --- a/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution +++ b/tutorials/multiphase/interDyMFoam/ras/DTCHull/system/fvSolution @@ -60,7 +60,7 @@ solvers cacheAgglomeration true; tolerance 5e-8; - relTol 0.001; + relTol 0; }; p_rghFinal @@ -86,8 +86,8 @@ PIMPLE { momentumPredictor no; - nOuterCorrectors 1; - nCorrectors 3; + nOuterCorrectors 3; + nCorrectors 1; nNonOrthogonalCorrectors 0; correctPhi yes; diff --git a/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict b/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict index 68f291d7762b8f7e65b310b8999f4530cfdad3fa..090aa63d16be1aaffe93f0b4c86e12f9a5a857da 100644 --- a/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict +++ b/tutorials/multiphase/interDyMFoam/ras/floatingObject/constant/dynamicMeshDict @@ -59,7 +59,12 @@ sixDoFRigidBodyMotionCoeffs }; report on; - accelerationRelaxation 0.3; + accelerationRelaxation 0.7; + + solver + { + type Newmark; + } constraints { diff --git a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution index bc278a91fa4a94e324454178c182a4a28e7bc0db..ae6f2be175d2ac5f2f149077cd96c80ac550582d 100644 --- a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution +++ b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution @@ -117,7 +117,7 @@ solvers PIMPLE { momentumPredictor no; - nOuterCorrectors 5; + nOuterCorrectors 3; nCorrectors 1; nNonOrthogonalCorrectors 0; correctPhi yes; diff --git a/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict b/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict deleted file mode 100644 index 7e7f4969af5564d68e522fad2b52d023ff5db4f9..0000000000000000000000000000000000000000 --- a/tutorials/multiphase/interFoam/ras/DTCHull/constant/dynamicMeshDict +++ /dev/null @@ -1,67 +0,0 @@ -/*--------------------------------*- C++ -*----------------------------------*\ -| ========= | | -| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | -| \\ / O peration | Version: dev | -| \\ / A nd | Web: www.OpenFOAM.org | -| \\/ M anipulation | | -\*---------------------------------------------------------------------------*/ -FoamFile -{ - version 2.0; - format ascii; - class dictionary; - object dynamicMeshDict; -} -// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - -dynamicFvMesh dynamicMotionSolverFvMesh; - -motionSolverLibs ("libsixDoFRigidBodyMotion.so"); - -solver sixDoFRigidBodyMotion; - -sixDoFRigidBodyMotionCoeffs -{ - patches (hull); - innerDistance 0.3; - outerDistance 1; - - centreOfMass (2.929541 0 0.2); - mass 412.73; - momentOfInertia (40 921 921); - rhoInf 1; - report on; - accelerationRelaxation 0.3; - value uniform (0 0 0); - - constraints - { - zAxis - { - sixDoFRigidBodyMotionConstraint line; - direction (0 0 1); - } - yPlane - { - sixDoFRigidBodyMotionConstraint axis; - axis (0 1 0); - } - } - - restraints - { - translationDamper - { - sixDoFRigidBodyMotionRestraint linearDamper; - coeff 8596; - } - rotationDamper - { - sixDoFRigidBodyMotionRestraint sphericalAngularDamper; - coeff 11586; - } - } -} - - -// ************************************************************************* // diff --git a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh index de5ac9e4379285e390239add1bd4e2ac6fe8c51e..d78578a44083638332e50f1175d5aaf6fcaabf7d 100644 --- a/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh +++ b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh @@ -22,19 +22,22 @@ boundaryField { inlet { - type fixedFluxPressure; - value $internalField; + type fixedFluxPressure; + value $internalField; } outlet { - type prghPressure; - p $internalField; - value $internalField; + type prghTotalPressure; + p0 $internalField; + U U.air; + phi phi.air; + rho rho.air; + value $internalField; } walls { - type fixedFluxPressure; - value $internalField; + type fixedFluxPressure; + value $internalField; } } diff --git a/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution b/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution index 6cf90aa61f3d9656df00ad8bd3f1d185e9cc46bc..2486938d7b52da8bc030e59d25cfc0262da3de61 100644 --- a/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution +++ b/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution @@ -87,6 +87,11 @@ PIMPLE relaxationFactors { + fields + { + iDmdt 1; + } + equations { ".*" 1;