# On the Lagged Diffusivity Method for the Solution of Nonlinear Finite Difference Systems
## 2 Discretization  
### Space Discretization  
For simplicity, let us assume \(\Omega\) to be a 2D bounded rectangular domain. Let us superimpose a grid \(\overline{\Omega}_{h}=\Omega_{h}\cup\partial\Omega_{h}\) to \(\overline{\Omega}\). For simplicity, we use uniform mesh size \(h\) along \(x\) and \(y\). Thus, \(\Omega_{h}\) is made of the mesh points \((x_{i},y_{j})\), \(i=1,...,N\), \(j=1,...,M\) with:  
\[x_{i+1}=x_{i}+h,\quad i=0,...,N\qquad\qquad y_{j+1}=y_{j}+h\quad j=0,...,M.\]  
Using this grid, \(u_{ij}(t):\overline{\Omega}_{h}\rightarrow\mathbb{R}\) denotes the grid function which approximates the solution \(u(x_{i},y_{j},t)\) at the mesh points \((x_{i},y_{j})\) of \(\overline{\Omega}_{h}\), \(i=0,...,N+1\), \(j=0,...,M+1\). This grid function satisfies the initial condition for \((x_{i},y_{j})\in\overline{\Omega}_{h}\) and \(t=t_{0}\) and the Dirichlet boundary conditions for \((x_{i},y_{j})\in\partial\Omega_{h}\), \(t>t_{0}\). We can now approximate the derivatives of Equation (1) by finite difference methods. In this regard, several choices are possible. In particular, it is interesting to analyze the space discretizations obtained by central finite differences and by the upwind scheme. Indeed, this choice influences the properties of the matrix of the discretized system, as we will see in the following. Regarding the upwind scheme, we consider, for instance, the case of positive \(\tilde{\sigma}\). If \(\tilde{\sigma}<\mathbf{0}\) the backward difference quotients which approximate the first-order derivatives are simply replaced by forward difference quotients. Thus, denoting the forward finite differences (in \(x\) and in \(y\) respectively) by \(\Delta_{x}\) and \(\Delta_{y}\), the backward finite differences by \(\
abla_{x}\), \(\
abla_{y}\) and the central finite differences by \(\delta_{x}\), \(\delta_{y}\), we have:  
* Discretization by central finite differences (error \(O(h^{2})\)) \[ \begin{split}&\frac{\mathrm{d}u_{i,j}(t)}{\mathrm{d}t}=\Delta_{x}[ \sigma(x_{i},y_{j},u_{i,j}(t))\
abla_{x}u_{i,j}(t)]+\Delta_{y}[\sigma(x_{i},y_{ j},u_{i,j}(t))\
abla_{y}u_{i,j}(t)]-\\ &-\tilde{v}_{1}\delta_{x}u_{i,j}(t)-\tilde{v}_{2}\delta_{y}u_{i, j}(t)-\alpha_{i,j}u_{i,j}(t)-g(x_{i},y_{j},u_{i,j}(t))+s(x_{i},y_{j},t).\end{split} \] (2)
* Discretization by upwind scheme (error \(O(h)\)), case \(\tilde{\sigma}>\mathbf{0}\) \[ \begin{split}&\frac{\mathrm{d}u_{i,j}(t)}{\mathrm{d}t}=\Delta_{x}[ \sigma_{i,j}(u_{i,j}(t))\
abla_{x}u_{i,j}(t)]+\Delta_{y}[\sigma_{i,j}(u_{i,j}(t ))\
abla_{y}u_{i,j}(t)]-\\ &-\tilde{v}_{1}\
abla_{x}u_{i,j}(t)-\tilde{v}_{2}\
abla_{y}u_{i,j}(t)-\alpha_{i,j}u_{i,j}(t)-g_{i,j}(u_{i,j}(t))+s_{i,j}(t).\end{split} \] (3)  
In both cases, we can write the discretization in a more convenient way by defining  
\[ \begin{split}&\text{Central FD}\qquad\text{Upwind}\\ & L_{i,j}=\frac{1}{h^{2}}\sigma(x_{i},y_{j},u_{i,j}(t))\qquad \qquad\qquad\tilde{L}_{i,j}=\frac{\tilde{v}_{1}}{2h}\qquad\qquad\tilde{L}_{i, j}=\frac{\tilde{v}_{1}}{h}\\ & B_{i,j}=\frac{1}{h^{2}}\sigma(x_{i},y_{j},u_{i,j}(t))\qquad \qquad\tilde{B}_{i,j}=\frac{\tilde{v}_{2}}{2h}\qquad\qquad\tilde{B}_{i,j}= \frac{\tilde{v}_{2}}{h}\\ & R_{i,j}=\frac{1}{h^{2}}\sigma(x_{i+1,}y_{j},u_{i+1,j}(t))\qquad \qquad\tilde{R}_{i,j}=-\tilde{L}_{i,j}\qquad\qquad\tilde{R}_{i,j}=0\\ & T_{i,j}=\frac{1}{h^{2}}\sigma(x_{i},y_{j+1},u_{i,j+1}(t))\qquad \qquad\tilde{T}_{i,j}=-\tilde{B}_{i,j}\qquad\qquad\tilde{T}_{i,j}=0\\ & D_{i,j}=L_{i,j}+B_{i,j}+R_{i,j}+T_{i,j}\qquad\qquad\tilde{D}_{i, j}=\alpha(x_{i},y_{j})\quad\tilde{D}_{i,j}=\alpha(x_{i},y_{j})+\tilde{L}_{i,j}+ \tilde{B}_{i,j}\end{split}