ENH: CHT solvers - updated equation sources

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/chtMultiRegionSimpleFoam.C

@@ -36,7 +36,6 @@ Description
 #include "regionProperties.H"
 #include "solidThermo.H"
 #include "radiationModel.H"
-#include "IOporosityModelList.H"
 #include "IObasicSourceList.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/fluid/EEqn.H

@@ -16,6 +16,9 @@
     );
 
     EEqn.relax();
+
+    sources.constrain(EEqn);
+
     EEqn.solve();
 
     thermo.correct();

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/fluid/UEqn.H

@@ -3,13 +3,13 @@
     (
         fvm::div(phi, U)
       + turb.divDevRhoReff(U)
-      + sources(rho, U)
+     ==
+        sources(rho, U)
     );
 
     UEqn().relax();
 
-    // Add porous resistance
-    pZones.addResistance(UEqn());
+    sources.constrain(UEqn());
 
     solve
     (

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/fluid/createFluidFields.H

@@ -18,7 +18,6 @@
     PtrList<dimensionedScalar> rhoMin(fluidRegions.size());
 
     PtrList<IObasicSourceList> heatSources(fluidRegions.size());
-    PtrList<IOporosityModelList> porousZonesFluid(fluidRegions.size());
 
     // Populate fluid field pointer lists
     forAll(fluidRegions, i)
@@ -200,13 +199,6 @@
             i,
             new IObasicSourceList(fluidRegions[i])
         );
-
-        Info<< "    Adding porous zones\n" << endl;
-        porousZonesFluid.set
-        (
-            i,
-            new IOporosityModelList(fluidRegions[i])
-        );
     }
 
 

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/fluid/pEqn.H

@@ -23,6 +23,8 @@
 
     phiHbyA += phig;
 
+    sources.relativeFlux(fvc::interpolate(rho), phiHbyA);
+
     dimensionedScalar compressibility = fvc::domainIntegrate(psi);
     bool compressible = (compressibility.value() > SMALL);
 
@@ -54,6 +56,7 @@
             // calculated from the relaxed pressure
             U = HbyA + rAU*fvc::reconstruct((phig - p_rghEqn.flux())/rhorAUf);
             U.correctBoundaryConditions();
+            sources.correct(U);
         }
     }
 

applications/solvers/heatTransfer/chtMultiRegionFoam/chtMultiRegionSimpleFoam/fluid/setRegionFluidFields.H

@@ -14,8 +14,6 @@
 
     IObasicSourceList& sources = heatSources[i];
 
-    const IOporosityModelList& pZones = porousZonesFluid[i];
-
     const dimensionedScalar initialMass
     (
         "initialMass",

applications/solvers/heatTransfer/chtMultiRegionFoam/fluid/pEqn.H

@@ -24,6 +24,8 @@
       + phig
     );
 
+    sources.relativeFlux(fvc::interpolate(rho), phiHbyA);
+
     {
         fvScalarMatrix p_rghDDtEqn
         (
@@ -64,6 +66,7 @@
                 U = HbyA
                   + rAU*fvc::reconstruct((phig + p_rghEqn.flux())/rhorAUf);
                 U.correctBoundaryConditions();
+                sources.correct(U);
                 K = 0.5*magSqr(U);
             }
         }

applications/solvers/heatTransfer/chtMultiRegionFoam/solid/solveSolid.H

@@ -10,9 +10,14 @@ if (finalIter)
         (
             fvm::ddt(betav*rho, h)
           - fvm::laplacian(betav*alpha, h, "laplacian(alpha,h)")
-          + sources(rho, h)
+          ==
+            sources(rho, h)
         );
+
         hEqn().relax();
+
+        sources.constrain(hEqn());
+
         hEqn().solve(mesh.solver(h.select(finalIter)));
     }
 }