From 664da2773b2010ef784abf242ca67ede46d383ee Mon Sep 17 00:00:00 2001
From: sergio <sergio>
Date: Fri, 5 Apr 2019 14:15:17 -0700
Subject: [PATCH] ENH: Correcting order of the compressibleContErr.H in comp
 solvers. Adding pMin,Pmax pressure control to buoyantPimple and
 chtMultiReagion

---
 .../solvers/compressible/rhoPimpleFoam/pEqn.H |  6 +++--
 .../buoyantPimpleFoam/buoyantPimpleFoam.C     |  1 +
 .../buoyantPimpleFoam/createFields.H          |  5 ++++
 .../heatTransfer/buoyantPimpleFoam/pEqn.H     | 11 +++++----
 .../chtMultiRegionFoam/chtMultiRegionFoam.C   |  1 +
 .../fluid/createFluidFields.H                 | 23 +++++++++++++++++++
 .../chtMultiRegionFoam/fluid/pEqn.H           | 18 ++++++---------
 .../fluid/setRegionFluidFields.H              |  5 ++++
 8 files changed, 53 insertions(+), 17 deletions(-)

diff --git a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
index f2f2c5c71c7..8e2c3d836d3 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
@@ -85,8 +85,6 @@ else
     }
 }
 
-#include "rhoEqn.H"
-#include "compressibleContinuityErrs.H"
 
 // Explicitly relax pressure for momentum corrector
 p.relax();
@@ -102,6 +100,10 @@ if (pressureControl.limit(p))
 }
 
 thermo.correctRho(psi*p - psip0, rhoMin, rhoMax) ;
+
+#include "rhoEqn.H"
+#include "compressibleContinuityErrs.H"
+
 rho = thermo.rho();
 
 // Correct rhoUf if the mesh is moving
diff --git a/applications/solvers/heatTransfer/buoyantPimpleFoam/buoyantPimpleFoam.C b/applications/solvers/heatTransfer/buoyantPimpleFoam/buoyantPimpleFoam.C
index 289c447233f..c83f59780a8 100644
--- a/applications/solvers/heatTransfer/buoyantPimpleFoam/buoyantPimpleFoam.C
+++ b/applications/solvers/heatTransfer/buoyantPimpleFoam/buoyantPimpleFoam.C
@@ -44,6 +44,7 @@ Description
 #include "radiationModel.H"
 #include "fvOptions.H"
 #include "pimpleControl.H"
+#include "pressureControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
diff --git a/applications/solvers/heatTransfer/buoyantPimpleFoam/createFields.H b/applications/solvers/heatTransfer/buoyantPimpleFoam/createFields.H
index 82f63c93fbb..22dc129a873 100644
--- a/applications/solvers/heatTransfer/buoyantPimpleFoam/createFields.H
+++ b/applications/solvers/heatTransfer/buoyantPimpleFoam/createFields.H
@@ -35,6 +35,7 @@ volVectorField U
 
 #include "compressibleCreatePhi.H"
 
+pressureControl pressureControl(p, rho, pimple.dict(), false);
 
 Info<< "Creating turbulence model\n" << endl;
 autoPtr<compressible::turbulenceModel> turbulence
@@ -104,3 +105,7 @@ dimensionedScalar initialMass("initialMass", fvc::domainIntegrate(rho));
 #include "createMRF.H"
 #include "createRadiationModel.H"
 #include "createFvOptions.H"
+
+
+const dimensionedScalar rhoMax("rhoMax", dimDensity, GREAT, pimple.dict());
+const dimensionedScalar rhoMin("rhoMin", dimDensity, Zero, pimple.dict());
diff --git a/applications/solvers/heatTransfer/buoyantPimpleFoam/pEqn.H b/applications/solvers/heatTransfer/buoyantPimpleFoam/pEqn.H
index 9a2553a7e98..c36e096dbd5 100644
--- a/applications/solvers/heatTransfer/buoyantPimpleFoam/pEqn.H
+++ b/applications/solvers/heatTransfer/buoyantPimpleFoam/pEqn.H
@@ -72,8 +72,7 @@ while (pimple.correctNonOrthogonal())
 
 p = p_rgh + rho*gh;
 
-#include "rhoEqn.H"
-#include "compressibleContinuityErrs.H"
+pressureControl.limit(p);
 
 if (p_rgh.needReference())
 {
@@ -90,16 +89,20 @@ if (p_rgh.needReference())
     {
         p += (initialMass - fvc::domainIntegrate(psi*p))
             /compressibility;
-        thermo.correctRho(psi*p - psip0);
+        thermo.correctRho(psi*p - psip0, rhoMin, rhoMax);
         rho = thermo.rho();
         p_rgh = p - rho*gh;
+        p_rgh.correctBoundaryConditions();
     }
 }
 else
 {
-    thermo.correctRho(psi*p - psip0);
+    thermo.correctRho(psi*p - psip0, rhoMin, rhoMax);
 }
 
+#include "rhoEqn.H"
+#include "compressibleContinuityErrs.H"
+
 rho = thermo.rho();
 
 if (thermo.dpdt())
diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionFoam.C b/applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionFoam.C
index ae3954a6891..0b38d827759 100644
--- a/applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionFoam.C
+++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionFoam.C
@@ -51,6 +51,7 @@ Description
 #include "fvOptions.H"
 #include "coordinateSystem.H"
 #include "loopControl.H"
+#include "pressureControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/createFluidFields.H b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/createFluidFields.H
index d0092f800dc..1c9db4f887b 100644
--- a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/createFluidFields.H
+++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/createFluidFields.H
@@ -25,6 +25,11 @@ PtrList<fv::options> fluidFvOptions(fluidRegions.size());
 List<label> pRefCellFluid(fluidRegions.size());
 List<scalar> pRefValueFluid(fluidRegions.size());
 
+PtrList<dimensionedScalar> rhoMinFluid(fluidRegions.size());
+PtrList<dimensionedScalar> rhoMaxFluid(fluidRegions.size());
+
+PtrList<pressureControl> pressureControls(fluidRegions.size());
+
 const uniformDimensionedVectorField& g = meshObjects::gravity::New(runTime);
 
 // Populate fluid field pointer lists
@@ -256,6 +261,24 @@ forAll(fluidRegions, i)
         fluidRegions[i].solutionDict().subDict("PIMPLE");
     pimpleDict.readIfPresent("frozenFlow", frozenFlowFluid[i]);
 
+    rhoMaxFluid.set
+    (
+        i,
+        new dimensionedScalar("rhoMax", dimDensity, GREAT, pimpleDict)
+    );
+
+    rhoMinFluid.set
+    (
+        i,
+        new dimensionedScalar("rhoMin", dimDensity, Zero, pimpleDict)
+    );
+
+    pressureControls.set
+    (
+        i,
+        new pressureControl(thermoFluid[i].p(), rhoFluid[i], pimpleDict, false)
+    );
+
     Info<< "    Adding MRF\n" << endl;
     MRFfluid.set
     (
diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/pEqn.H b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/pEqn.H
index 9bc978b2540..206c32eeb2b 100644
--- a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/pEqn.H
+++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/pEqn.H
@@ -81,17 +81,9 @@ constrainPressure(p_rgh, rho, U, phiHbyA, rhorAUf, MRF);
 
     p = p_rgh + rho*gh;
 
-    // Thermodynamic density update
-    //thermo.correctRho(psi*p - psip0);
 }
 
-
-
-// Solve continuity
-#include "rhoEqn.H"
-
-// Update continuity errors
-#include "compressibleContinuityErrors.H"
+pressureControl.limit(p);
 
 // For closed-volume cases adjust the pressure and density levels
 // to obey overall mass continuity
@@ -110,16 +102,20 @@ if (closedVolume)
     {
         p += (initialMass - fvc::domainIntegrate(psi*p))
             /compressibility;
-        thermo.correctRho(psi*p - psip0);
+        thermo.correctRho(psi*p - psip0, rhoMin, rhoMax);
         rho = thermo.rho();
         p_rgh = p - rho*gh;
+        p_rgh.correctBoundaryConditions();
     }
 }
 else
 {
-    thermo.correctRho(psi*p - psip0);
+    thermo.correctRho(psi*p - psip0,  rhoMin, rhoMax);
 }
 
+#include "rhoEqn.H"
+#include "compressibleContinuityErrors.H"
+
 rho = thermo.rho();
 
 // Update pressure time derivative if needed
diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/setRegionFluidFields.H b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/setRegionFluidFields.H
index 42cfe36e8e9..50935d3490e 100644
--- a/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/setRegionFluidFields.H
+++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/setRegionFluidFields.H
@@ -61,3 +61,8 @@
     const label pRefCell = pRefCellFluid[i];
     const scalar pRefValue = pRefValueFluid[i];
 
+    const dimensionedScalar rhoMax = rhoMaxFluid[i];
+    const dimensionedScalar rhoMin = rhoMinFluid[i];
+
+    const pressureControl& pressureControl = pressureControls[i];
+
-- 
GitLab