I want to solve PDE system, which consists of Poisson equation, continuity equations for electron and hole with drift-diffusion equation numerically, by using method called Newton's method.
This method is explained in several literatures, including
https://nanohub.org/resources/1565/download/ddmodel_solution_details_word.pdf
(But I'm not sure I understand this method at all)
$$ \begin{gathered} \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}V} \right) = - q\left( {p - n + N_D^ + - N_A^ - } \right) \hfill \\ \frac{\partial }{{\partial x}}{J_n} = q\left( {\frac{{\partial n}}{{\partial t}} + \frac{{n - {n_0}}}{{{\tau _n}}} - {G_{op}}} \right) \hfill \\ \frac{\partial }{{\partial x}}{J_p} = - q\left( {\frac{{\partial p}}{{\partial t}} + \frac{{p - {p_0}}}{{{\tau _p}}} - {G_{op}}} \right) \hfill \\ {J_n} = - qn{\mu _n}\frac{\partial }{{\partial x}}V + q{D_n}\frac{\partial }{{\partial x}}n \hfill \\ {J_p} = - qp{\mu _p}\frac{\partial }{{\partial x}}V - q{D_p}\frac{\partial }{{\partial x}}p \hfill \\ \end{gathered} $$
Arranging these, I got $$ \begin{gathered} 0 = \frac{\partial }{{\partial x}} \cdot \left( {\varepsilon \frac{\partial }{{\partial x}}V} \right) + q\left( {p - n + N_D^ + - N_A^ - } \right) \hfill \\ 0 = \frac{{\partial n}}{{\partial t}} + {\mu _n}\frac{{\partial n}}{{\partial x}}\frac{{\partial V}}{{\partial x}} + {\mu _n}n\frac{{{\partial ^2}V}}{{\partial {x^2}}} - {D_n}\frac{{{\partial ^2}n}}{{\partial {x^2}}} + \frac{{n - {n_0}}}{{{\tau _n}}} - {G_{op}} \hfill \\ 0 = \frac{{\partial p}}{{\partial t}} - {\mu _p}\frac{{\partial p}}{{\partial x}}\frac{{\partial V}}{{\partial x}} - {\mu _p}p\frac{{{\partial ^2}V}}{{\partial {x^2}}} + {D_p}\frac{{{\partial ^2}p}}{{\partial {x^2}}} + \frac{{p - {p_0}}}{{{\tau _p}}} - {G_{op}} \hfill \\ \end{gathered} $$
Now, what I understand from general Newton's method for solving equation system is following; I have system of equation $$ \begin{aligned} f_1(x,y,z)=0 \\ f_2(x,y,z)=0 \\ f_3(x,y,z)=0 \end{aligned} \qquad\rightarrow\qquad \mathbf{F}(x,y,z)=0 $$ and initial value (starting value) $$ \mathbf{X}_0=(x_0,y_0,z_0) $$ Then I create Jacobian Matrix $$ \mathbf{J}(x,y,z)= \left[ \begin{array}{ccc} \partial_xf_1 & \partial_yf_1 & \partial_zf_1\\ \partial_xf_2 & \partial_yf_2 & \partial_zf_2\\ \partial_xf_3 & \partial_yf_3 & \partial_zf_3 \end{array} \right] $$ The $\mathbf{X}$ that give $\mathbf{F}=0$ is $$ \mathbf{X}_{n+1}=\mathbf{X}_n-\mathbf{J}^{-1}(\mathbf{X}_n)\mathbf{F}(\mathbf{X}_n) \qquad (n=0,1,2...) $$
I got this idea, and Now I want to apply this into the poisson, continuity, drift-diffusion eq. system. First, I have my system of equation, with independent variables $V, n, p$.
$$ \begin{aligned} 0 =& f_1(V,n,p)= \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}V} \right) + q\left( {p - n + N_D^ + - N_A^ - } \right) \\ 0 =& f_2(V,n,p)= \frac{{\partial n}}{{\partial t}} + {\mu _n}\frac{{\partial n}}{{\partial x}}\frac{{\partial V}}{{\partial x}} + {\mu _n}n\frac{{{\partial ^2}V}}{{\partial {x^2}}} - {D_n}\frac{{{\partial ^2}n}}{{\partial {x^2}}} + \frac{{n - {n_0}}}{{{\tau _n}}} - {G_{op}} \\ 0 =& f_3(V,n,p)= \frac{{\partial p}}{{\partial t}} - {\mu _p}\frac{{\partial p}}{{\partial x}}\frac{{\partial V}}{{\partial x}} - {\mu _p}p\frac{{{\partial ^2}V}}{{\partial {x^2}}} + {D_p}\frac{{{\partial ^2}p}}{{\partial {x^2}}} + \frac{{p - {p_0}}}{{{\tau _p}}} - {G_{op}} \end{aligned} $$ I assumed that the initial condition might be $$ V(x)=\left\{ \begin{aligned} 0 \qquad & (x_1<x<0)\\ k_BT\ln\left(\frac{N_AN_D}{n_i^2}\right)\qquad &(0<x<x_2)\\ \end{aligned}\right. \qquad n(x)=\left\{ \begin{aligned} N_D \qquad & (x_1<x<0)\\ \frac{n_i^2}{N_A} \qquad &(0<x<x_2)\\ \end{aligned}\right. \qquad p(x)=\left\{ \begin{aligned} \frac{n_i^2}{N_D} \qquad & (x_1<x<0)\\ N_A \qquad &(0<x<x_2)\\ \end{aligned}\right. $$
I got stuck when I was creating Jacobian $ \mathbf{J}_{11}$ $$ \partial_Vf_1=\frac{\partial }{\partial V}\left[ \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}V} \right) + q\left( {p - n + N_D^ + - N_A^ - } \right) \right] \stackrel{?}{=} \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}} \right) $$ Apparently, I got an Jacobian element which is just operator! Am I proceeding correctly? If I may go further, I have my Jacobian that looks like $$ \mathbf{J}(V,n,p)= \left[ \begin{array}{ccc} \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}} \right) & -q & q \\ \mu_n\frac{\partial^2 n}{\partial x^2}+mu_n n\frac{\partial^2 }{\partial x^2} & \frac{{\partial }}{{\partial t}} + {\mu _n}\frac{{\partial^2 V}}{{\partial x^2}} + {\mu _n}\frac{{{\partial ^2}V}}{{\partial {x^2}}} - {D_n}\frac{{{\partial ^2}}}{{\partial {x^2}}} + \frac{1}{{{\tau _n}}} & 0 \\ -\mu_p\frac{\partial^2 p}{\partial x^2}-mu_p p\frac{\partial^2 }{\partial x^2} & 0 & \frac{{\partial }}{{\partial t}} - {\mu _p}\frac{{\partial^2 V}}{{\partial x^2}} - {\mu _p}\frac{{{\partial ^2}V}}{{\partial {x^2}}} + {D_p}\frac{{{\partial ^2}}}{{\partial {x^2}}} + \frac{1}{{{\tau _p}}} \end{array} \right] $$
$$ = \left[ \begin{array}{ccc} \frac{\partial }{{\partial x}} \left( {\varepsilon \frac{\partial }{{\partial x}}} \right) & -q & q \\ \mu_n\frac{\partial^2 n}{\partial x^2}+\mu_n n\frac{\partial^2 }{\partial x^2} & \frac{{\partial }}{{\partial t}} + 2{\mu _n}\frac{{\partial^2 V}}{{\partial x^2}} - {D_n}\frac{{{\partial ^2}}}{{\partial {x^2}}} + \frac{1}{{{\tau _n}}} & 0 \\ -\mu_p\frac{\partial^2 p}{\partial x^2}-\mu_p p\frac{\partial^2 }{\partial x^2} & 0 & \frac{{\partial }}{{\partial t}} - 2{\mu _p}\frac{{{\partial ^2}V}}{{\partial {x^2}}} + {D_p}\frac{{{\partial ^2}}}{{\partial {x^2}}} + \frac{1}{{{\tau _p}}} \end{array} \right] $$ I don't have any idea how to apply inverse matrix of this "operator" matrix
Is this correct result?
