From cfc6c89487876aba3c82387d2ee5a7b0f050bd83 Mon Sep 17 00:00:00 2001
From: Kutalmis Bercin <kutalmis.bercin@esi-group.com>
Date: Thu, 23 Feb 2023 13:42:31 +0000
Subject: [PATCH] WIP: electricPotential: add space charge density hook

---
 .../electricPotential/electricPotential.C     | 19 ++++++++++++++++++-
 .../electricPotential/electricPotential.H     | 13 ++++++++++++-
 2 files changed, 30 insertions(+), 2 deletions(-)

diff --git a/src/functionObjects/solvers/electricPotential/electricPotential.C b/src/functionObjects/solvers/electricPotential/electricPotential.C
index c88bb2686c2..c25569518b6 100644
--- a/src/functionObjects/solvers/electricPotential/electricPotential.C
+++ b/src/functionObjects/solvers/electricPotential/electricPotential.C
@@ -210,8 +210,10 @@ Foam::functionObjects::electricPotential::electricPotential
             IOobject::scopedName(typeName, "V")
         )
     ),
+    erhoName_(),
     nCorr_(1),
-    writeDerivedFields_(false)
+    writeDerivedFields_(false),
+    spaceChargeDensity_(false)
 {
     read(dict);
 
@@ -246,6 +248,11 @@ bool Foam::functionObjects::electricPotential::read(const dictionary& dict)
     dict.readIfPresent("nCorr", nCorr_);
     dict.readIfPresent("writeDerivedFields", writeDerivedFields_);
 
+    if (dict.readIfPresent("spaceChargeDensity", spaceChargeDensity_))
+    {
+        dict.readEntry("erho", erhoName_);
+    }
+
     // If flow is multiphase
     if (!phasesDict_.empty())
     {
@@ -329,6 +336,9 @@ bool Foam::functionObjects::electricPotential::execute()
     tmp<volScalarField> tsigma = this->sigma();
     const volScalarField& sigma = tsigma();
 
+    tmp<volScalarField> tepsilon = this->epsilonm();
+    const auto& epsilon = tepsilon();
+
     volScalarField& eV = getOrReadField(fieldName_);
 
     for (label i = 1; i <= nCorr_; ++i)
@@ -340,6 +350,13 @@ bool Foam::functionObjects::electricPotential::execute()
             fvOptions_(eV)
         );
 
+        if (spaceChargeDensity_)
+        {
+            volScalarField& erho = getOrReadField(erhoName_);
+
+            eVEqn += erho/epsilon;
+        }
+
         eVEqn.relax();
 
         fvOptions_.constrain(eVEqn);
diff --git a/src/functionObjects/solvers/electricPotential/electricPotential.H b/src/functionObjects/solvers/electricPotential/electricPotential.H
index 4f015532deb..e8906684322 100644
--- a/src/functionObjects/solvers/electricPotential/electricPotential.H
+++ b/src/functionObjects/solvers/electricPotential/electricPotential.H
@@ -37,7 +37,7 @@ Description
     The steady-state equation of the charge conservation:
 
     \f[
-        \nabla \cdot \left( \sigma \nabla V \right) = 0
+        \nabla \cdot \left(\sigma \nabla V \right) = 0
     \f]
 
     where
@@ -121,6 +121,8 @@ Usage
         nCorr                 <label>;
         writeDerivedFields    <bool>;
         fieldName             <word>;
+        spaceChargeDensity    <bool>;
+        erhoName              <word>;
 
         // Inherited entries
         ...
@@ -137,6 +139,9 @@ Usage
       nCorr     | Number of corrector iterations     | label | no   | 1
       writeDerivedFields | Flag to write extra fields | bool | no   | false
       fieldName | Name of operand field   | word | no | electricPotential:V
+      spaceChargeDensity | Flag to include space charge density | bool | no <!--
+                -->                                                 | false
+      erhoName  | Name of space charge density field  | word | no   | -
     \endtable
 
     The inherited entries are elaborated in:
@@ -206,6 +211,9 @@ class electricPotential
         //- Name of the operand field
         word fieldName_;
 
+        //- Name of the space charge density field
+        word erhoName_;
+
         //- Number of corrector iterations
         label nCorr_;
 
@@ -213,6 +221,9 @@ class electricPotential
         //- electric field, current density and free-charge density
         bool writeDerivedFields_;
 
+        //- Flag to add the space charge density to the rhs of the potential eq
+        bool spaceChargeDensity_;
+
 
     // Private Member Functions
 
-- 
GitLab