From 9c8e2fa709c1b33a1cd26e4c087b6bf08ec49d31 Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Thu, 26 Jan 2012 16:23:34 +0000
Subject: [PATCH] ENH: Updated use of field sources in simple/pimple solvers

---
 .../solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H | 4 +++-
 .../solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H    | 3 ++-
 applications/solvers/incompressible/pimpleFoam/UEqn.H         | 4 +---
 applications/solvers/incompressible/pimpleFoam/pEqn.H         | 3 ++-
 .../solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H    | 4 +++-
 .../solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H    | 3 ++-
 .../solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H    | 4 +---
 .../solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H    | 3 ++-
 .../solvers/incompressible/simpleFoam/MRFSimpleFoam/pEqn.H    | 1 +
 .../solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H    | 1 +
 applications/solvers/incompressible/simpleFoam/pEqn.H         | 1 +
 applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H  | 3 +--
 applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H  | 3 ++-
 applications/solvers/lagrangian/coalChemistryFoam/UEqn.H      | 3 +--
 applications/solvers/lagrangian/coalChemistryFoam/pEqn.H      | 4 +++-
 .../lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H  | 3 +--
 .../lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H  | 4 +++-
 17 files changed, 30 insertions(+), 21 deletions(-)

diff --git a/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H b/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H
index a86b142276b..8b496a6aed5 100644
--- a/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/UrelEqn.H
@@ -9,4 +9,6 @@
 
     UrelEqn().relax();
 
-    solve(UrelEqn() == -fvc::grad(p));
+    sources.constrain(UrelEqn());
+
+    solve(UrelEqn() == -fvc::grad(p) + sources(Urel));
diff --git a/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H b/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H
index df3b09adc64..72c0d662909 100644
--- a/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/SRFPimpleFoam/pEqn.H
@@ -1,5 +1,5 @@
 volScalarField rAUrel(1.0/UrelEqn().A());
-Urel = rAUrel*UrelEqn().H();
+Urel = rAUrel*(UrelEqn() == sources(Urel))().H();
 
 if (pimple.nCorrPISO() <= 1)
 {
@@ -37,3 +37,4 @@ p.relax();
 // Momentum corrector
 Urel -= rAUrel*fvc::grad(p);
 Urel.correctBoundaryConditions();
+sources.correct(Urel);
diff --git a/applications/solvers/incompressible/pimpleFoam/UEqn.H b/applications/solvers/incompressible/pimpleFoam/UEqn.H
index d7be010a9f9..a6aef51dafb 100644
--- a/applications/solvers/incompressible/pimpleFoam/UEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/UEqn.H
@@ -5,8 +5,6 @@ tmp<fvVectorMatrix> UEqn
     fvm::ddt(U)
   + fvm::div(phi, U)
   + turbulence->divDevReff(U)
-  ==
-    sources(U)
 );
 
 UEqn().relax();
@@ -17,5 +15,5 @@ volScalarField rAU(1.0/UEqn().A());
 
 if (pimple.momentumPredictor())
 {
-    solve(UEqn() == -fvc::grad(p));
+    solve(UEqn() == -fvc::grad(p) + sources(U));
 }
diff --git a/applications/solvers/incompressible/pimpleFoam/pEqn.H b/applications/solvers/incompressible/pimpleFoam/pEqn.H
index 61acd7d6978..55062637c21 100644
--- a/applications/solvers/incompressible/pimpleFoam/pEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/pEqn.H
@@ -1,4 +1,4 @@
-U = rAU*UEqn().H();
+U = rAU*(UEqn() == sources(U))().H();
 
 if (pimple.nCorrPISO() <= 1)
 {
@@ -36,3 +36,4 @@ p.relax();
 
 U -= rAU*fvc::grad(p);
 U.correctBoundaryConditions();
+sources.correct(U);
diff --git a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
index 5a211caa9c2..618816f5fe1 100644
--- a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H
@@ -9,9 +9,11 @@ tmp<fvVectorMatrix> UEqn
 
 UEqn().relax();
 
+sources.constrain(UEqn());
+
 rAU = 1.0/UEqn().A();
 
 if (pimple.momentumPredictor())
 {
-    solve(UEqn() == -fvc::grad(p));
+    solve(UEqn() == -fvc::grad(p) + sources(U));
 }
diff --git a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H
index aeef26a8377..dd90826e7d0 100644
--- a/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H
+++ b/applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/pEqn.H
@@ -1,4 +1,4 @@
-U = rAU*UEqn().H();
+U = rAU*(UEqn() == sources(U))().H();
 
 if (pimple.nCorrPISO() <= 1)
 {
@@ -46,3 +46,4 @@ fvc::makeRelative(phi, U);
 
 U -= rAU*fvc::grad(p);
 U.correctBoundaryConditions();
+sources.correct(U);
diff --git a/applications/solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H b/applications/solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H
index c869df39a23..74e84d9bfc9 100644
--- a/applications/solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H
+++ b/applications/solvers/incompressible/potentialFreeSurfaceFoam/UEqn.H
@@ -3,8 +3,6 @@ tmp<fvVectorMatrix> UEqn
     fvm::ddt(U)
   + fvm::div(phi, U)
   + turbulence->divDevReff(U)
-  ==
-    sources(U)
 );
 
 
@@ -14,5 +12,5 @@ sources.constrain(UEqn());
 
 if (pimple.momentumPredictor())
 {
-    solve(UEqn() == -fvc::grad(p_gh));
+    solve(UEqn() == -fvc::grad(p_gh) + sources(U));
 }
diff --git a/applications/solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H b/applications/solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H
index a5adeb9e914..d7cf6cd5b17 100644
--- a/applications/solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H
+++ b/applications/solvers/incompressible/potentialFreeSurfaceFoam/pEqn.H
@@ -1,7 +1,7 @@
 volScalarField rAU(1.0/UEqn().A());
 surfaceScalarField rAUf(rAU.name() + 'f', fvc::interpolate(rAU));
 
-U = rAU*UEqn().H();
+U = rAU*(UEqn() == sources(U))().H();
 
 if (pimple.nCorrPISO() <= 1)
 {
@@ -41,3 +41,4 @@ p = p_gh + (g & (mesh.C() + zeta - refLevel));
 
 U -= rAU*fvc::grad(p_gh);
 U.correctBoundaryConditions();
+sources.correct(U);
diff --git a/applications/solvers/incompressible/simpleFoam/MRFSimpleFoam/pEqn.H b/applications/solvers/incompressible/simpleFoam/MRFSimpleFoam/pEqn.H
index 2d2540d3038..93d6537b6a9 100644
--- a/applications/solvers/incompressible/simpleFoam/MRFSimpleFoam/pEqn.H
+++ b/applications/solvers/incompressible/simpleFoam/MRFSimpleFoam/pEqn.H
@@ -34,4 +34,5 @@
     // Momentum corrector
     U -= rAU*fvc::grad(p);
     U.correctBoundaryConditions();
+    sources.correct(U);
 }
diff --git a/applications/solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H b/applications/solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H
index 54eab9cba7b..b9a9079b92c 100644
--- a/applications/solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H
+++ b/applications/solvers/incompressible/simpleFoam/SRFSimpleFoam/pEqn.H
@@ -34,4 +34,5 @@
     // Momentum corrector
     Urel -= rAUrel*fvc::grad(p);
     Urel.correctBoundaryConditions();
+    sources.correct(Urel);
 }
diff --git a/applications/solvers/incompressible/simpleFoam/pEqn.H b/applications/solvers/incompressible/simpleFoam/pEqn.H
index 966dedb6168..78dd40500ba 100644
--- a/applications/solvers/incompressible/simpleFoam/pEqn.H
+++ b/applications/solvers/incompressible/simpleFoam/pEqn.H
@@ -34,4 +34,5 @@
     // Momentum corrector
     U -= rAU*fvc::grad(p);
     U.correctBoundaryConditions();
+    sources.correct(U);
 }
diff --git a/applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H b/applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H
index 70f75201906..b381b2901d3 100644
--- a/applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H
@@ -7,7 +7,6 @@
      ==
         rho.dimensionedInternalField()*g
       + parcels.SU(U)
-      + sources(rho, U)
     );
 
     sources.constrain(UEqn);
@@ -16,5 +15,5 @@
 
     if (pimple.momentumPredictor())
     {
-        solve(UEqn == -fvc::grad(p));
+        solve(UEqn == -fvc::grad(p) + sources(rho, U));
     }
diff --git a/applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H b/applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H
index 141faa30bf2..8a77a2656e4 100644
--- a/applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H
@@ -6,7 +6,7 @@
     thermo.rho() -= psi*p;
 
     volScalarField rAU(1.0/UEqn.A());
-    U = rAU*UEqn.H();
+    U = rAU*(UEqn == sources(rho, U))().H();
 
     if (pZones.size() > 0)
     {
@@ -60,6 +60,7 @@
 
     U -= rAU*fvc::grad(p);
     U.correctBoundaryConditions();
+    sources.correct(U);
 
     rho = thermo.rho();
     rho = max(rho, rhoMin);
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
index fe7b047d9a2..24f56e83d23 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
@@ -7,7 +7,6 @@
         rho.dimensionedInternalField()*g
       + coalParcels.SU(U)
       + limestoneParcels.SU(U)
-      + sources(rho, U)
     );
 
     UEqn.relax();
@@ -16,6 +15,6 @@
 
     if (pimple.momentumPredictor())
     {
-        solve(UEqn == -fvc::grad(p));
+        solve(UEqn == -fvc::grad(p) + sources(rho, U));
         K = 0.5*magSqr(U);
     }
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
index 56435d31ce1..cb8553bbfbb 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
@@ -1,7 +1,7 @@
 rho = thermo.rho();
 
 volScalarField rAU(1.0/UEqn.A());
-U = rAU*UEqn.H();
+U = rAU*(UEqn == sources(rho, U))().H();
 
 if (pimple.transonic())
 {
@@ -74,6 +74,8 @@ else
 
 U -= rAU*fvc::grad(p);
 U.correctBoundaryConditions();
+sources.correct(U);
+
 K = 0.5*magSqr(U);
 
 dpdt = fvc::ddt(p);
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
index 28f6f81d9f0..ce25309532b 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
@@ -7,7 +7,6 @@
      ==
         rho.dimensionedInternalField()*g
       + parcels.SU(U)
-      + sources(rho, U)
     );
 
     UEqn.relax();
@@ -18,7 +17,7 @@
 
     if (pimple.momentumPredictor())
     {
-        solve(UEqn == -fvc::grad(p));
+        solve(UEqn == -fvc::grad(p) + sources(rho, U));
         K = 0.5*magSqr(U);
     }
 
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
index 080faeccedd..6ed7affe77c 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
@@ -6,7 +6,7 @@
     thermo.rho() -= psi*p;
 
     volScalarField rAU(1.0/UEqn.A());
-    U = rAU*UEqn.H();
+    U = rAU*(UEqn == sources(rho, U))().H();
 
     if (pZones.size() > 0)
     {
@@ -59,6 +59,8 @@
 
     U -= rAU*fvc::grad(p);
     U.correctBoundaryConditions();
+    sources.correct(U);
+
     K = 0.5*magSqr(U);
 
     dpdt = fvc::ddt(p);
-- 
GitLab