I am trying to solve this equation: $$ \frac{\partial^2 A}{\partial r^2} + \frac{1}{r}\frac{\partial A}{\partial r} - \frac{A}{r^2} + \frac{\partial^2 A}{\partial z^2} = -g $$
This is basically the radial (or angular) equation of the vector laplacian in cylindrical coordinates for a system with azimuthal symmetry, that is, $A = A(r, z)$. In summary, poisson equation.
Using a grid where $r\in [0, R]$ and $z\in[-L, L]$, with equal spacing $h$, I can discretize it:
$$ \frac{A_{i,j+1} + A_{i,j-1} - 2A_{i,j}}{h^2} + \frac{1}{r}\frac{A_{i,j+1} - A_{i,j-1}}{2h} - \frac{1}{r^2}A_{ij} + \frac{A_{i+1,j} + A_{i-1,j} - 2A_{i,j}}{h^2} = -g $$
Multiplying everything by $h^2$: $$ A_{i,j+1} + A_{i,j-1} - 2A_{i,j} + \frac{h}{2r}\left(A_{i,j+1} - A_{i,j-1}\right) - \frac{h^2}{r^2}A_{ij} + A_{i+1,j} + A_{i-1,j} - 2A_{i,j} = -h^2g $$
Re-arranging: $$ \left(-4-\frac{h^2}{r^2}\right)A_{i,j} + A_{i,j+1} + A_{i,j-1} + \frac{h}{2r}\left(A_{i,j+1} - A_{i,j-1}\right) + A_{i+1,j} + A_{i-1,j}= -h^2g $$
Therefore: $$ A_{ij} = \frac{1}{4 + \frac{h^2}{r^2}}\left[ A_{i,j+1} + A_{i,j-1} + A_{i+1,j} + A_{i-1,j} + \frac{h}{2r}\left(A_{i,j+1} - A_{i,j-1}\right) + h^2 g\right] $$
I've been using jacobi iterations with neumann conditions with zero derivative at the boundary $z=-L$, $z=L$ and $r=R$. Furthermore, by symmetry, it can be argued that $A(r=0, z) = 0$, therefore, I've been imposing dirichlet boundary $A=0$ at $r=0$.
I have the analytical solution for a physical system that obeys this equation (the magnetic potential of a current ring of radius $a$). In this case, $g = \frac{I}{a^2}\delta(r-a)\delta(z)$. $$A(r,z) = \frac{Ia}{4\pi}\int_0^{2\pi}\frac{\cos(x) dx}{\sqrt{r^2 + z^2 + a^2 - 2ra\cos(x)}}$$
I've been trying to validate it, but so far I am completely stuck. I've been trying to make this work for two days now, and nothing. I am not sure what I am doing wrong, and I am starting to doubt even if the discretization itself is correct.
If you wish, here are the codes I am using:
Any kind of help is appreciated. Thanks.
- Is this numerical scheme correct?
- Is there something wrong with the code? (if you are so kind take a look at it..). Also, both jacobi and seidel iterations seemed to be giving different results.. aaaa!
