From 7d9231d3d113c8bfcee297efe1ed0236cd9c23a6 Mon Sep 17 00:00:00 2001
From: "vinugholap@gmail.com" <Vinugitlab$123>
Date: Wed, 17 Feb 2021 11:12:57 +0100
Subject: [PATCH] Correction in mlaw, compute & store em field source

---
 .../materialLaw/mlawLinearElecMagTherMech.cpp       | 13 +++++++++----
 .../materialLaw/mlawLinearElecMagTherMech.h         |  1 +
 dG3D/src/dG3DMaterialLaw.cpp                        |  6 ++++--
 3 files changed, 14 insertions(+), 6 deletions(-)

diff --git a/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.cpp b/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.cpp
index 7050eec59..2a1febb72 100644
--- a/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.cpp
+++ b/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.cpp
@@ -203,6 +203,7 @@ void mlawLinearElecMagTherMech::constitutive(
                                             SVector3 & dmechSourcedB,
                                             SVector3 & sourceVectorField, // Magnetic source vector field (T, grad V, F)
                                             SVector3 & dsourceVectorFielddt,
+                                            double & emFieldSource,
                                             STensor3 & dHdA,
                                             STensor3 & dHdB,
                                             STensor33 & dPdB,
@@ -350,15 +351,17 @@ void mlawLinearElecMagTherMech::constitutive(
 
                     // derivatives of sourceVectorField
                     dsourceVectorFielddB(i, j) += 0.0;
-                    dsourceVectorFielddgradV(i, j) += _l0(i, j);
-                    dsourceVectorFielddgradT(i, j) += _seebeck * _l0(i ,j);
+                    dsourceVectorFielddgradV(i, j) += -_l0(i, j);
+                    dsourceVectorFielddgradT(i, j) += -_seebeck * _l0(i ,j);
 
                     dsourceVectorFielddV(i) += 0.0;
-                    dsourceVectorFielddT(i) += dseebeckdT * _l0(i, j) * gradT(j);
-                    dsourceVectorFielddA(i, j) += _l0(i, j) * invdh;
+                    dsourceVectorFielddT(i) += -dseebeckdT * _l0(i, j) * gradT(j);
+                    dsourceVectorFielddA(i, j) += -_l0(i, j) * invdh;
                 }
             }
             sourceVectorField(i) = -fluxD(i) - js0(i);
+            // w_AV = j_e \cdot e
+            emFieldSource += sourceVectorField(i) * (-gradV(i) - (An(i) - A0(i)) * invdh);
         }
 
         // Need to check contribution to jy from source js0
@@ -396,6 +399,8 @@ void mlawLinearElecMagTherMech::constitutive(
                 }
             }
             sourceVectorField(i) = -fluxD(i) - js0(i);
+            // w_AV = j_e \cdot e
+            emFieldSource += sourceVectorField(i) * (-gradV(i) - (An(i) - A0(i)) * invdh);
         }
 
         // For test with unit cube to check
diff --git a/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.h b/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.h
index 56c1c123a..170006824 100644
--- a/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.h
+++ b/NonLinearSolver/materialLaw/mlawLinearElecMagTherMech.h
@@ -150,6 +150,7 @@ class mlawLinearElecMagTherMech : public mlawLinearElecTherMech
                             SVector3 & dmechSourcedB,
                             SVector3 & sourceVectorField, // Magnetic source vector field (T, grad V, F)
                             SVector3 & dsourceVectorFielddt,
+                            double & emFieldSource,
                             STensor3 & dHdA,
                             STensor3 & dHdB,
                             STensor33 & dPdB,
diff --git a/dG3D/src/dG3DMaterialLaw.cpp b/dG3D/src/dG3DMaterialLaw.cpp
index 9cd840818..926ca4056 100644
--- a/dG3D/src/dG3DMaterialLaw.cpp
+++ b/dG3D/src/dG3DMaterialLaw.cpp
@@ -5724,6 +5724,7 @@ const bool checkfrac, const bool dTangent)
     SVector3 & dmechSourcedB = ipvcur->getRefTodMechanicalSourcedMagneticVectorCurl();
     SVector3 & sourceVectorField = ipvcur->getRefTosourceVectorField(); // Magnetic source vector field (T, grad V, F)
     SVector3 & dsourceVectorFielddt = ipvcur->getRefTodsourceVectorFielddt();
+    double & emFieldSource = ipvcur->getRefToemFieldSource();
     STensor3 & dHdA = ipvcur->getRefTodMagneticFielddMagneticVectorPotential();
     STensor3 & dHdB = ipvcur->getRefTodMagneticFielddMagneticVectorCurl();
     STensor33 & dPdB = ipvcur->getRefTodPdMagneticVectorCurl();
@@ -5753,7 +5754,7 @@ const bool checkfrac, const bool dTangent)
     dl20dv,dk10dv,dk20dv,dl10dF,dl20dF,dk10dF,dk20dF,fluxjy,djydV,djydgradV,djydT,djydgradT,
     djydF,djydA, djydB, djydgradVdT,djydgradVdV,djydgradVdF,djydgradTdT,djydgradTdV,djydgradTdF,stiff,
     A0, An, B, H, js0, dBdT, dBdGradT, dBdV, dBdGradV, dBdF, dBdA, dfluxTdA, dfluxDdA, dfluxTdB, dfluxDdB, dw_TdA, dw_TdB, dmechSourcedB, sourceVectorField,
-    dsourceVectorFielddt, dHdA, dHdB, dPdB, dHdF, dsourceVectorFielddF, mechFieldSource, dmechFieldSourcedB, dmechFieldSourcedF,
+    dsourceVectorFielddt, emFieldSource, dHdA, dHdB, dPdB, dHdF, dsourceVectorFielddF, mechFieldSource, dmechFieldSourcedB, dmechFieldSourcedF,
     dHdT, dHdgradT, dHdV, dHdgradV, dsourceVectorFielddT, dsourceVectorFielddgradT, dsourceVectorFielddV, dsourceVectorFielddgradV, dsourceVectorFielddA,
     dsourceVectorFielddB, dmechFieldSourcedT, dmechFieldSourcedgradT, dmechFieldSourcedV, dmechFieldSourcedgradV);
 
@@ -5970,6 +5971,7 @@ const bool checkfrac, const bool dTangent)
   SVector3 & dmechSourcedB = ipvcur->getRefTodMechanicalSourcedMagneticVectorCurl();
   SVector3 & sourceVectorField = ipvcur->getRefTosourceVectorField(); // Magnetic source vector field (T, grad V, F)
   SVector3 & dsourceVectorFielddt = ipvcur->getRefTodsourceVectorFielddt();
+  double & emFieldSource = ipvcur->getRefToemFieldSource();
   STensor3 & dHdA = ipvcur->getRefTodMagneticFielddMagneticVectorPotential();
   STensor3 & dHdB = ipvcur->getRefTodMagneticFielddMagneticVectorCurl();
   STensor33 & dPdB = ipvcur->getRefTodPdMagneticVectorCurl();
@@ -5999,7 +6001,7 @@ const bool checkfrac, const bool dTangent)
     dl20dv,dk10dv,dk20dv,dl10dF,dl20dF,dk10dF,dk20dF,fluxjy,djydV,djydgradV,djydT,djydgradT,
     djydF, djydA, djydB, djydgradVdT,djydgradVdV,djydgradVdF,djydgradTdT,djydgradTdV,djydgradTdF,stiff,
     A0, An, B, H, js0, dBdT, dBdGradT, dBdV, dBdGradV, dBdF, dBdA, dfluxTdA, dfluxDdA, dfluxTdB, dfluxDdB, dw_TdA, dw_TdB, dmechSourcedB,
-    sourceVectorField, dsourceVectorFielddt, dHdA, dHdB, dPdB, dHdF, dsourceVectorFielddF, mechFieldSource, dmechFieldSourcedB, dmechFieldSourcedF,
+    sourceVectorField, dsourceVectorFielddt, emFieldSource, dHdA, dHdB, dPdB, dHdF, dsourceVectorFielddF, mechFieldSource, dmechFieldSourcedB, dmechFieldSourcedF,
     dHdT, dHdgradT, dHdV, dHdgradV, dsourceVectorFielddT, dsourceVectorFielddgradT, dsourceVectorFielddV, dsourceVectorFielddgradV, dsourceVectorFielddA,
     dsourceVectorFielddB, dmechFieldSourcedT, dmechFieldSourcedgradT, dmechFieldSourcedV, dmechFieldSourcedgradV);
 
-- 
GitLab