ENH: moved sources to eqn rhs

applications/solvers/incompressible/pimpleFoam/UEqn.H

@@ -5,7 +5,8 @@ tmp<fvVectorMatrix> UEqn
     fvm::ddt(U)
   + fvm::div(phi, U)
   + turbulence->divDevReff(U)
-  + sources(U)
+  ==
+    sources(U)
 );
 
 UEqn().relax();

applications/solvers/incompressible/simpleFoam/MRFSimpleFoam/UEqn.H

@@ -4,7 +4,8 @@
     (
         fvm::div(phi, U)
       + turbulence->divDevReff(U)
-      + sources(U)
+      ==
+        sources(U)
     );
 
     mrfZones.addCoriolis(UEqn());

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

@@ -5,7 +5,8 @@
         fvm::div(phi, Urel)
       + turbulence->divDevReff(Urel)
       + SRF->Su()
-      + sources(Urel)
+     ==
+        sources(Urel)
     );
 
     UrelEqn().relax();

applications/solvers/incompressible/simpleFoam/UEqn.H

@@ -4,7 +4,8 @@
     (
         fvm::div(phi, U)
       + turbulence->divDevReff(U)
-      + sources(U)
+      ==
+        sources(U)
     );
 
     UEqn().relax();

applications/solvers/incompressible/simpleFoam/porousSimpleFoam/UEqn.H

@@ -4,7 +4,8 @@
     (
         fvm::div(phi, U)
       + turbulence->divDevReff(U)
-      + sources(U)
+      ==
+        sources(U)
     );
 
     UEqn().relax();

applications/solvers/lagrangian/LTSReactingParcelFoam/UEqn.H

@@ -4,10 +4,10 @@
         fvm::ddt(rho, U)
       + fvm::div(phi, U)
       + turbulence->divDevRhoReff(U)
-      + sources(rho, U)
      ==
         rho.dimensionedInternalField()*g
       + parcels.SU(U)
+      + sources(rho, U)
     );
 
     sources.constrain(UEqn);

applications/solvers/lagrangian/LTSReactingParcelFoam/YEqn.H

@@ -28,10 +28,10 @@ if (solveSpecies)
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
-              + sources(Yi)
              ==
                 parcels.SYi(i, Yi)
               + combustion->R(Yi)
+              + sources(Yi)
             );
 
             sources.constrain(YEqn);

applications/solvers/lagrangian/LTSReactingParcelFoam/hsEqn.H

@@ -4,12 +4,12 @@
         fvm::ddt(rho, hs)
       + mvConvection->fvmDiv(phi, hs)
       - fvm::laplacian(turbulence->alphaEff(), hs)
-      + sources(rho, hs)
      ==
         DpDt
       + parcels.Sh(hs)
       + radiation->Shs(thermo)
       + combustion->Sh()
+      + sources(rho, hs)
     );
 
     sources.constrain(hsEqn);

applications/solvers/lagrangian/LTSReactingParcelFoam/pEqn.H

@@ -27,9 +27,9 @@
     (
         fvc::ddt(rho) + psi*correction(fvm::ddt(p))
       + fvc::div(phi)
-      + sources(psi, p, rho.name())
      ==
         parcels.Srho()
+      + sources(psi, p, rho.name())
     );
 
     while (pimple.correctNonOrthogonal())

applications/solvers/lagrangian/LTSReactingParcelFoam/rhoEqn.H

@@ -34,9 +34,9 @@ Description
     (
         fvm::ddt(rho)
       + fvc::div(phi)
-      + sources(rho)
       ==
         parcels.Srho(rho)
+      + sources(rho)
     );
 
     sources.constrain(rhoEqn);

applications/solvers/lagrangian/coalChemistryFoam/UEqn.H

@@ -3,11 +3,11 @@
         fvm::ddt(rho, U)
       + fvm::div(phi, U)
       + turbulence->divDevRhoReff(U)
-      + sources(rho, U)
      ==
         rho.dimensionedInternalField()*g
       + coalParcels.SU(U)
       + limestoneParcels.SU(U)
+      + sources(rho, U)
     );
 
     UEqn.relax();

applications/solvers/lagrangian/coalChemistryFoam/YEqn.H

@@ -27,10 +27,10 @@ tmp<fv::convectionScheme<scalar> > mvConvection
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
-              + sources(rho, Yi)
               ==
                 coalParcels.SYi(i, Yi)
               + combustion->R(Yi)
+              + sources(rho, Yi)
             );
 
             YiEqn.relax();

applications/solvers/lagrangian/coalChemistryFoam/hsEqn.H

@@ -4,13 +4,13 @@
         fvm::ddt(rho, hs)
       + mvConvection->fvmDiv(phi, hs)
       - fvm::laplacian(turbulence->alphaEff(), hs)
-      + sources(rho, hs)
      ==
         DpDt
       + combustion->Sh()
       + coalParcels.Sh(hs)
       + limestoneParcels.Sh(hs)
       + radiation->Shs(thermo)
+      + sources(rho, hs)
     );
 
     hsEqn.relax();

applications/solvers/lagrangian/coalChemistryFoam/pEqn.H

@@ -22,9 +22,9 @@ if (pimple.transonic())
             fvm::ddt(psi, p)
           + fvm::div(phid, p)
           - fvm::laplacian(rho*rAU, p)
-          + sources(psi, p, rho.name())
          ==
             coalParcels.Srho()
+          + sources(psi, p, rho.name())
         );
 
         sources.constrain(pEqn, rho.name());
@@ -53,9 +53,9 @@ else
             fvm::ddt(psi, p)
           + fvc::div(phi)
           - fvm::laplacian(rho*rAU, p)
-          + sources(psi, p, rho.name())
          ==
             coalParcels.Srho()
+          + sources(psi, p, rho.name())
         );
 
         sources.constrain(pEqn, rho.name());

applications/solvers/lagrangian/coalChemistryFoam/rhoEqn.H

@@ -34,9 +34,9 @@ Description
     (
         fvm::ddt(rho)
       + fvc::div(phi)
-      + sources(rho)
       ==
         coalParcels.Srho(rho)
+      + sources(rho)
     );
 
     sources.constrain(rhoEqn);

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H

@@ -4,10 +4,10 @@
         fvm::ddt(rho, U)
       + fvm::div(phi, U)
       + turbulence->divDevRhoReff(U)
-      + sources(rho, U)
      ==
         rho.dimensionedInternalField()*g
       + parcels.SU(U)
+      + sources(rho, U)
     );
 
     sources.constrain(UEqn);

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/YEqn.H

@@ -29,10 +29,10 @@ if (solveSpecies)
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
-              + sources(rho, Yi)
               ==
                 parcels.SYi(i, Yi)
               + combustion->R(Yi)
+              + sources(rho, Yi)
             );
 
             sources.constrain(YiEqn);

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/hsEqn.H

@@ -34,12 +34,12 @@
             fvm::ddt(rho, hs)
           + mvConvection->fvmDiv(phi, hs)
           - fvm::laplacian(turbulence->alphaEff(), hs)
-          + sources(rho, hs)
          ==
             pWork()
           + parcels.Sh(hs)
           + radiation->Shs(thermo)
           + combustion->Sh()
+          + sources(rho, hs)
         );
 
         sources.constrain(hsEqn);

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H

@@ -27,9 +27,9 @@
     (
         fvc::ddt(rho) + psi*correction(fvm::ddt(p))
       + fvc::div(phi)
-      + sources(psi, p, rho.name())
      ==
         parcels.Srho()
+      + sources(psi, p, rho.name())
     );
 
 

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/rhoEqn.H

@@ -34,9 +34,9 @@ Description
     (
         fvm::ddt(rho)
       + fvc::div(phi)
-      + sources(rho)
       ==
         parcels.Srho(rho)
+      + sources(rho)
     );
 
     sources.constrain(rhoEqn);

src/finiteVolume/cfdTools/general/fieldSources/basicSource/actuationDiskSource/actuationDiskSourceTemplates.C

@@ -55,7 +55,7 @@ void Foam::actuationDiskSource::addActuationDiskAxialInertialResistance
     }
     forAll(cells, i)
     {
-        Usource[cells[i]] += ((Vcells[cells[i]]/V())*T[i]*E) & U[cells[i]];
+        Usource[cells[i]] -= ((Vcells[cells[i]]/V())*T[i]*E) & U[cells[i]];
     }
 }
 

src/finiteVolume/cfdTools/general/fieldSources/basicSource/basicSource/basicSource.H

@@ -33,6 +33,8 @@ Description
         duration        1000.0;                 // duration
         selectionMode   cellSet;                // cellSet // points //cellZone
 
+    On evaluation, source expects to be added to the rhs of the equation
+
 
 SourceFiles
     basicSource.C

src/finiteVolume/cfdTools/general/fieldSources/basicSource/explicitSetValue/ExplicitSetValue.C

@@ -81,32 +81,9 @@ void Foam::ExplicitSetValue<Type>::setValue
             << ">::setValue for source " << name_ << endl;
     }
 
-    DimensionedField<Type, volMesh> rhs
-    (
-        IOobject
-        (
-            name_ + fieldNames_[fieldI] + "rhs",
-            eqn.psi().mesh().time().timeName(),
-            eqn.psi().mesh(),
-            IOobject::NO_READ,
-            IOobject::NO_WRITE,
-            false
-        ),
-        eqn.psi().mesh(),
-        dimensioned<Type>
-        (
-            "zero",
-            dimless,
-            pTraits<Type>::zero
-        )
-    );
-
     List<Type> values(cells_.size());
 
-    forAll(values, i)
+    UIndirectList<Type>(values, cells_) = injectionRate_[fieldI];
-    {
-        values[i] = injectionRate_[fieldI];
-    }
 
     eqn.setValues(cells_, values);
 }

src/finiteVolume/cfdTools/general/fieldSources/basicSource/explicitSource/ExplicitSource.C

@@ -137,7 +137,7 @@ void Foam::ExplicitSource<Type>::addSup
     const label fieldI
 )
 {
-//    if (debug)
+    if (debug)
     {
         Info<< "ExplicitSource<"<< pTraits<Type>::typeName
             << ">::addSup for source " << name_ << endl;
@@ -163,12 +163,9 @@ void Foam::ExplicitSource<Type>::addSup
         false
     );
 
-    forAll(cells_, i)
+    UIndirectList<Type>(Su, cells_) = injectionRate_[fieldI]/VDash_;
-    {
-        Su[cells_[i]] = injectionRate_[fieldI]/VDash_;
-    }
 
-    eqn -= Su;
+    eqn += Su;
 }
 
 

src/finiteVolume/cfdTools/general/fieldSources/basicSource/radialActuationDiskSource/radialActuationDiskSourceTemplates.C

@@ -79,7 +79,7 @@ addRadialActuationDiskAxialInertialResistance
 
     forAll(cells, i)
     {
-        Usource[cells[i]] += ((Vcells[cells[i]]/V_)*Tr[i]*E) & U[cells[i]];
+        Usource[cells[i]] -= ((Vcells[cells[i]]/V_)*Tr[i]*E) & U[cells[i]];
     }
 
     if (debug)

src/finiteVolume/cfdTools/general/fieldSources/basicSource/rotorDiskSource/rotorDiskSource.C

@@ -399,7 +399,7 @@ Foam::rotorDiskSource::~rotorDiskSource()
 
 void Foam::rotorDiskSource::addSup(fvMatrix<vector>& eqn, const label fieldI)
 {
-    // add source to lhs of eqn
+    // add source to rhs of eqn
 
     const volVectorField& U = eqn.psi();
 
@@ -410,7 +410,7 @@ void Foam::rotorDiskSource::addSup(fvMatrix<vector>& eqn, const label fieldI)
         const volScalarField& rho =
             mesh_.lookupObject<volScalarField>(rhoName_);
 
-        eqn += calculateForces
+        eqn -= calculateForces
             (
                 rho.internalField(),
                 inflowVelocity(U),
@@ -419,7 +419,7 @@ void Foam::rotorDiskSource::addSup(fvMatrix<vector>& eqn, const label fieldI)
     }
     else
     {
-        eqn += calculateForces
+        eqn -= calculateForces
             (
                 oneField(),
                 inflowVelocity(U),