100'  NAME--DE-2OR
110'
120'  DESCRIPTION--SOLVES THE INITIAL VALUE PROBLEM
130'     Y'' = F(X,Y,Y')
140'     Y(X0) = Y0
150'     Y'(X0)= Y0'
160'
170'  SOURCE--RUNGA KUTTA METHOD (2ND ORDER ACCURACY)
180'
190'  INSTRUCTIONS--THE INTEGRATION STEP IS H AND THE SOLUTION IS
200'  TABULATED ON THE INTERVAL X0<= X <=B IS STEPS OF SIZE L.
210'
220'     430 DEF FNF(X,Y,Z) = .........
230'     440 DATA X0,Y0,Y0',B
240'     450 DATA H,L
250'
260'  THE VALUES OF Y' MAY BE INCLUDED IN THE TABULATION BY
270'  TYPING THE FOLLOWING LINES:
280'     540 PRINT"VALUE OF X","VALUE OF Y","VALUE OF Y'"
290'     560 PRINT X,Y,Z
300'     830 PRINT Q,Y2,Z2
310'
320'  FOR MODIFICATIONS OF THE PRINT ROUTINE, FOR EXAMPLE TO PLOT
330'  THE SOLUTION RATHER THAN TO TABULATE ITS VALUES, THE FOLLOWING
340'  LINES SHOULD BE MODIFIED:
350'     540 PRINT "VALUE OF X","VALUE OF Y"
360'     550 PRINT
370'     560 PRINT X,Y,Z
380'     830 PRINT Q,Y2,Z2
390'
400'
410'  *  *  *  *  *  *  *  MAIN PROGRAM  *  *  *  *  *  *  *  *  *  *
420'
430 DEF FNF(X,Y,Z)=-Y-Z/X
440 DATA .01,1,0,20
450 DATA .01,.25
460 READ X0,Y,Z,B,H,L
470 LET X=X0
480 LET E=1/L
490 IF (B-X0)*H>0 THEN 510
500 LET H=-H
510 IF (B-X0)*L>0 THEN 530
520 LET L=-L
530 REM
540 PRINT "VALUE OF X","VALUE OF Y"
550 PRINT
560 PRINT X,Y
570 GOTO 630
580 LET P1=Y+H*Z
590 LET Q1=Z+H*FNF(X,Y,Z)
600 LET Y=(Y+P1)/2+.5*H*Q1
610 LET Z=(Z+Q1)/2+.5*H*FNF(X+H,P1,Q1)
620 LET X=X+H
630 IF (X+H-B)*SGN(L)+1E-7 >= 0 THEN 660
640 LET R=X+H
650 GOTO 670
660 LET R=B
670 LET A=INT(E*R)/E + L*(SGN(H)-1)/2
680 IF (A+L-R)*SGN(L)-1E-7 >= 0 THEN 700
690 LET A=A+L
700 IF (X-A)*SGN(L)+1E-7 >= 0 THEN 730
710 LET Q=A
720 GOSUB 790
730 IF R=B THEN 750
740 GOTO 580
750 IF ABS(R-A)<.5*1E-6 THEN 780
760 LET Q=R
770 GOSUB 790
780 STOP
790 LET P2=Y+(Q-X)*Z
800 LET Q2=Z+(Q-X)*FNF(X,Y,Z)
810 LET Y2=(Y+P2)/2+.5*(Q-X)*Q2
820 LET Z2=(Z+Q2)/2+.5*(Q-X)*FNF(Q,P2,Q2)
830 PRINT Q,Y2
840 RETURN
850 END