Scheme for solving Laplace's (Poisson's) equation, with Jacobi's method

We assume now that we have an estimate for the unknown functions \( u_{11} \), \( u_{12} \), \( u_{21} \) and \( u_{22} \). We will call this the zeroth value and label it as \( u^{(0)}_{11} \), \( u^{(0)}_{12} \), \( u^{(0)}_{21} \) and \( u^{(0)}_{22} \). We can then set up an iterative scheme where the next solution is defined in terms of the previous one as $$ \begin{align} u^{(1)}_{11} =&\frac{1}{4}(b_1-u^{(0)}_{12} -u^{(0)}_{21}) \nonumber \\ u^{(1)}_{12} =&\frac{1}{4}(b_2-u^{(0)}_{11}-u^{(0)}_{22}) \nonumber \\ u^{(1)}_{21} =&\frac{1}{4}(b_3-u^{(0)}_{11}-u^{(0)}_{22}) \nonumber \\ u^{(1)}_{22}=&\frac{1}{4}(b_4-u^{(0)}_{12}-u^{(0)}_{21}), \nonumber \end{align} $$ where we have defined the vector $$ \begin{equation*} \mathbf{b}= \begin{bmatrix} u_{01}+u_{10}-\tilde{\rho}_{11}\\ u_{13}+u_{02}-\tilde{\rho}_{12}\\ u_{31}+u_{20}-\tilde{\rho}_{21} \\ u_{32}+u_{23}-\tilde{\rho}_{22}\\ \end{bmatrix}. \end{equation*} $$