100 '  Schroedinger equation
110 '  energy dependence of u(0,e)
120 '  y''=z'=(V(x)-e)*y  y'=z,  y:wave function
130 CONSOLE 0,25,0,1:SCREEN 3,0,0,1
140 VIEW (0,0)-(639,399):CLS 3
150 DIM E(39),PH(29),YYD(29),EDG(100)
160 DIM Y(1001)
170 VIEW (130,10)-(620,380),,6
180 WINDOW (0,-20)-(1500,20)
190 DEF FNP(X)=X*X  :REM potential function
200 DEF FNF(Z)=Z
210 DEF FNG(X,Y,E)=(X*X-E)*Y
220 CLS 3 :LINE (500,0)-(1500,0)
230 LOCATE 30,3:PRINT "potential"
240 LOCATE 30,22 :PRINT "energy E"
250 LOCATE 20,12 :PRINT "eigen function"
260 X2=-10:H=.01
270 FOR I=0 TO 1500   :REM display potential
280 PSET (I,-FNP(X2)),7
290 X2=X2+H
300 NEXT I
310 LINE (1000,0)-(600,20):E=.5
320 FOR KK=1 TO 20
330 C=KK MOD 4 :IF C>0 THEN C=6 ELSE C=1
340 LOCATE 1,1 :PRINT USING "##.#";E
350 X=-SQR(E)-2.3
360 Y=.000001:Z=.0001:H=.01:IMAX=-X/H
370 GOSUB *RUNGE
380 EDG(KK)=Y :E=E+.5
390 LINE (1000-4*KK,-1000*Y+.2*KK)-(1000-4*KK,.2*KK)
400 NEXT KK
410 END
420 *RUNGE         :REM runge kutta integral formula
430 FOR I=0 TO IMAX
440 PSET(1000-IMAX+I-4*KK,-1000*Y+.2*KK),C
450 K1=H*FNF(Z)    :L1=H*FNG(X,Y,E)
460 K2=H*FNF(Z+L1/2):L2=H*FNG(X+H/2,Y+K1/2,E)
470 K3=H*FNF(Z+L2/2):L3=H*FNG(X+H/2,Y+K2/2,E)
480 K4=H*FNF(Z+L3):L4=H*FNG(X+H,Y+K3,E)
490 X=X+H
500 Y=Y+(K1+2*K2+2*K3+K4)/6  :Z=Z+(L1+2*L2+2*L3+L4)/6
510 Y2=Y2+Y*Y :Y(I)=Y
520 NEXT I
530 RETURN
540 END
