From edb77356aa1e7eca456e0080469c8996b3332f32 Mon Sep 17 00:00:00 2001 From: Vuko Vukcevic Date: Mon, 2 Jan 2017 16:56:52 +0100 Subject: [PATCH] Run-time bugfixes regarding time consistency --- .../pimpleControl/pimpleControl.C | 17 +++++----- .../pimpleControl/pimpleControl.H | 18 +++++----- .../solutionControl/solutionControl.C | 34 +++++++++++-------- .../solutionControl/solutionControl.H | 19 ++++++----- .../CoEulerDdtScheme/CoEulerDdtScheme.C | 2 +- .../EulerDdtScheme/EulerDdtScheme.C | 2 +- .../ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C | 2 +- .../backwardDdtScheme/backwardDdtScheme.C | 4 +-- .../steadyInertialDdtScheme.C | 2 +- 9 files changed, 54 insertions(+), 46 deletions(-) diff --git a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C index 63b60c4b1..52fc7c7a5 100644 --- a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C +++ b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.C @@ -175,16 +175,17 @@ const Foam::dimensionedScalar Foam::pimpleControl::relaxFactor const Foam::volVectorField& U ) const { - return - dimensionedScalar - ( - "alphaU", - dimless, - mesh_.solutionDict().equationRelaxationFactor + scalar urf = 1; + + if (mesh_.solutionDict().relaxEquation(U.select(finalIter()))) + { + urf = mesh_.solutionDict().equationRelaxationFactor ( U.select(finalIter()) - ) - ); + ); + } + + return dimensionedScalar("alphaU", dimless, urf); } diff --git a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.H b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.H index 8e559c4f2..0584c7d66 100644 --- a/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.H +++ b/src/finiteVolume/cfdTools/general/solutionControl/pimpleControl/pimpleControl.H @@ -91,6 +91,15 @@ protected: virtual bool criteriaSatisfied(); + // Time and under-relaxation consistency + + //- Get relaxation factor for velocity field + virtual const dimensionedScalar relaxFactor + ( + const volVectorField& U + ) const; + + public: // Static Data Members @@ -123,15 +132,6 @@ public: inline label corrPISO() const; - // Time and under-relaxation consistency - - //- Get relaxation factor for velocity field - virtual const dimensionedScalar relaxFactor - ( - const volVectorField& U - ) const; - - // Solution control //- PIMPLE loop diff --git a/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.C b/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.C index 301aaaa9b..ade2b7dd8 100644 --- a/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.C +++ b/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.C @@ -222,6 +222,22 @@ Foam::scalar Foam::solutionControl::maxResidual } +const Foam::dimensionedScalar Foam::solutionControl::relaxFactor +( + const Foam::volVectorField& U +) const +{ + scalar urf = 1; + + if (mesh_.solutionDict().relaxEquation(U.name())) + { + urf = mesh_.solutionDict().equationRelaxationFactor(U.name()); + } + + return dimensionedScalar("alphaU", dimless, urf); +} + + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // Foam::solutionControl::solutionControl(fvMesh& mesh, const word& algorithmName) @@ -331,6 +347,9 @@ void Foam::solutionControl::calcTimeConsistentFlux // Interpolate original rAU on the faces const surfaceScalarField rAUf = fvc::interpolate(rAU); + // Store previous iteration for the correct handling of under-relaxation + phi.storePrevIter(); + // Calculate the ordinary part of the flux (H/A) phi = (faceU & mesh_.Sf()); @@ -471,19 +490,4 @@ const Foam::volScalarField& Foam::solutionControl::aCoeff() const } -const Foam::dimensionedScalar Foam::solutionControl::relaxFactor -( - const Foam::volVectorField& U -) const -{ - return - dimensionedScalar - ( - "alphaU", - dimless, - mesh_.solutionDict().equationRelaxationFactor(U.name()) - ); -} - - // ************************************************************************* // diff --git a/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.H b/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.H index ce4878afa..e1cf6ddce 100644 --- a/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.H +++ b/src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.H @@ -145,6 +145,17 @@ protected: ) const; + // Time and under-relaxation consistency + + //- Get relaxation factor for velocity field. Overriden in + // pimpleControl since different relaxation factor may be used for + // final iteration. + virtual const dimensionedScalar relaxFactor + ( + const volVectorField& U + ) const; + + private: // Private data @@ -251,14 +262,6 @@ public: //- Const access to aCoeff (needed for pressure equation) const volScalarField& aCoeff() const; - //- Get relaxation factor for velocity field. Overriden in - // pimpleControl since different relaxation factor may be used for - // final iteration. - virtual const dimensionedScalar relaxFactor - ( - const volVectorField& U - ) const; - // Evolution diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C index 91ee1e522..737a6e4bd 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C @@ -660,7 +660,7 @@ CoEulerDdtScheme::fvcDdtConsistentPhiCorr const surfaceScalarField& rAUf ) const { - return (mesh().Sf() & faceU.oldTime())/rAUf*CofrDeltaT(); + return (mesh().Sf() & faceU.oldTime())*rAUf*CofrDeltaT(); } diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C index d39884280..3565a4799 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C @@ -522,7 +522,7 @@ EulerDdtScheme::fvcDdtConsistentPhiCorr const surfaceScalarField& rAUf ) const { - return (mesh().Sf() & faceU.oldTime())/rAUf/mesh().time().deltaT(); + return (mesh().Sf() & faceU.oldTime())*rAUf/mesh().time().deltaT(); } diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C index c6b205992..a9df6eefb 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C @@ -665,7 +665,7 @@ SLTSDdtScheme::fvcDdtConsistentPhiCorr const surfaceScalarField& rAUf ) const { - return (mesh().Sf() & faceU.oldTime())/rAUf*fvc::interpolate(SLrDeltaT()); + return (mesh().Sf() & faceU.oldTime())*rAUf*fvc::interpolate(SLrDeltaT()); } diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C index f77c6ac4e..b7c6e0137 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C @@ -732,8 +732,8 @@ backwardDdtScheme::fvcDdtConsistentPhiCorr const scalar rDeltaT = 1.0/deltaT; - const dimensionedScalar beta("beta", dimTime, coefft0*rDeltaT); - const dimensionedScalar gamma("gamma", dimTime, -coefft00*rDeltaT); + const dimensionedScalar beta("beta", dimless/dimTime, coefft0*rDeltaT); + const dimensionedScalar gamma("gamma", dimless/dimTime, -coefft00*rDeltaT); return rAUf* diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C index 562db3d6e..c5bad65b4 100644 --- a/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C +++ b/src/finiteVolume/finiteVolume/ddtSchemes/steadyInertialDdtScheme/steadyInertialDdtScheme.C @@ -666,7 +666,7 @@ steadyInertialDdtScheme::fvcDdtConsistentPhiCorr const surfaceScalarField& rAUf ) const { - return (mesh().Sf() & faceU.oldTime())/rAUf*CofrDeltaT(); + return (mesh().Sf() & faceU.oldTime())*rAUf*CofrDeltaT(); }