From d3e5ea4914a78152c17be6cc438540d525b9ba5c Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Mon, 12 Oct 2015 21:26:40 +0100 Subject: [PATCH 01/16] fvSchemes: setFluxRequired now adds entry quietly --- src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C b/src/finiteVolume/finiteVolume/fvSchemes/fvSchemes.C index a3351652f9..2841fb85e3 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); } From 9427d134f99b5d9bce823b96144b2a961df789c5 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Mon, 12 Oct 2015 21:27:42 +0100 Subject: [PATCH 02/16] sixDoFRigidBodyMotion: Time integration now switches between symplectic and Crank-Nicolson For explicit motion (and the first iteration of iterative motion correction) the 2nd-order symplectic motion integrator is used. For iterative correction a form of lagged Crank-Nicolson is used in which the current time-step values correspond to the current iteration. This converges to a 2nd-order implicit solution. --- ...gidBodyDisplacementPointPatchVectorField.C | 19 +++-- ...gidBodyDisplacementPointPatchVectorField.C | 10 ++- .../sixDoFRigidBodyMotion.C | 60 +++++++++++---- .../sixDoFRigidBodyMotion.H | 19 +++-- .../sixDoFRigidBodyMotionSolver.C | 76 ++++++++++++------- .../sixDoFRigidBodyMotionSolver.H | 4 + 6 files changed, 130 insertions(+), 58 deletions(-) diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/sixDoFRigidBodyDisplacement/sixDoFRigidBodyDisplacementPointPatchVectorField.C index dc731f060a..a22e4e5865 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); @@ -228,6 +224,11 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() f.calcForcesMoment(); + // 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(firstIter, t.deltaTValue(), t.deltaT0Value()); + // Get the forces on the patch faces at the current positions if (lookupGravity_ == 1) @@ -243,9 +244,11 @@ void sixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() motion_.updateAcceleration ( + firstIter, ramp*(f.forceEff() + motion_.mass()*g_), ramp*(f.momentEff() + motion_.mass()*(motion_.momentArm() ^ g_)), - t.deltaTValue() + t.deltaTValue(), + t.deltaT0Value() ); Field::operator= diff --git a/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C b/src/sixDoFRigidBodyMotion/pointPatchFields/derived/uncoupledSixDoFRigidBodyDisplacement/uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField.C index fb832ca6c8..03bdfcd779 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,13 +152,15 @@ 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()); + motion_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value()); vector gravity = vector::zero; @@ -173,9 +175,11 @@ void uncoupledSixDoFRigidBodyDisplacementPointPatchVectorField::updateCoeffs() // Do not modify the accelerations, except with gravity, where available motion_.updateAcceleration ( + firstIter, gravity*motion_.mass(), vector::zero, - t.deltaTValue() + t.deltaTValue(), + t.deltaT0Value() ); Field::operator= diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C index 8155e218dd..9bb33ecbf9 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.C @@ -259,22 +259,38 @@ void Foam::sixDoFRigidBodyMotion::addConstraints void Foam::sixDoFRigidBodyMotion::updatePosition ( + bool firstIter, scalar deltaT, scalar deltaT0 ) { - // First leapfrog velocity adjust and motion part, required before - // force calculation - if (Pstream::master()) { - v() = tConstraints_ & (v0() + aDamp_*0.5*deltaT0*a()); - pi() = rConstraints_ & (pi0() + aDamp_*0.5*deltaT0*tau()); + if (firstIter) + { + // First simplectic step: + // Half-step for linear and angular velocities + // Update position and orientation - // Leapfrog move part - centreOfRotation() = centreOfRotation0() + deltaT*v(); + v() = tConstraints_ & (v0() + aDamp_*0.5*deltaT0*a()); + pi() = rConstraints_ & (pi0() + aDamp_*0.5*deltaT0*tau()); - // Leapfrog orientation adjustment + centreOfRotation() = centreOfRotation0() + deltaT*v(); + } + else + { + // For subsequent iterations use Crank-Nicolson + + v() = tConstraints_ + & (v0() + aDamp_*0.5*deltaT*(a() + motionState0_.a())); + pi() = rConstraints_ + & (pi0() + aDamp_*0.5*deltaT*(tau() + motionState0_.tau())); + + centreOfRotation() = + centreOfRotation0() + 0.5*deltaT*(v() + motionState0_.v()); + } + + // Correct orientation Tuple2 Qpi = rotate(Q0(), pi(), deltaT); Q() = Qpi.first(); pi() = rConstraints_ & Qpi.second(); @@ -286,16 +302,15 @@ void Foam::sixDoFRigidBodyMotion::updatePosition void Foam::sixDoFRigidBodyMotion::updateAcceleration ( + 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 @@ -315,9 +330,23 @@ void Foam::sixDoFRigidBodyMotion::updateAcceleration } first = false; - // Correct velocities - v() += tConstraints_ & aDamp_*0.5*deltaT*a(); - pi() += rConstraints_ & aDamp_*0.5*deltaT*tau(); + if (firstIter) + { + // Second simplectic step: + // Complete update of linear and angular velocities + + v() += tConstraints_ & aDamp_*0.5*deltaT*a(); + pi() += rConstraints_ & aDamp_*0.5*deltaT*tau(); + } + else + { + // For subsequent iterations use Crank-Nicolson + + v() = tConstraints_ + & (v0() + aDamp_*0.5*deltaT*(a() + motionState0_.a())); + pi() = rConstraints_ + & (pi0() + aDamp_*0.5*deltaT*(tau() + motionState0_.tau())); + } if (report_) { @@ -329,6 +358,7 @@ void Foam::sixDoFRigidBodyMotion::updateAcceleration } + void Foam::sixDoFRigidBodyMotion::status() const { Info<< "6-DoF rigid body motion" << nl diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion/sixDoFRigidBodyMotion.H index b08c7d011f..70f2e7674f 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 @@ -42,6 +42,9 @@ Description url = {http://link.aip.org/link/?JCP/107/5840/1}, doi = {10.1063/1.474310} + Changes to Crank-Nicolson integration if subsequent iterations to handle + implicit coupling between pressure and motion. + Can add restraints (e.g. a spring) and constraints (e.g. motion may only be on a plane). @@ -298,18 +301,22 @@ public: // Update state - //- First leapfrog velocity adjust and motion part, required + //- First symplectic velocity adjust and motion part, required // before force calculation. Takes old timestep for variable // timestep cases. - void updatePosition(scalar deltaT, scalar deltaT0); + // Changes to Crank-Nicolson integration for subsequent iterations. + void updatePosition(bool firstIter, scalar deltaT, scalar deltaT0); - //- Second leapfrog velocity adjust part - // required after motion and force calculation + //- Second symplectic velocity adjust part + // required after motion and force calculation. + // Changes to Crank-Nicolson integration for subsequent iterations. void updateAcceleration ( + 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/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.C index 15bcc0df93..39453bb9d7 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("test", false)), rhoInf_(1.0), rhoName_(coeffDict().lookupOrDefault("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("g")) @@ -218,12 +203,51 @@ 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_.updatePosition(firstIter, t.deltaTValue(), t.deltaT0Value()); + + motion_.updateAcceleration + ( + 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(); + + // 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(firstIter, t.deltaTValue(), t.deltaT0Value()); + + motion_.updateAcceleration + ( + 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 4f4e1ca31c..5bf52d69cf 100644 --- a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H +++ b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H @@ -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_; From 7e8a4e1fe88250058fdd85939d2d20cf83e60f75 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Mon, 12 Oct 2015 21:33:56 +0100 Subject: [PATCH 03/16] Update header --- .../sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H b/src/sixDoFRigidBodyMotion/sixDoFRigidBodyMotionSolver/sixDoFRigidBodyMotionSolver.H index 5bf52d69cf..2ad4921aaa 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 From 1dead33a8908c1366d36f852f974433cf208260f Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Tue, 13 Oct 2015 22:28:26 +0100 Subject: [PATCH 04/16] uniformFixedValueFvPatchField: Remove the inconsistent optional "value" read in the construction from dictionary. It is important that the initial value is obtained from the table provided to avoid the user having to evaluate a consistent one or risk the code crashing from a very sudden change in the value. --- .../uniformFixedValue/uniformFixedValueFvPatchField.C | 11 ++--------- 1 file changed, 2 insertions(+), 9 deletions(-) diff --git a/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C index 11c8a88685..e587038890 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::uniformFixedValueFvPatchField fixedValueFvPatchField(p, iF), uniformValue_(DataEntry::New("uniformValue", dict)) { - if (dict.found("value")) - { - fvPatchField::operator==(Field("value", dict, p.size())); - } - else - { - const scalar t = this->db().time().timeOutputValue(); - fvPatchField::operator==(uniformValue_->value(t)); - } + const scalar t = this->db().time().timeOutputValue(); + fvPatchField::operator==(uniformValue_->value(t)); } From d3b8af85474751f153daf40fbda9030558fc88f7 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Tue, 13 Oct 2015 22:31:03 +0100 Subject: [PATCH 05/16] tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection: Add missing relaxation entry --- .../laminar/steamInjection/system/fvSolution | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution b/tutorials/multiphase/reactingTwoPhaseEulerFoam/laminar/steamInjection/system/fvSolution index 6cf90aa61f..2486938d7b 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; From f3d4e512424153ed631129adc552d192514423b3 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Wed, 14 Oct 2015 13:15:17 +0100 Subject: [PATCH 06/16] prghTotalPressureFvPatchScalarField: Total pressure BC for p_rgh Resolves some stability issues with the outlet of multiphase problems. --- src/finiteVolume/Make/files | 1 + .../prghTotalPressureFvPatchScalarField.C | 213 +++++++++++++++ .../prghTotalPressureFvPatchScalarField.H | 243 ++++++++++++++++++ .../laminar/bubbleColumn/0/p_rgh | 17 +- 4 files changed, 467 insertions(+), 7 deletions(-) create mode 100644 src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.C create mode 100644 src/finiteVolume/fields/fvPatchFields/derived/prghTotalPressure/prghTotalPressureFvPatchScalarField.H diff --git a/src/finiteVolume/Make/files b/src/finiteVolume/Make/files index c7100518ab..1b9c242765 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 0000000000..7b8789df44 --- /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 . + +\*---------------------------------------------------------------------------*/ + +#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& iF +) +: + fixedValueFvPatchScalarField(p, iF), + UName_("U"), + phiName_("phi"), + rhoName_("rho"), + p0_(p.size(), 0.0) +{} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedValueFvPatchScalarField(p, iF), + UName_(dict.lookupOrDefault("U", "U")), + phiName_(dict.lookupOrDefault("phi", "phi")), + rhoName_(dict.lookupOrDefault("rhoName", "rho")), + p0_("p0", dict, p.size()) +{ + if (dict.found("value")) + { + fvPatchScalarField::operator= + ( + scalarField("value", dict, p.size()) + ); + } + else + { + fvPatchField::operator=(p0_); + } +} + + +Foam::prghTotalPressureFvPatchScalarField:: +prghTotalPressureFvPatchScalarField +( + const prghTotalPressureFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& 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& 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(ptf); + + p0_.rmap(tiptf.p0_, addr); +} + + +void Foam::prghTotalPressureFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + const scalarField& rhop = + patch().lookupPatchField(rhoName_); + + const scalarField& phip = + patch().lookupPatchField(phiName_); + + const vectorField& Up = + patch().lookupPatchField(UName_); + + const uniformDimensionedVectorField& g = + db().lookupObject("g"); + + const uniformDimensionedScalarField& hRef = + db().lookupObject("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(os, "U", "U", UName_); + writeEntryIfDifferent(os, "phi", "phi", phiName_); + writeEntryIfDifferent(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 0000000000..f809029261 --- /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 . + +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& + ); + + //- Construct from patch, internal field and dictionary + prghTotalPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given + // prghTotalPressureFvPatchScalarField onto a new patch + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct as copy + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new prghTotalPressureFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + prghTotalPressureFvPatchScalarField + ( + const prghTotalPressureFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + 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/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh b/tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumn/0/p_rgh index de5ac9e437..d78578a440 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; } } From 79f20242484ae06ad5a94c2414910c8b6983c9d5 Mon Sep 17 00:00:00 2001 From: Chris Greenshields Date: Thu, 15 Oct 2015 12:10:12 +0100 Subject: [PATCH 07/16] template cases: corrected BirdCarreauCoeffs names (m -> k) --- etc/templates/axisymmetricJet/constant/transportProperties | 2 +- etc/templates/closedVolumeRotating/constant/transportProperties | 2 +- etc/templates/inflowOutflow/constant/transportProperties | 2 +- .../inflowOutflowRotating/constant/transportProperties | 2 +- 4 files changed, 4 insertions(+), 4 deletions(-) diff --git a/etc/templates/axisymmetricJet/constant/transportProperties b/etc/templates/axisymmetricJet/constant/transportProperties index 455ea9b0bf..2a62bb7f6c 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 455ea9b0bf..2a62bb7f6c 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 455ea9b0bf..2a62bb7f6c 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 455ea9b0bf..2a62bb7f6c 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; } From eb1080c933785abc47942851cdcad5f8a43146d6 Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sat, 17 Oct 2015 13:56:34 +0100 Subject: [PATCH 08/16] checkMesh: Provide the number of geometric and solution directions. Resolves bug-report http://www.openfoam.org/mantisbt/view.php?id=1872 --- .../mesh/manipulation/checkMesh/checkGeometry.C | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C b/applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C index 987e26d2e3..68c316f24d 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