simpleFoam/SRFSimpleFoam: Added support for SIMPLEC

SIMPLEC (SIMPLE-consistent) is selected by setting "consistent" option true/yes:

SIMPLE
{
    nNonOrthogonalCorrectors 0;
    consistent yes;
}

which relaxes the pressure in a "consistent" manner and additional
relaxation of the pressure is not generally necessary.  In addition
convergence of the p-U system is better and reliable with less
aggressive relaxation of the momentum equation, e.g. for the motorbike
tutorial:

relaxationFactors
{
    equations
    {
        U               0.9;
        k               0.7;
        omega           0.7;
    }
}

The cost per iteration is marginally higher but the convergence rate is
better so the number of iterations can be reduced.

The SIMPLEC algorithm also provides benefit for cases with large
body-forces, e.g. SRF, see tutorials/incompressible/SRFSimpleFoam/mixer
and feature request http://www.openfoam.org/mantisbt/view.php?id=1714

applications/solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H

@@ -2,17 +2,28 @@
     volScalarField rAUrel(1.0/UrelEqn().A());
     volVectorField HbyA("HbyA", Urel);
     HbyA = rAUrel*UrelEqn().H();
-    UrelEqn.clear();
 
     surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(HbyA) & mesh.Sf());
     adjustPhi(phiHbyA, Urel, p);
 
+    tmp<volScalarField> rAtUrel(rAUrel);
+
+    if (simple.consistent())
+    {
+        rAtUrel = 1.0/(1.0/rAUrel - UrelEqn().H1());
+        phiHbyA +=
+            fvc::interpolate(rAtUrel() - rAUrel)*fvc::snGrad(p)*mesh.magSf();
+        HbyA -= (rAUrel - rAtUrel())*fvc::grad(p);
+    }
+
+    UrelEqn.clear();
+
     // Non-orthogonal pressure corrector loop
     while (simple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
-            fvm::laplacian(rAUrel, p) == fvc::div(phiHbyA)
+            fvm::laplacian(rAtUrel(), p) == fvc::div(phiHbyA)
         );
 
         pEqn.setReference(pRefCell, pRefValue);
@@ -31,7 +42,7 @@
     p.relax();
 
     // Momentum corrector
-    Urel = HbyA - rAUrel*fvc::grad(p);
+    Urel = HbyA - rAtUrel()*fvc::grad(p);
     Urel.correctBoundaryConditions();
     fvOptions.correct(Urel);
 }

applications/solvers/incompressible/simpleFoam/pEqn.H

@@ -2,20 +2,29 @@
     volScalarField rAU(1.0/UEqn().A());
     volVectorField HbyA("HbyA", U);
     HbyA = rAU*UEqn().H();
-    UEqn.clear();
 
     surfaceScalarField phiHbyA("phiHbyA", fvc::interpolate(HbyA) & mesh.Sf());
-
     fvOptions.makeRelative(phiHbyA);
-
     adjustPhi(phiHbyA, U, p);
 
+    tmp<volScalarField> rAtU(rAU);
+
+    if (simple.consistent())
+    {
+        rAtU = 1.0/(1.0/rAU - UEqn().H1());
+        phiHbyA +=
+            fvc::interpolate(rAtU() - rAU)*fvc::snGrad(p)*mesh.magSf();
+        HbyA -= (rAU - rAtU())*fvc::grad(p);
+    }
+
+    UEqn.clear();
+
     // Non-orthogonal pressure corrector loop
     while (simple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
-            fvm::laplacian(rAU, p) == fvc::div(phiHbyA)
+            fvm::laplacian(rAtU(), p) == fvc::div(phiHbyA)
         );
 
         pEqn.setReference(pRefCell, pRefValue);
@@ -34,7 +43,7 @@
     p.relax();
 
     // Momentum corrector
-    U = HbyA - rAU*fvc::grad(p);
+    U = HbyA - rAtU()*fvc::grad(p);
     U.correctBoundaryConditions();
     fvOptions.correct(U);
 }

src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -37,16 +37,21 @@ namespace Foam
 
 void Foam::solutionControl::read(const bool absTolOnly)
 {
-    const dictionary& solnDict = this->dict();
+    const dictionary& solutionDict = this->dict();
 
     // Read solution controls
     nNonOrthCorr_ =
-        solnDict.lookupOrDefault<label>("nNonOrthogonalCorrectors", 0);
+        solutionDict.lookupOrDefault<label>("nNonOrthogonalCorrectors", 0);
-    momentumPredictor_ = solnDict.lookupOrDefault("momentumPredictor", true);
+    momentumPredictor_ =
-    transonic_ = solnDict.lookupOrDefault("transonic", false);
+        solutionDict.lookupOrDefault("momentumPredictor", true);
+    transonic_ = solutionDict.lookupOrDefault("transonic", false);
+    consistent_ = solutionDict.lookupOrDefault("consistent", false);
 
     // Read residual information
-    const dictionary residualDict(solnDict.subOrEmptyDict("residualControl"));
+    const dictionary residualDict
+    (
+        solutionDict.subOrEmptyDict("residualControl")
+    );
 
     DynamicList<fieldData> data(residualControl_);
 
@@ -177,6 +182,7 @@ Foam::solutionControl::solutionControl(fvMesh& mesh, const word& algorithmName)
     nNonOrthCorr_(0),
     momentumPredictor_(true),
     transonic_(false),
+    consistent_(false),
     corr_(0),
     corrNonOrtho_(0)
 {}

src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControl.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -83,6 +83,10 @@ protected:
             //- Flag to indicate to solve using transonic algorithm
             bool transonic_;
 
+            //- Flag to indicate to relax pressure using the
+            //  "consistent" approach of SIMPLEC
+            bool consistent_;
+
 
         // Evolution
 
@@ -169,6 +173,10 @@ public:
             //- Flag to indicate to solve using transonic algorithm
             inline bool transonic() const;
 
+            //- Flag to indicate to relax pressure using the
+            //  "consistent" approach of SIMPLEC
+            inline bool consistent() const;
+
 
         // Evolution
 

src/finiteVolume/cfdTools/general/solutionControl/solutionControl/solutionControlI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -67,6 +67,12 @@ inline bool Foam::solutionControl::transonic() const
 }
 
 
+inline bool Foam::solutionControl::consistent() const
+{
+    return consistent_;
+}
+
+
 inline bool Foam::solutionControl::correctNonOrthogonal()
 {
     corrNonOrtho_++;

tutorials/incompressible/SRFSimpleFoam/mixer/system/fvSolution

@@ -37,22 +37,14 @@ solvers
 SIMPLE
 {
     nNonOrthogonalCorrectors 0;
+    consistent          yes;
 }
 
 relaxationFactors
 {
-    fields
-    {
-        p               0.3;
-    }
     equations
     {
-        Urel            0.7;
+        ".*"            0.9;
-        k               0.7;
-        epsilon         0.7;
-        omega           0.7;
-        R               0.7;
-        nuTilda         0.7;
     }
 }
 

tutorials/incompressible/SRFSimpleFoam/mixer/turbulenceProperties

@@ -1,30 +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;
-    location    "constant";
-    object      turbulenceProperties;
-}
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-simulationType RAS;
-
-RAS
-{
-    RASModel        kOmegaSST;
-
-    turbulence      on;
-
-    printCoeffs     on;
-}
-
-
-// ************************************************************************* //

tutorials/incompressible/simpleFoam/motorBike/system/fvSolution

@@ -66,6 +66,7 @@ solvers
 SIMPLE
 {
     nNonOrthogonalCorrectors 0;
+    consistent yes;
 }
 
 potentialFlow
@@ -75,13 +76,9 @@ potentialFlow
 
 relaxationFactors
 {
-    fields
-    {
-        p               0.3;
-    }
     equations
     {
-        U               0.7;
+        U               0.9;
         k               0.7;
         omega           0.7;
     }