From 3d72ced0d422ff32c681553e1be6365474c59ce2 Mon Sep 17 00:00:00 2001
From: ignis <ignis@postech.ac.kr>
Date: Mon, 14 Nov 2016 03:25:07 +0900
Subject: [PATCH] plasmaReactingFoam write electron velocities and fvm::div
 instead of mvConvection->fvmDiv

---
 .../combustion/plasmaReactingFoam/EEqn.H      |  2 +-
 .../combustion/plasmaReactingFoam/YEqn.H      | 12 ++++----
 .../plasmaReactingFoam/createFields.H         | 28 +++++++++++++++++--
 3 files changed, 33 insertions(+), 9 deletions(-)

diff --git a/applications/solvers/combustion/plasmaReactingFoam/EEqn.H b/applications/solvers/combustion/plasmaReactingFoam/EEqn.H
index 9267c9a9..9262b913 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/EEqn.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/EEqn.H
@@ -3,7 +3,7 @@
 
     fvScalarMatrix EEqn
     (
-        fvm::ddt(rho, he) + mvConvection->fvmDiv(phi, he)
+        fvm::ddt(rho, he) + fvm::div(phi, he, "div(phi,Yi_h)")
       + fvc::ddt(rho, K) + fvc::div(phi, K)
       + (
             he.name() == "e"
diff --git a/applications/solvers/combustion/plasmaReactingFoam/YEqn.H b/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
index 3d515c89..639bb567 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
@@ -132,9 +132,9 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
         if (Y[i].name() == electronSpecie)
         {
-            Udrift = - linearInterpolate
-                ((mue/ng)*E + ((De/ng/Te)*fvc::grad(Te)));
-            ve = (Udrift & mesh.Sf()) + q;
+            Udrift = -(mue/ng)*E;
+            Uthermal = -((De/ng/Te)*fvc::grad(Te));
+            ve = (linearInterpolate(Udrift+Uthermal) & mesh.Sf()) + q;
 
             // Wall electron flux correction
             forAll (wallPatcheIDs, pidx)
@@ -216,10 +216,10 @@ tmp<fv::convectionScheme<scalar> > mvConvection
                 fvm::ddt(rho, Yi)
               +
               ( nCharge != 0
-              ? mvConvection->fvmDiv(phi_drift, Yi)
+              ? fvm::div(phi_drift, Yi, "div(phi,Yi_h)")
               : ( neutrals.contains(Y[i].name())
-                  ? mvConvection->fvmDiv(phi_neutral, Yi)
-                  : mvConvection->fvmDiv(phi, Yi)
+                  ? fvm::div(phi_neutral, Yi, "div(phi,Yi_h)")
+                  : fvm::div(phi, Yi, "div(phi,Yi_h)")
                 )
               )
               // - fvm::laplacian(turbulence->muEff(), Yi)
diff --git a/applications/solvers/combustion/plasmaReactingFoam/createFields.H b/applications/solvers/combustion/plasmaReactingFoam/createFields.H
index fbdee77d..e7ec79a8 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/createFields.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/createFields.H
@@ -330,7 +330,31 @@ bolos.presolve();
 */
 
 Info<< "Calculating face flux field ve\n" << endl;
-surfaceVectorField Udrift ("Udrift", - linearInterpolate(mue*E/ng));
+volVectorField Udrift
+(
+    IOobject
+    (
+        "Udrift",
+        runTime.timeName(),
+        mesh,
+        IOobject::NO_READ,
+        IOobject::AUTO_WRITE
+    ),
+    - (mue*E/ng)
+);
+
+volVectorField Uthermal
+(
+    IOobject
+    (
+        "Uthermal",
+        runTime.timeName(),
+        mesh,
+        IOobject::NO_READ,
+        IOobject::AUTO_WRITE
+    ),
+    - ((De/ng/Te)*fvc::grad(Te))
+);
 
 surfaceScalarField ve
 (
@@ -342,7 +366,7 @@ surfaceScalarField ve
         IOobject::NO_READ,
         IOobject::NO_WRITE
     ),
-    Udrift & mesh.Sf()
+    linearInterpolate(Udrift+Uthermal) & mesh.Sf()
 );