(*This seems to be incorrect.*)
Needs["NDSolve`FEM`"]
a1 = 0.1; a2 = 0.2; a3 = 0.3; a4 = 0.4; b = 0.7; sigma0 = 37*10^6;
omega = 2*10^3*Pi; mu0 = 4*10^-7*Pi;
sigma = If[a1 <= r <= a2 || a3 <= r <= a4, sigma0, 0];
bm = ToBoundaryMesh[
   "Coordinates" -> {{0}, {a1}, {a2}, {a3}, {a4}, {b}},
   "BoundaryElements" -> {PointElement[{{1}, {2}, {3}, {4}, {5}, {6}}]}];
bm["Wireframe"]
rm = ToElementMesh[bm, "MeshOrder" -> 2, MaxCellMeasure -> 0.0001];
L = psi''[r] + psi'[r]/r - (1/r^2 + I*omega*sigma*mu0)*psi[r];
B = DirichletCondition[psi[r] == 0, True];
{vals, funs} = NDEigensystem[{L,B}, psi[r], {r} [Element] rm, 10];
lambda = b*Sqrt[-vals];
Sort[lambda, Re@#1 < Re@#2 &] // Chop

(*This seems to be resonable.*)
a1 = 0.1; a2 = 0.2; a3 = 0.3; a4 = 0.4; b = 0.7; sigma0 =37*10^6;
omega = 2*10^3*Pi; mu0 = 4*10^-7*Pi;
sigma = If[a1 <= r <= a2 || a3 <= r <= a4, sigma0, 0];
L = psi''[r] + psi'[r]/r - (1/r^2 + I*omega*sigma*mu0)*psi[r];
B = DirichletCondition[psi[r] == 0, True];
{vals, funs} = NDEigensystem[{L,B}, psi[r], {r, 0, b}, 10, Method ->
{"PDEDiscretization"-> {"FiniteElement", "MeshOptions"->
{"MaxCellMeasure"-> 0.0001,"MeshOrder"-> 2}}}];
lambda = b*Sqrt[-vals];
Sort[lambda, Re@#1 < Re@#2 &] // Chop

{7.383877183370188 + 0.031093165122796473 I,
 14.641300240427059 + 0.0630628542925003 I,
 21.5557952622858 + 0.5427490911070447 I,
 21.92517611601835 + 0.09489510794076031 I,
 26.471413185406863 + 0.34211325835858836 I,
 29.216261799663346 + 0.1267126418808621 I,
 36.510328314631025 + 0.15855627905576822 I,
 42.90614817349877 + 1.09274892501598 I,
 43.80593433114688 + 0.19044807206663805 I,
 48.46846289924706 + 0.6276371418796622 I}