Example for a Poisson equation in axisymmetric cylindrical coordinates

[Jul3k]

I am trying to model an axis symmetric poisson equation as a two dimensional problem (the magnetic field of a solenoid). The weak formulation of the should be: $\int_\Omega \nabla u \cdot \nabla v \, r \, dx$
I specify the bilinear form by the following code based on a Fenics implementation

@BilinearForm
def claplace(u, v, w):
    r = abs(w.x[1])    
    return dot(grad(u),grad(v)) *r

I am setting the magnetic vector potential to zero on the boundary faces:

boundary = V.get_dofs()
XD = V.zeros()
XD[boundary] = 0

I am unsure if this is the correct approach to implement a Poisson equation in axisymmetric cylindrical coordinates. I am encountering an issue where the solution along the symmetry axis does not appear to be smooth. Furthermore, the results I am obtaining do not align with the analytical solution for an ideal solenoid.

[kinnala]

I think you should have Neumann condition at the axisymmetric boundary?

Something like:

from skfem import *
from skfem.helpers import *

@BilinearForm
def claplace(u, v, w):
    r = abs(w.x[1])    
    return dot(grad(u),grad(v)) *r

@LinearForm
def load(v, w):
    return 1. * v

m = MeshTri().refined(4)
basis = Basis(m, ElementTriP2())

A = asm(claplace, basis)
b = asm(load, basis)

y = solve(*condense(A, b, D=basis.get_dofs({'left', 'right', 'top'})))

basis.plot(y, shading='gouraud').show()

[Jul3k]

For the solenoid a positive current density on one side of the symmetry axis should get canceled by the negative current density on the opposite side. That's why I set the magnetic vector potential to zero on the symmetry axis. I am starting to believe that the problem is the weak formulation itself. In two dimensions it can be treated like a poisson equation but this seems not to be valid for the rotational symmetry.