I am trying to find the eigensystem of a system of 3 coupled ODEs. Analytically, the system spectrum should have a gap at [-2M,2M] (except for the degenerate state with E=0, strictly at the middle of the gap, which is of no interest). However, NDEigensystem with Finite Elements Method gives unwanted solutions inside the gap. Why is it so? How can I tune the procedure to get rid of the parasitic solutions? The MWE is below. Thank you!
Needs["NDSolve`FEM`"];
M = 0.01;
B = 50;
A = 5;
R = 100;
{vals, funs} = NDEigensystem[({
{-2 B u''[x] - (2 B)/x u'[x] + (2 B/x^2 + 2 M) u[x] -
2 A ( v'[x] + 2/x v[x])},
{A ( u'[x] - 1/x u[x]) - A ( z'[x] + 3/x z[x])},
{2 A ( v'[x] - 2/x v[x]) +
2 B z''[x] + (2 B)/x z'[x] + (-18 B/x^2 - 2 M) z[x]}
}), {u[x], v[x], z[x]}, {x, 0, R}, 1,
Method -> {"SpatialDiscretization" -> {"FiniteElement", \
{"MeshOptions" -> {"MaxCellMeasure" -> 0.001}}},
"Eigensystem" -> {"Arnoldi", "BasisSize" -> 100,
"MaxIterations" -> 10000, "Shift" -> -M}}];
vals
{-0.00974081}
upd: I suspect the parasitic states arise because we are trying to solve the problem within a finite region, mesh being also of finite coarseness. If I take a finer mesh, the spatial region of $x \in (0,d)$ where parasitic solutions do exist shrinks, i.e. $d$ decreases, but never disappears completely. It is vital for me however to get rid of these solution in the spectrum.
upd2. I tried to use the graded mesh,
mesh = ToGradedMesh[{Line[{{0}, {R}}], <|"Alignment" -> "Left",
"MinimalDistance" -> 0.00001, "ElementCount" -> 1000|>}];
But the problem persists.
upd3. The numerical scheme is stabilized by detuning the band dispersion slightly, in addition to the graded mesh. This gives rise to the appearance of v''[x] in the second equation. The unwanted solutions in the gap (half of the gap in this case) disappear.
{u[x], v[x], z[x]}
is the right order? It matters. What are the right-hand sides of your ODEs? Do you really want to use the defaultNeumannValue
boundary conditions, especially sincev
enters only as first order in derivatives? The ODEs are singular atx = 0
. The eigenfunction solution forv
is very noisy. $\endgroup$x = 0
to see that they are well-behaved. You also should perform a series expansion of the ODEs atx = 0
, which will uniquely determine the boundary conditions there. For the ODEs in your questions, the lowest order nonzero terms areu'[0]
andz'''[0]
. $\endgroup$