On 2/21/09, Barton Willis <willisb at unk.edu> wrote:
> I saw the glitch too. Did you compare the numerical values near 80? There
> is a tiny(!) chance that the glitch is due to gnuplot, not in the data. I'd guess there
> is a greater chance that the stopping criteria for eigens_by_jacobi is subideal.
>
> Also, I think Maxima has alternative numerical eigenvalue code (translated
> from fortran). Maybe you could try that code to see if the glitch vanishes.
I tried the program shown in Michel's message of Feb 20 with different
values for s (5, 10, or 20, I have a slow computer) and the bump is present
for the values I tried. Also, I replaced eigens_by_jacobi by dgeev (from
LAPACK) and the bump is the same.
I think the next thing to consider is the construction of the matrix M.
Maybe Mathematica and Maxima are computing different values for
some radicals or something; just guessing.
FWIW
Robert Dodier
PS. Here is the program after modifications.
showtime:true;
s:20;
numer:true;
hbar: 1/s;
kappa:-sqrt(2)/2;
/* this has to do with quantum oscillators and spin problems */
H[i,j]:= if i=j then -kappa*hbar*(2*s-2*i)
else if j = i-1 then hbar^(3/2) *i*sqrt(2*s+1-i)
else if j = i+1 then hbar^(3/2) * (i+1)*sqrt(2*s-i)
else 0;
M: genmatrix(H,2*s+1,2*s+1)$
/* M is symmetric*/
load (lapack);
t: dgeev (M, false, false)$
eig:part(t,1);
eig: sort(-eig);
u:create_list(i,i,1,2*s+1);
plot2d([discrete,u,eig], [plot_format, gnuplot]);