I'm trying to recreate an example for Discretization of the Laplace equation from this side: https://folk.ntnu.no/leifh/teaching/tkt4140/._main055.html
to a more general case.
I want to use is it for a general grid with a right-hand-side (like Neumann boundary condition) and to use it for a non-uniformal grid later. But I'm having difficulties, because the solutions of the example and my general non-uniform dicretized grid are not the same.
First Laplace Operator: $\nabla^2 u = 0. \;$
$\frac{\partial^2 u}{\partial x^2} = \frac{u_{i+1,j} - 2u_{i,j} + u_{i-1,j}}{\Delta x^2}, \, \frac{\partial^2 u}{\partial y^2} = \frac{u_{i,j+1} - 2u_{i,j} + u_{i,j-1}}{\Delta y^2}$.
$\Rightarrow-2u_{i,j} (\Delta x^2 + \Delta y^2) + (u_{i+1,j} + u_{i-1,j})\Delta y^2 + (u_{i,j+1} + u_{i,j-1})\Delta x^2 = 0 $.
In this pictur we can see the nodes $T_{ij}, \, i=1,2,3, \, j=1,2,3,4,5 $ (above $T$ is $u$ for temperature) with the boundaries. I used in my case the distance $1$, not $0.25$, like in this example. So $\Delta x = 1, \Delta y = 1$.
I now present $3$ different types of nodes an calculate the rows for the matrix from the coefficients:
$\cdot$ Node $T_{11}$ (corner):
$-2 T_{11} (\Delta y^2 + \Delta x^2) + (T_{21} + 0)\Delta y^2 + (T_{1,2} + 0)\Delta x^2 = 0$
$\Rightarrow -4 \cdot T_{11} + 1 \cdot T_{21}\Delta y^2 + 1 \cdot T_{12}\Delta x^2 = 0$
Which gives me the same first row as in the example.
$\cdot$ Node $T_{12}$ (boundary, not corner):
$-2 T_{12} (\Delta y^2 + \Delta x^2) + (0 + T_{22})\Delta y^2 + (T_{11} + T_{13})\Delta x^2 = 0 $
$\Rightarrow -6 \cdot T_{12} + 1 \cdot T_{22} + 1 \cdot T_{11} + 1 \cdot T_{13} = 0$
This is the second row, but it's not the same, as in the example before. There it is $-4$ instead of $-6$.
$\cdot$ Node $T_{22}$ (not boundary):
$-2 T_{22} (\Delta y^2 + \Delta x^2) + (T_{12} + T_{32})\Delta y^2 + (T_{21} + T_{23})\Delta x^2 = 0$
$\Rightarrow -8 \cdot T_{12} + 1 \cdot T_{12} + 1 \cdot T_{32} + 1 \cdot T_{21} + 1 \cdot T_{23} = 0$
This is the $7$th row in the matrix. And it is again not the same as in the example ($-4$ instead of $-8$).
What I get is the same matrix (with same location of $1$s and $0$s), but with a different diagonal as in the example:
$\text{diag} (A) = (-4, -6, -6, -6, -4, -6, -8, -8, -8, -6, -4, -6, -6, -6, -4)$ for my $A$
$\text{diag} (A) = (-4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4, -4)$ for example $A$.
And $b = (0, 0, 0, 0, -100, 0, 0, 0, 0, -100, 0, 0, 0, 0, -100)$, stays the same.
When I now solve this system $Ax = b$, I get a kind of different solution, which physically makes no sense (the solution is ordered like the nodes in the picture) (I solved this system in python with np.linalg.solve(A, b)):
$\underbrace{\begin{bmatrix}
T_{15} & T_{25} & T_{35} \\
T_{14} & T_{24} & T_{34} \\
T_{13} & T_{23} & T_{33} \\
T_{12} & T_{22} & T_{32} \\
T_{11} & T_{21} & T_{31}
\end{bmatrix}}_{T} :
\underbrace{\begin{bmatrix}
33.915 & 28.878 & 33.915 \\
6.784 & 5.436 & 6.784 \\
1.35 & 1.043 & 1.35 \\
0.273 & 0.206 & 0.273 \\
0.084 & 0.062 & 0.084
\end{bmatrix}}_{\text{my solution}} \;
\underbrace{\begin{bmatrix}
43.193 & 53.154 & 43.193 \\
19.62 & 26.228 & 19.62 \\
9.057 & 12.518 & 9.057 \\
4.092 & 5.731 & 4.092 \\
1.578 & 2.222 & 1.578
\end{bmatrix}}_{\text{example solution}}$
The values are a bit off, but i think that's not a big problem.
"My solution" makes no sense because as we can see for example in the top line of the solution ($33.915, 28.878, 33.915$) doesn't make sense for the temperatures. The middle one should be the highest because of the boundarys (see picture), just like in the example.
What am I missing or doing wrong here?
Thank's for any help or advise.