# [Maxima] Now sure what to make of this

Richard Hennessy rvh2007 at comcast.net
Sun May 11 12:58:24 CDT 2008

```This eigenvalue of Shrodinger's eq

-hbar^2/(2*m)*diff(f(x),x, 2) + V(x)*f(x) = E * f(x)

E = 8.03423854737058519307573738131977856841318856113874712108445697431664974582b0

is accurate to 72 decimals for hbar = 1, k = 2/3, m = 2 and V(x) = 1/2 * m * k^2 * x^4.  Not bad for off the cuff hacking.

I stayed up all night for many nights and worked on this for many days but I got it and it is right.

FYI.

Rich

------------Original Message------------
From: "Richard Hennessy"<rvh2007 at comcast.net>
To: "Maxima List" <maxima at math.utexas.edu>
Date: Sat, May-10-2008 2:02 PM
Subject: Re: [Maxima] Not sure what to make of this

Sorry,

Never mind.  I will figure it out.  The potential function should read

1/2*m*k^2*x^4

Sorry about that, I forgot the 1/2.

Rich

------------Original Message------------
From: "Richard Hennessy"<rvh2007 at comcast.net>
To: "Maxima List" <maxima at math.utexas.edu>
Date: Fri, May-9-2008 4:37 PM
Subject: [Maxima] Not sure what to make of this

I have the following two functions.  The first is a solution to Shrodinger for the m*k^2*x^4 potential.  The second tries to find an eigenvalue to the first given two initial guess's.  It is not elegant but it seems to work.  I keep getting the same answer for the eigenvalues accurate to 32 number of digits.  The accurate part I don't know for sure.  All I know is I get the same 32 digits all the time.  How accurate is my answer?  The second function prints out the answer two 32 decimal places and it gives that same exact answer for different parameters and number of iterations as long as it is high enough.  Can anyone tell me if this is right, the method seems good but when I change a digit in the answer is has no effect on my 'tests' for correctness.  I suspect the first 20 digits or so may be right, but how to tell is a problem.  Any ideas

ratprint:false;

ratepsilon:10^-100;

terms:1600;

fpp:500;

f(x,y,terms):=
block
([c,p,n,hbar,k,m],
c[-4]:0,c[-3]:0,c[-2]:0,c[-1]:0,c[0]:1,c[1]:0,hbar:1,m:2,k:2/3,
for n: 2 thru terms do (c[n]:(k^2 * m^2 * c[n-6] - 2 * m * c[n-2] * y)/(n * (n - 1) * hbar^2)),
p:c[0],
for n: 1 thru terms do (p:p + c[n]*x^n),
p
)\$

fsolve(py1, py2, xg, iterations):=
block
(
fpprec:500,i:0, ii:0, j:0, px:xg, pyb:bfloat(py1), pyt:bfloat(py2),
while ii < iterations do
block
(
ft:f(px,pyt,terms),
st:block(dfx:(f(px+1/100,pyt,terms)-ft)*100,
dfy:(f(px,pyt+(pyt-pyb)/1000000,terms)-ft)/((pyt-pyb)/1000000), -dfx/dfy),
fb:f(px,pyb,terms),
sb:block(dfx:(f(px+1/100,pyb,terms)-fb)*100,
dfy:(f(px,pyb+(pyt-pyb)/1000000,terms)-fb)/((pyt-pyb)/1000000), -dfx/dfy),
if st=sb then go(endloop),
p2:(pyb*st-pyt*sb)/(st-sb),
pyt:(pyt+p2*999)/1000,
pyb:(pyb+p2*999)/1000,
ii:ii+1,
i:ii,
go (again),
endloop,
i:ii,
ii:10000,
again,
9
),
fpprec:32,
/* return a list of num iterations and high and low value */
[i,bfloat(pyt),bfloat(pyb)]
)\$

With this

fsolve(8, 9, 6, 15);

->  [15,8.0342385473705851930757373813198b0,8.0342385473705851930757373813198b0]

Thanks,

Rich

PS You can reduce the terms to 1200 with no effect on answer.  You can change 6 to 7 with no effect at 1600 terms.

_______________________________________________
Maxima mailing list
Maxima at math.utexas.edu
http://www.math.utexas.edu/mailman/listinfo/maxima

_______________________________________________
Maxima mailing list
Maxima at math.utexas.edu
http://www.math.utexas.edu/mailman/listinfo/maxima

```