From 005aa102df4406b6b09e4c1281b5ceedbd898541 Mon Sep 17 00:00:00 2001
From: ignis <ignis@postech.ac.kr>
Date: Tue, 4 Apr 2017 10:14:18 +0900
Subject: [PATCH] electron and ion drift relaxation

---
 applications/solvers/combustion/plasmaReactingFoam/YEqn.H    | 5 +++++
 .../solvers/combustion/plasmaReactingFoam/createFields.H     | 5 +++++
 applications/solvers/combustion/plasmaReactingFoam/neEqn.H   | 4 ++++
 3 files changed, 14 insertions(+)

diff --git a/applications/solvers/combustion/plasmaReactingFoam/YEqn.H b/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
index 9f23be3b..ac1e073a 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/YEqn.H
@@ -119,6 +119,11 @@ tmp<fv::convectionScheme<scalar> > mvConvection
             {
                 phis[i] = phi;
                 phis[i] += fvc::interpolate((rho*Di/T*(eCharge*z/kB))*E) & mesh.Sf();
+
+                if (relaxDrift < 1.0 && relaxDrift > 0.0)
+                {
+                    phis[i] *= relaxDrift;
+                }
             }
 
             if (ions.contains(Y[i].name()))
diff --git a/applications/solvers/combustion/plasmaReactingFoam/createFields.H b/applications/solvers/combustion/plasmaReactingFoam/createFields.H
index d48e25cc..1bfe5605 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/createFields.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/createFields.H
@@ -12,6 +12,11 @@ IOdictionary physicalProperties
     )
 );
 
+scalar relaxDrift
+(
+    physicalProperties.lookupOrDefault("relaxDrift", 1.0)
+);
+
 dimensionedScalar epsilon0
 (
     physicalProperties.lookup("epsilon0")
diff --git a/applications/solvers/combustion/plasmaReactingFoam/neEqn.H b/applications/solvers/combustion/plasmaReactingFoam/neEqn.H
index 5c551986..dd1978f5 100644
--- a/applications/solvers/combustion/plasmaReactingFoam/neEqn.H
+++ b/applications/solvers/combustion/plasmaReactingFoam/neEqn.H
@@ -33,6 +33,10 @@
 
 
     Udrift = -(mue/ng)*E;
+    if (relaxDrift < 1.0 && relaxDrift > 0.0)
+    {
+        Udrift *= relaxDrift;
+    }
     Uthermal = -((De/ng/Te)*fvc::grad(Te));
     ve = (linearInterpolate(Udrift+Uthermal) & mesh.Sf()) + q;