home *** CD-ROM | disk | FTP | other *** search
/ Best Objectech Shareware Selections / UNTITLED.iso / boss / educ / math / 023 / mathf.bas (.txt) < prev    next >
Encoding:
GW-BASIC  |  1989-11-24  |  18.3 KB  |  327 lines

  1. 9  GOTO 40
  2. 15                                                                                                                                                    :
  3. 20  + - * / ^ = SQR LOG EXP ABS SIN COS TAN INT RETURN
  4. 25  '
  5. 30  '  << MATHEMATICAL  FUNCTIONS  ---  MATHF.BAS --- by  Dr Russell Langley >>
  6. 35  '
  7. 40  CLEAR:DEF SEG:KEY OFF:ON ERROR GOTO 1610:DEFINT I-N:DIM X(401),Y(401),Z(401),T$(6),TC$(6),U$(9),UC$(9):HD$=" M A T H E M A T I C A L   F U N C T I O N S ":VER$="(RL,8)"
  8. 45  DAT$=MID$(DATE$,4,2)+" "+MID$("JanFebMarAprMayJunJulAugSepOctNovDec",-2+3*VAL(LEFT$(DATE$,2)),3)+" "+RIGHT$(DATE$,4)
  9. 50  F1$="(##)   ######.# ######.###":F2$=F1$+" ######.###":W$="WORKING."
  10. 55  CLS:SCREEN 0,0,0,0:IT=0
  11. 60  PRINT DAT$;TAB(40-LEN(HD$)\2);:COLOR 0,7:PRINT HD$;:COLOR 7,0:PRINT TAB(73)VER$
  12. 65  LOCATE 3,8:PRINT "Press the NUMBER of your choice, and then press <ENTER> to run it.":PRINT TAB(22)"If trouble, try <Ctrl-Break> & GOTO 9.":PRINT STRING$(80,196):LOCATE 7,1:K=14
  13. 70  PRINT:PRINT TAB(K)"1.  Derivative of any Function"CHR$(10)TAB(K)"2.  Exponential Curve Fitting, Y = K + A * (B ^ X)"CHR$(10)TAB(K)"3.  Factorials, Permutations, & Combinations"CHR$(10)TAB(K)"4.  Function Plotter, Y = f(X)"
  14. 75  PRINT TAB(K)"5.  Integration"CHR$(10)TAB(K)"6.  Interpolation, linear"CHR$(10)TAB(K)"7.      ditto    , linear reciprocal"CHR$(10)TAB(K)"8.      ditto    , curvilinear"CHR$(10)TAB(K)"9.  Polynomial Roots"
  15. 80  PRINT TAB(K-1)"10.  Simultaneous Linear Equations, solving"CHR$(10)TAB(K-1)"11.  Return to Main Menu":PRINT
  16. 85  PRINT"     Note:  Data entry only from keyboard in this program.":PRINT
  17. 90  LOCATE ,K-1:INPUT"===> Option (1-11) ";OP:ON OP GOTO 880,630,1520,215,915,455,510,570,790,1080,1630:PRINT "Silly!":GOTO 90
  18. 95  '
  19. 100  ' -- SUBROUTINES
  20. 105  '
  21. 110  ZZ$=STRING$(37-LEN(Z$)\2,177):LOCATE 1,1:PRINT ZZ$"  ";:COLOR 15,0:PRINT Z$;:COLOR 7,0:PRINT "  "ZZ$:RETURN
  22. 115  INPUT" (Y/N)";YN$:IF YN$="" THEN YN$="N":RETURN ELSE YN$=CHR$(ASC(YN$) AND 95):IF YN$="Y" OR YN$="N" THEN RETURN ELSE PRINT"WHAT?  ";:GOTO 115
  23. 120  IF OP<>4 THEN PRINT TAB(46)"Press ANY KEY to continue";:Z$=INKEY$ ELSE PRINT TAB(68)"Press a key";:Z$=INKEY$
  24. 125  IF INKEY$ =""THEN 125 ELSE IF OP<>4 THEN LOCATE ,39:PRINT SPACE$(35):RETURN ELSE LOCATE 25,68:PRINT SPACE$(11);:RETURN
  25. 130  ' -- Read row I of A(I,J) for J=1 to M, in FF.
  26. 135  INPUT X$:IF X$="" THEN PRINT"WHAT";:GOTO 135 ELSE K=0
  27. 140  FOR J=1 TO M:Y$="":QP=0
  28. 145  K=K+1:IF K=LEN(X$) THEN Y$=Y$+RIGHT$(X$,1):IF J=M THEN 155 ELSE J=M:FLAG=1:GOTO 160
  29. 150  Z$=MID$(X$,K,1):IF Z$>" " THEN Y$=Y$+Z$:QP=1:GOTO 145 ELSE IF QP=0 THEN 145
  30. 155  A(I,J)=VAL(Y$)
  31. 160  NEXT J:IF FLAG=1 THEN FLAG=0:BEEP:PRINT"Not enough values.  Please re-do whole line.":GOTO 135 ELSE RETURN
  32. 165  '  -- Line Counter
  33. 170  L=L+1:IF L>20 THEN GOSUB 120:L=1
  34. 175  RETURN
  35. 180  ' -- Permutation Sub
  36. 185  M=N1:E=0:IF X=0 THEN RETURN
  37. 190  FOR J=X TO 1 STEP -1:IF M>=NT THEN M=M/NT:E=E+T
  38. 195  M=M*N:N=N-N1:NEXT J:RETURN
  39. 200  '
  40. 205  ' -- FUNCTION PLOTTER
  41. 210  '
  42. 215  CLS:Z$="F U N C T I O N   P L O T T E R":GOSUB 110
  43. 220  PRINT:PRINT:PRINT"After the `?', enter required function Y=f(X) in BASIC, starting with `Y='":PRINT"Any fractions must have 1 or more digits before the decimal point.":PRINT
  44. 225  IF T$(1)="+" THEN 255
  45. 230  DATA "+","-","*","/","^","=",SQR,LOG,EXP,ABS,SIN,COS,TAN,INT,RETURN
  46. 235  FOR I=1 TO 6:READ T$(I):NEXT I:FOR I=1 TO 9:READ U$(I):NEXT I
  47. 240  START=PEEK(&H30)+256*PEEK(&H31):AD1=PEEK(START)+256*PEEK(START+1):AD2=PEEK(AD1)+256*PEEK(AD1+1)+4:FOR I=1 TO 6:IF PEEK(AD2)<>32 THEN TC$(I)=TC$(I)+CHR$(PEEK(AD2)) ELSE I=I-1
  48. 245  AD2=AD2+1:NEXT I:AD2=AD2+1:FOR I=1 TO 9:IF PEEK(AD2)<>32 THEN UC$(I)=UC$(I)+CHR$(PEEK(AD2)):I=I-1
  49. 250  AD2=AD2+1:NEXT I
  50. 255  INPUT Z$:IF Z$="" THEN IF Z2$>"" THEN PRINT "Ok, I will use: "Z2$:PRINT:GOTO 285 ELSE BEEP:GOTO 255
  51. 260  FOR I=1 TO LEN(Z$):A=ASC(MID$(Z$,I,1)):IF A>96 THEN A=A-32:MID$(Z$,I,1)=CHR$(A)
  52. 265  NEXT I:PRINT:Z2$=Z$
  53. 270  FU$=Z$:AD=AD1+4:FOR I=1 TO LEN(FU$):FOR J=1 TO 6:IF MID$(FU$,I,1)=T$(J) THEN MID$(FU$,I,1)=TC$(J)
  54. 275  NEXT J:NEXT I:FOR J=1 TO 8:U=INSTR(FU$,U$(J)):IF U>0 THEN FU$=LEFT$(FU$,U-1)+UC$(J)+MID$(FU$,U+3)
  55. 280  NEXT J:FU$=FU$+":"+UC$(9)+CHR$(0):FOR I=1 TO LEN(FU$):POKE AD,ASC(MID$(FU$,I,1)):AD=AD+1:NEXT I
  56. 285  IF OP<>4 THEN RETURN
  57. 290  IT=1:INPUT"Smallest X ";SX:INPUT"Biggest  X ";BX:IF BX<=SX THEN PRINT"Silly":GOTO 290
  58. 295  PRINT:PRINT W$:XI=(BX-SX)/400:X=SX:GOSUB 15:SY=Y:BY=SY:FOR I=0 TO 400:X(I)=SX+I*XI:X=X(I):GOSUB 15:Y(I)=Y:IF Y(I)<SY THEN SY=Y(I) ELSE IF Y(I)>BY THEN BY=Y(I)
  59. 300  NEXT I:B=BX:S=SX:H=5:IF ABS(BY-SY)<9.999E-06 THEN PRINT"Y = Constant,"BY:GOTO 435
  60. 305  '
  61. 310  '  -- Select Scales  (H = No. of Axis Labels)
  62. 315  '
  63. 320  Q=1:KT=0
  64. 325  R=B-S:C=S
  65. 330  IF R<=1 THEN KT=KT+1:R=R*10:GOTO 330
  66. 335  IF R>10 THEN KT=KT-1:R=R/10:GOTO 335
  67. 340  IF Q>2 THEN 350 ELSE C=C*10^KT:IF C<0 AND C<>INT(C) THEN C=C-1
  68. 345  C=INT(C)/10^KT:R=(B-C)/H:KT=0:Q=Q+2:GOTO 330
  69. 350  F=INT(R):IF F<>R THEN F=F+1
  70. 355  IF R<0.5 THEN F=F-0.5
  71. 360  F=F/10^KT:IF Q=4 THEN 365 ELSE IF(B-S)/(H*F)<=0.8 THEN KT=1:Q=2:GOTO 325
  72. 365  S=C:D=F*INT(C/F):IF D<0 AND D<>C THEN D=D-F
  73. 370  IF D+H*F>=B THEN S=D
  74. 375  IT=IT+1:IF IT=2 THEN X=S:U=F:H=4:B=BY:S=SY:GOTO 320 ELSE Y=S:V=F
  75. 380  '
  76. 385  ' -- Display Axes (X,Y = Smallest Labels;  U,V = Gapes betw Labels)
  77. 390  '
  78. 395  SCREEN 2:CLS:LINE(70,1)-(70,160):LINE(70,160)-(580,160):FOR I=0 TO 160 STEP 40:LINE (65,I)-(70,I):NEXT I:FOR J=70 TO 570 STEP 100:LINE(J,160)-(J,164):NEXT J
  79. 400  LOCATE 4,3:PRINT"Y";:Q=-2:L(1)=0:L(2)=384:L(3)=832:FOR I=1 TO 3:Q=Q+2:YL(I)=Y+(H-Q)*V:LOCATE 5*Q+1,3:PRINT USING"####.#";YL(I);:NEXT I
  80. 405  LOCATE 22,3:PRINT;:FOR J=1 TO 6:XL(J)=SX+(J-1)*U:PRINT USING"########.#  ";XL(J);:NEXT J:LOCATE 23,64:PRINT"X";
  81. 410  '
  82. 415  '  -- Display Graph
  83. 420  '
  84. 425  X=14:Y=40
  85. 430  XI=(XL(6)-XL(1))/100:YI=(YL(1)-YL(3))/40:FOR I=0 TO 400:X0=X:Y0=Y:X=14+(X(I)-XL(1))/XI:Y=40-(Y(I)-YL(3))/YI:LINE(5*X0,4*Y0)-(5*X,4*Y):NEXT I:LOCATE 25,1:PRINT"Function: "Z2$;:GOSUB 120
  86. 435  LOCATE 25,1:PRINT SPACE$(79);:LOCATE 25,1:PRINT"Another graph of this function";:GOSUB 115:IF YN$="Y" THEN SCREEN 0:GOTO 285 ELSE 55
  87. 440  '
  88. 445  ' -- LINEAR INTERPOLATION
  89. 450  '
  90. 455  CLS:Z$="L I N E A R   I N T E R P O L A T I O N":GOSUB 110
  91. 460  PRINT:PRINT "First enter the X and Y COORDINATES of 2 known points which surround the":PRINT"required interpolate:":PRINT
  92. 465  INPUT;"X(1) ";X1:PRINT SPACE$(6);:INPUT "Y(1) ";Y1:INPUT;"X(2) ";X2:PRINT SPACE$(6);:INPUT "Y(2) ";Y2:PRINT
  93. 470  INPUT"What INTERPOLATE, X ";X
  94. 475  Y=Y1+(X-X1)*(Y2-Y1)/(X2-X1):Y=INT(Y*1000+0.5)/1000
  95. 480  LOCATE CSRLIN-1,40:PRINT "Y = ";Y
  96. 485  PRINT:PRINT"Another point on this line";:GOSUB 115:IF YN$="Y" THEN LOCATE CSRLIN-1,1:PRINT SPACE$(50);:LOCATE ,1:GOTO 470
  97. 490  PRINT:PRINT"Another line";:GOSUB 115:IF YN$="Y" THEN 455 ELSE 55
  98. 495  '
  99. 500  ' -- LINEAR RECIPROCAL INTERPOLATION
  100. 505  '
  101. 510  CLS:Z$="L I N E A R   R E C I P R O C A L   I N T E R P O L A T I O N":GOSUB 110
  102. 515  PRINT:PRINT "First enter the X and Y COORDINATES of 2 known points which surround the":PRINT"required interpolate:":PRINT
  103. 520  INPUT;"X(1) ";X1:PRINT SPACE$(6);:INPUT "Y(1) ";Y1:INPUT;"X(2) ";X2:PRINT SPACE$(6);:INPUT "Y(2) ";Y2:PRINT
  104. 525  INPUT"What INTERPOLATE, X ";X
  105. 530  IF X2>1E+06 THEN P=(X-X1)/X ELSE P=(X1-X)*X2/(X1-X2)/X
  106. 535  Y=Y1+P*(Y2-Y1):Y=INT(Y*1000+0.5)/1000
  107. 540  LOCATE CSRLIN-1,40:PRINT "Y = ";Y
  108. 545  PRINT:PRINT"Another point on this line";:GOSUB 115:IF YN$="Y" THEN LOCATE CSRLIN-1,1:PRINT SPACE$(50);:LOCATE ,1:GOTO 525
  109. 550  PRINT:PRINT"Another line";:GOSUB 115:IF YN$="Y" THEN 510 ELSE 55
  110. 555  '
  111. 560  ' -- CURVILINEAR INTERPOLATION
  112. 565  '
  113. 570  CLS:Z$="C U R V I L I N E A R   I N T E R P O L A T I O N   (LAGRANGE)":GOSUB 110:PRINT
  114. 575  INPUT"Number of known points (3 or more) ";N:IF N<3 THEN 575
  115. 580  PRINT:FOR I=1 TO N:PRINT"X("MID$(STR$(I),2)") ";:INPUT;X(I):PRINT SPACE$(6);:PRINT "Y("MID$(STR$(I),2)") ";:INPUT Y(I):NEXT I:PRINT
  116. 585  INPUT"What INTERPOLATE, X ";X
  117. 590  Y=0:FOR J=1 TO N:T=1:FOR I=1 TO N:IF I=J THEN 595 ELSE T=T*(X-X(I))/(X(J)-X(I))
  118. 595  NEXT I:Y=Y+T*Y(J):NEXT J
  119. 600  LOCATE CSRLIN-1,40:PRINT"Y =";INT(Y*10000+0.5)/10000
  120. 605  PRINT:PRINT"Another point on this curve";:GOSUB 115:IF YN$="Y"THEN LOCATE CSRLIN-1,1:PRINT SPACE$(50);:LOCATE ,1:GOTO 585
  121. 610  PRINT:PRINT"Another curve";:GOSUB 115:IF YN$="Y"THEN 570 ELSE 55
  122. 615  '
  123. 620  '  -- EXPONENTIAL CURVES
  124. 625  '
  125. 630  CLS:Z$="E X P O N E N T I A L   C U R V E   F I T T I N G":GOSUB 110:MS$="METHOD OF SEMIAVERAGES":PRINT STRING$(37-LEN(MS$)\2,32)MS$:PRINT:PRINT
  126. 635  INPUT"No. of equally spaced X-Values (any multiple of 3) ";N:Z=N/3:IF N<3 OR INT(Z)<>Z THEN 635 ELSE NN=N:M=Z
  127. 640  INPUT"First X-Value ";X(1):IF X(1)>100 THEN PRINT"Too big. Code X-Values by subtraction so that all X's are < 100.":GOTO 640 ELSE INPUT"Interval between successive X-Values ";XI:FOR I=2 TO N:X(I)=X(I-1)+XI:NEXT I
  128. 645  PRINT:IF N=3 THEN PRINT"Enter 3 Y-Values read from your hand-drawn smooth curve:"ELSE PRINT"ENTER";N;"Observed Y-Values corresponding to observed X-Values:"
  129. 650  FOR I=1 TO N
  130. 655  PRINT"Y(";MID$(STR$(I),2);") ";:INPUT Y$:IF Y$>" " THEN Y(I)=VAL(Y$):NEXT I  ELSE BEEP:GOTO 655
  131. 660  K=0:FOR I=1 TO 3:S(I)=0:FOR J=1 TO M:K=K+1:S(I)=S(I)+Y(K):NEXT J:NEXT I:BM=(S(2)-S(3))/(S(1)-S(2)):B=BM^(1/M):A=(S(1)-S(2))*(1-B)/(1-BM)^2:K!=(S(1)*S(3)-S(2)^2)/(M*(S(1)+S(3)-2*S(2))):K!=INT(K!*10000+0.5)/10000
  132. 665  BM=(S(2)-S(3))/(S(1)-S(2)):B=BM^(1/M):A=(S(1)-S(2))*(1-B)/(1-BM)^2:K!=(S(1)*S(3)-S(2)^2)/(M*(S(1)+S(3)-2*S(2))):K!=INT(K!*10000+0.5)/10000
  133. 670  A=A/B^(X(1)/XI):AA=INT(A*10000+0.5)/10000:B=B^(1/XI):BB=INT(B*10000+0.5)/10000
  134. 675  PRINT:PRINT"FUNCTION:  Y ="K!;:IF A>0 THEN PRINT"+"AA;ELSE PRINT AA;
  135. 680  PRINT"* ("BB"^ X )":PRINT:PRINT "Estimates based on this equation:"
  136. 685  FOR I=1 TO N:Z(I)=K!+A*B^X(I):NEXT I:PRINT SPACE$(13)"X      EST Y":FOR I=1 TO N:PRINT USING F1$;I;X(I);Z(I):IF INT(I/15)=I/15 THEN GOSUB 120
  137. 690  NEXT I:PRINT
  138. 695  PRINT"Do you want more estimates based on this equation";:GOSUB 115:IF YN$="N" THEN 720ELSE PRINT"Ok, use a slash '/' to signal end."
  139. 700  FOR I=N+1 TO 100:INPUT"X ";X$:IF X$="/"THEN NN=I-1:I=100 ELSE X(I)=VAL(X$):Z(I)=K!+A*B^X(I):Z(I)=INT(Z(I)*1000+0.5)/1000:LOCATE CSRLIN-1,15:PRINT "Y ="Z(I)
  140. 705  NEXT I:PRINT
  141. 710  '
  142. 715  '         -- Printout
  143. 720  PRINT"Do you want printout of these results";:GOSUB 115:IF YN$="N" THEN 770
  144. 725  LINE INPUT "PROBLEM ID ? ";ID$
  145. 730  BEEP:LOCATE 24,30:LINE INPUT;"Press <ENTER> when PRINTER is ON.";Z$:LOCATE 24,30:PRINT SPACE$(50);
  146. 735  LPRINT CHR$(14)"EXPONENTIAL CURVE FITTING"CHR$(20) '= W I D E letters on/off with Epson Dot Matrix Printers.
  147. 740  LPRINT SPACE$(12)MS$
  148. 745  LPRINT" ":LPRINT"Problem: "ID$:LPRINT" ":LPRINT"DATA         X      0BS Y      EST Y"
  149. 750  FOR I=1 TO N:LPRINT USING F2$;I;X(I);Y(I);Z(I):NEXT I:LPRINT " "
  150. 755  LPRINT"FUNCTION: Y ="K!;:IF A>0 THEN LPRINT"+"AA;ELSE LPRINT AA;
  151. 760  LPRINT"* ("BB"^ X )":LPRINT " ":IF NN>N THEN LPRINT"ALSO         X      EST Y"
  152. 765  FOR I=N+1 TO NN:LPRINT USING"       ######.# ######.###";X(I);Z(I):NEXT I:LPRINT " ":LPRINT " "
  153. 770  PRINT"Do you want another Exponential Curve";:GOSUB 115:IF YN$="Y" THEN 630ELSE 55
  154. 775  '
  155. 780  '  -- POLYSOLVER (after B. P. Douglass)
  156. 785  '
  157. 790  CLS:CLEAR:DEFDBL A-Z:DEFINT I-N:Z$="P O L Y S O L V E R":GOSUB 110:PRINT:PRINT SPACE$(5)"This seeks N roots of:  Y = A + B(1)*X +B(2)*X^2 + ... + B(N)*X^N":PRINT SPACE$(10)"(Complex or unstable roots will cause non-convergence)":PRINT
  158. 795  INPUT"Order N of your polynomial (1-10) ";N:NN=N
  159. 800  INPUT"A ";A(0):FOR I=1 TO N:PRINT"B("MID$(STR$(I),2)") ";:INPUT A(I):NEXT I:FOR I=0 TO N:D(I)=A(I):NEXT I
  160. 805  PRINT:PRINT"Approx roots will first be found using NEWTON & HORNER'S method.":PRINT"Enter trial value for each root in turn:":FOR J=1 TO NN:K=0
  161. 810  PRINT"Trial Root #"J;:INPUT X:PRINT W$;:LOCATE ,1
  162. 815  GOSUB 860:IF ABS(B(0))<9.99E-07 THEN 825ELSE IF C(1)=0 THEN C(1)=9.9999E-05
  163. 820  IF K>29 THEN 850ELSE K=K+1:X=X-B(0)/C(1):GOTO 815
  164. 825  X(J)=X:FOR I=1 TO N:A(I-1)=B(I):NEXT I:N=N-1:PRINT"Rough Root #"J"= "CSNG(X)"  after "K" iterations.":NEXT J
  165. 830  PRINT:PRINT"Improved ANSWERS obtained by re-solving with above approx roots:":N=NN:FOR I=0 TO N:A(I)=D(I):NEXT I:FOR J=1 TO N:PRINT "WORKING";:LOCATE ,1:X=X(J)
  166. 835  GOSUB 860:IF ABS(B(0))<9.8E-08 THEN 845ELSE IF C(1)=0 THEN C(1)=9.999E-06
  167. 840  X=X-B(0)/C(1):GOTO 835
  168. 845  X(J)=X:PRINT"Root #"J"= "CSNG(X):NEXT:GOSUB 120:GOTO 9
  169. 850  PRINT"Incomplete convergence after 30 iterations.":PRINT"Rough Root = "X"  but differences still = "B(0)
  170. 855  PRINT"Is this acceptable";:GOSUB 115:IF YN$="Y"THEN 825ELSE K=0:GOTO 810
  171. 860  B(N)=A(N):FOR I=N-1 TO 0 STEP -1:B(I)=A(I)+B(I+1)*X:NEXT I:C(N)=B(N):FOR I=N-1 TO 1 STEP -1:C(I)=B(I)+C(I+1)*X:NEXT I:RETURN
  172. 865  '
  173. 870  ' -- DERIVATIVES
  174. 875  '
  175. 880  CLS:Z$="DERIVATIVE OF ANY FUNCTION":GOSUB 110:GOSUB 220:PRINT
  176. 885  INPUT"At what value of X do you want the derivative ";XD:D=0.000999999*XD:IF D=0 THEN D=9.99E-07 ELSE IF ABS(D)<9.99E-07 THEN D=9.99E-07*SGN(D)
  177. 890  X=XD+D:GOSUB 15:YD=Y:X=XD-D:GOSUB 15:PRINT"At X = "STR$(XD)", the derivative f'(X) = "(YD-Y)/2/D
  178. 895  PRINT:PRINT"Another derivative of this function";:GOSUB 115:IF YN$="Y"THEN LOCATE CSRLIN-1,1:GOTO 885ELSE 55
  179. 900  '
  180. 905  '  -- INTEGRATION by SIMPSON'S RULE
  181. 910  '
  182. 915  CLS:Z$="I N T E G R A T I O N":GOSUB 110:PRINT
  183. 920  INPUT"Simpson's rule or Romberg's algorithm (S/R) ";Z$:IF Z$="" THEN BEEP:GOTO 920 ELSE Z$=CHR$(ASC(Z$) AND 95):IF Z$="S"THEN PRINT ELSE IF Z$="R"THEN 1020 ELSE BEEP:GOTO 920
  184. 925  PRINT"Will you enter a FORMULA, or equally spaced Y VALUES (F/V) ";
  185. 930  INPUT KIND$:IF KIND$=""THEN BEEP:GOTO 925ELSE KIND$=CHR$(ASC(KIND$) AND 95):IF KIND$="F"THEN GOSUB 220 ELSE IF KIND$<>"V"THEN BEEP:GOTO 925
  186. 935  INPUT"Start integrating at what X value ";X0:INPUT"End integrating at what X value ";XN:IF XN<=X0 THEN PRINT"SILLY":GOTO 935
  187. 940  YE=0:YO=0:IF KIND$="V" THEN 970
  188. 945  INPUT"Enter an Even Number of X Intervals (2-100) ";N:IF N<2 OR N/2<>INT(N/2)THEN PRINT"That's not an even number.":BEEP:GOTO 945
  189. 950  PRINT:PRINT W$;:LOCATE ,1:W=(XN-X0)/N:X=X0:GOSUB 15:Y(0)=Y:X=XN:GOSUB 15:Y(N)=Y
  190. 955  FOR I=1 TO N-1:X=X0+I*W:GOSUB 15:Y(I)=Y
  191. 960  IF I/2=INT(I/2)THEN YE=YE+Y(I) ELSE YO=YO+Y(I)
  192. 965  NEXT I:GOTO 995
  193. 970  PRINT"Enter an Odd Number (3-101) of equally spaced X Values at which Y is known,":INPUT "including the Starting & Ending Values ";N:IF N<3 OR N/2=INT(N/2)THEN PRINT"That's not an odd number.":BEEP:GOTO 970
  194. 975  PRINT:N=N-1:W=(XN-X0)/N
  195. 980  FOR I=0 TO N:PRINT"At X = ";X0+I*W;TAB(20)"Y("MID$(STR$(I+1),2)") ";:INPUT Y(I)
  196. 985  IF I=0 OR I=N THEN 990 ELSE IF I/2=INT(I/2)THEN YE=YE+Y(I) ELSE YO=YO+Y(I)
  197. 990  NEXT I:PRINT
  198. 995  PRINT"INTEGRAL AREA = ";W/3*(Y(0)+2*YE+4*YO+Y(N))"  (SIMPSON)":PRINT
  199. 1000  PRINT"Another Integral of this function";:GOSUB 115:IF YN$="Y" THEN 935ELSE 55
  200. 1005  '
  201. 1010  ' -- INTEGRATION by ROMBERG'S ALGORITHM (after B. P. Douglass)
  202. 1015  '
  203. 1020  CLEAR:ON ERROR GOTO 1610:GOSUB 220
  204. 1025  INPUT"Single or Double precision (S/D) ";Z$:IF Z$="D" OR X$="d" THEN DEFDBL R,S,W,X
  205. 1030  INPUT"Start integrating at what X value ";X0:INPUT"End integrating at what X value ";XN:IF XN<=X0 THEN PRINT"SILLY":GOTO 1030
  206. 1035  INPUT"How many iterations (2-15) ";N:IF N<2 OR N>15 THEN BEEP:GOTO 1035
  207. 1040  DIM R(N,N):W=XN-X0:X=X0:GOSUB 15:Y0=Y:X=XN:GOSUB 15:R(1,1)=W*(Y0+Y)/2:L=1
  208. 1045  FOR I=2 TO N:W=W/2:L=L+L:S=0:FOR K=1 TO L-1 STEP 2:X=X0+W*K:GOSUB 15:S=S+Y:NEXT K:R(I,1)=R(I-1,1)/2+W*S:PRINT"R("MID$(STR$(I),2)",1) ="R(I,1);
  209. 1050  EM=1:FOR J=2 TO I:EM=EM*4:R(I,J)=R(I,J-1)+(R(I,J-1)-R(I-1,J-1))/(EM-1):IF J>= I THEN I$=MID$(STR$(I),2):PRINT"   R("I$","I$") ="R(I,J)
  210. 1055  NEXT J
  211. 1060  NEXT I:PRINT"FINAL INTEGRAL ="R(N,N)"  (ROMBERG)":PRINT:LINE INPUT "Press <ENTER> to return to Math Functions Menu. ";Z$:GOTO 9
  212. 1065  '
  213. 1070  ' -- SIMULTANEOUS EQUATIONS
  214. 1075  '
  215. 1080  CLS:Z$=" SOLVING  SIMULTANEOUS  LINEAR  EQUATIONS ":GOSUB 110
  216. 1085  PRINT:PRINT TAB(18)"B1.X1   +   B2.X2   +  ...  +   BN.XN  =  Y":PRINT
  217. 1090  PRINT"     Given values of X's & Y for N equations, this program solves":PRINT"     for the unknown B coefficients, if possible.":PRINT"     The computations can be done with:"
  218. 1095  PRINT"        Single precision (6 digit accuracy), or":PRINT"        Double precision (16 digit accuracy, but slower).":PRINT
  219. 1100  INPUT"Which do you want (S/D) ";SD$
  220. 1105  IF SD$="S" OR SD$="s" THEN A0=9.99E-07:A3=0.000999999 ELSE IF SD$="D" OR SD$="d" THEN DEFDBL A,B,C,D,S,P,Y:A0=9.9999E-13:A3=0.001 ELSE PRINT"Please reply 'S' for Single, or 'D' for Double precision.":GOTO 1100
  221. 1110  PRINT
  222. 1115  INPUT"How many unknown coefficients, N (2-10) ";N:IF N<2 OR N>10 THEN 1115 ELSE M=N+1
  223. 1120  DIM A(N,M),B(N),C(N,N),D(N,N),IP(N),KR(N),KC(N),P(N)
  224. 1125  F4$="####.##### ":F5$="#####.#### ":F8$="########.# ":R$="Row##: ":D$="Determinant of X Matrix = "
  225. 1130  '   --Enter data from Keyboard.
  226. 1135  PRINT:PRINT"Enter"N"rows of data, each with X1 to X"MID$(STR$(N),2)" & Y, in Free Format:":I=1
  227. 1140  PRINT"Row"I;:GOSUB 135:I=I+1:IF I<=N THEN 1140
  228. 1145  CLS:IF QR=0 THEN PRINT"Data was read as ---" ELSE PRINT"Revised data is---"
  229. 1150  FOR I=1 TO N:PRINT USING R$;I;:FOR J=1 TO M:PRINT A(I,J);:NEXT J:PRINT:NEXT I:PRINT
  230. 1155  PRINT"Do you want to make any CHANGES to those values";:GOSUB 115:IF YN$="N" THEN 1220 ELSE INPUT"Which row";I
  231. 1160  IF I<1 OR I>N THEN PRINT"WHAT?   ";:GOTO 1155
  232. 1165  '   --Edit a row.
  233. 1170  PRINT USING R$;I;:FOR J=1 TO M:PRINT A(I,J);:NEXT J:PRINT
  234. 1175  PRINT"Enter DATUM # (1-"MID$(STR$(M),2)") & RIGHT VALUE, in Free Format: ";
  235. 1180  INPUT "",X$:IF X$="" THEN 1155 ELSE QR=1:K=0
  236. 1185  FOR J=1 TO 2:Y$="":QP=0
  237. 1190  K=K+1:IF K=LEN(X$) THEN Y$=Y$+RIGHT$(X$,1):IF J=2 THEN 1200 ELSE J=2:FLAG=1:GOTO 1205
  238. 1195  Z$=MID$(X$,K,1):IF Z$>" " THEN Y$=Y$+Z$:QP=1:GOTO 1190 ELSE IF QP=0 THEN 1190
  239. 1200  Y(J)=VAL(Y$):IF Y(1)<1 OR Y(1)>M THEN PRINT"WHAT?   ";:GOTO 1175
  240. 1205  NEXT J:IF FLAG=1 THEN FLAG=0:PRINT "No, please enter 2 VALUES on the same line.":GOTO 1175
  241. 1210  J=Y(1):A(I,J)=Y(2):GOTO 1145
  242. 1215  '   -- Calc now --
  243. 1220  COLOR 23,0:PRINT:PRINT"Working.":COLOR 7,0
  244. 1225  '   -- Copy A into Square C of order N.
  245. 1230  FOR I=1 TO N:FOR J=1 TO N:C(I,J)=A(I,J):NEXT J:NEXT I
  246. 1235  '   -- Inverse of Diag C(n,n) --<UNK! {000A}>      Uses I-M,D,Q,Z,FLAG,IR,IC,IP(n),KR(n),KC(n),P(n)
  247. 1240  D=1:FLAG=0:FOR I=1 TO N:IP(I)=0:NEXT I:FOR I=1 TO N:Z=0:FOR J=1 TO N:IF IP(J)=1 THEN 1265
  248. 1245  FOR K=1 TO N:ON SGN(IP(K)-1)+2 GOTO 1255,1260,1250
  249. 1250  I=N:J=N:K=N:FLAG=1:GOTO 1260
  250. 1255  IF Z < ABS(C(J,K)) THEN IR=J:IC=K:Z=C(J,K)
  251. 1260  NEXT K
  252. 1265  NEXT J:IF FLAG=1 THEN D=0:GOTO 1295
  253. 1270  IP(IC)=IP(IC)+1:IF IR<>IC THEN D=-D:FOR L=1 TO N:Z=C(IR,L):C(IR,L)=C(IC,L):C(IC,L)=Z:NEXT L
  254. 1275  KR(I)=IR:KC(I)=IC:P(I)=C(IC,IC):D=D*P(I):IF ABS(P(I))<=9.99E-07 THEN D=0:I=N:GOTO 1295
  255. 1280  C(IC,IC)=1:FOR L=1 TO N:C(IC,L)=C(IC,L)/P(I):NEXT L
  256. 1285  FOR K=1 TO N:IF K<>IC THEN Z=C(K,IC):C(K,IC)=0:FOR L=1 TO N:C(K,L)=C(K,L)-C(IC,L)*Z:NEXT L
  257. 1290  NEXT K
  258. 1295  NEXT I:IF D=0 THEN 1360
  259. 1300  FOR I=1 TO N:Q=N-I+1:IF KR(Q)=KC(Q) THEN 1310 ELSE K=KR(Q):L=KC(Q)
  260. 1305  FOR J=1 TO N:Z=C(J,K):C(J,K)=C(J,L):C(J,L)=Z:NEXT J
  261. 1310  NEXT I
  262. 1315  '   -- Confirm that C*C(INV) = I.
  263. 1320  FOR I=1 TO N:FOR J=1 TO N:D(I,J)=0:FOR K=1 TO N:D(I,J)=D(I,J)+A(I,K)*C(K,J):NEXT K:NEXT J:NEXT I
  264. 1325  FOR I=1 TO N:FOR J=1 TO N:IF I=J AND ABS(D(I,I)-1)>A3 THEN J=N:I=N:FLAG=1 ELSE IF I<>J AND ABS(D(I,J))>A3 THEN J=N:I=N:FLAG=1
  265. 1330  NEXT J
  266. 1335  NEXT I
  267. 1340  IF FLAG=1 THEN D=0:GOTO 1360
  268. 1345  '  -- Compute Coefficients
  269. 1350  FOR I=1 TO N:B(I)=0:FOR J=1 TO N:B(I)=B(I)+C(I,J)*A(J,M):NEXT J:NEXT I
  270. 1355  '  -- Show Answers.
  271. 1360  CLS:PRINT D$;D:PRINT:IF D<>0 THEN 1380
  272. 1365  PRINT"These equations are NOT soluble.":IF SD$="S" OR SD$="s" THEN PRINT"Try Double Precision."
  273. 1370  PRINT:GOSUB 120:GOTO 1425
  274. 1375  '   -- Show B Coefficients
  275. 1380  PRINT"B COEFFICIENTS:":FOR I=1 TO N:PRINT SPACE$(5)"B("MID$(STR$(I),2)") = "B(I):NEXT I:PRINT:GOSUB 120
  276. 1385  '   -- Show Inverse X Matrix
  277. 1390  PRINT"INVERSE X MATRIX:":L=2:C=ABS(C(1,1)):IF C<10 THEN F$=F4$ ELSE IF C<1000 THEN F$=F5$ ELSE F$=F8$
  278. 1395  FOR I=1 TO N:PRINT USING R$;I;:IF N>5 THEN MM=5 ELSE MM=N
  279. 1400  FOR J=1 TO MM:PRINT USING F$;C(I,J);:NEXT J:PRINT:GOSUB 170
  280. 1405  IF N>5 THEN PRINT SPACE$(7);:FOR J=6 TO N:PRINT USING F$;C(I,J);:NEXT J:PRINT:GOSUB 170
  281. 1410  NEXT I:GOSUB 120
  282. 1415  '
  283. 1420  ' -- Printout
  284. 1425  PRINT"Do you want PRINTOUT of these results";:GOSUB 115:IF YN$="N" THEN 1500
  285. 1430  LINE INPUT "Problem ID ? ";ID$
  286. 1435  BEEP:LOCATE 24,30:LINE INPUT;"Press <ENTER> when PRINTER is ON.";Z$:LOCATE 24,30:PRINT SPACE$(50);
  287. 1440  LPRINT CHR$(14)"SOLVING LINEAR EQUATIONS"CHR$(20) '= W I D E letters on/off with Epson Dot Matrix Printers.
  288. 1445  LPRINT" ":LPRINT"PROBLEM: ";ID$:LPRINT STRING$(7,"="):LPRINT" "
  289. 1450  LPRINT"DATA:":FOR I=1 TO N:LPRINT USING R$;I;:FOR J=1 TO M:LPRINT A(I,J);:NEXT J:LPRINT" ":NEXT I:LPRINT" "
  290. 1455  LPRINT D$;D:LPRINT" ":IF D<>0 THEN 1470
  291. 1460  LPRINT"These equations are NOT soluble.":IF SD$="S" OR SD$="s" THEN LPRINT"Try Double Precision."
  292. 1465  LPRINT" ":LPRINT" ":GOTO 1500
  293. 1470  LPRINT"B COEFFICIENTS:":FOR I=1 TO N:LPRINT SPACE$(5)"B("MID$(STR$(I),2)") = "B(I):NEXT I:LPRINT" "
  294. 1475  LPRINT"INVERSE X MATRIX:":L=2:C=ABS(C(1,1)):IF C<10 THEN F$=F4$ ELSE IF C<1000 THEN F$=F5$ ELSE F$=F8$
  295. 1480  FOR I=1 TO N:LPRINT USING R$;I;:IF N>5 THEN MM=5 ELSE MM=N
  296. 1485  FOR J=1 TO MM:LPRINT USING F$;C(I,J);:NEXT J:LPRINT" "
  297. 1490  IF N>5 THEN FOR J=6 TO N:LPRINT USING F$;C(I,J);:NEXT J:LPRINT" "
  298. 1495  NEXT I:LPRINT" ":LPRINT" "
  299. 1500  GOTO 9
  300. 1505  '
  301. 1510  ' -- FACTORIALS etc.
  302. 1515  '
  303. 1520  CLS:DEFDBL L,M,N:N1=1:NT=1E+10:T=10:LN=2.30259:LG=0.434294
  304. 1525  Z$="FACTORIALS, PERMUTATIONS, & COMBINATIONS":GOSUB 110:PRINT:PRINT
  305. 1530  PRINT TAB(10)"This program computes Double Precision  n!,  nPx,  or  nCx,":PRINT TAB(10)"using a separate exponent which enables products to greatly":PRINT TAB(10)"exceed the normal upper limit of 1E38.":PRINT
  306. 1535  INPUT"Factorial, Permutation, Combination, or Menu (F/P/C/M)";X$:IF X$="" THEN 9 ELSE X$=CHR$(ASC(X$) AND 95):IF X$="M" THEN 9 ELSE IF INSTR("FPC",X$)=0 THEN BEEP:GOTO 1535
  307. 1540  INPUT"n";N:IF N<2 THEN BEEP:PRINT"n must be 2 or more!":GOTO 1540 ELSE IF X$<>"F" THEN 1550
  308. 1545  PRINT"n! =  ";:X=N:GOSUB 185:GOTO 1565' n! = nPn
  309. 1550  INPUT"x";X:IF X<0 OR X>N THEN BEEP:PRINT"x must be > 0 and < n.":GOTO 1550 ELSE GOSUB 185
  310. 1555  IF X$="P" THEN PRINT"nPx = ";:GOTO 1565' nPx
  311. 1560  PRINT"nCx = ";:MM=M:EE=E:N=VAL(STR$(X)):GOSUB 185:M=MM/M:E=EE-E' nCx = nPx/xPx
  312. 1565  FOR J=1 TO 12:IF E=0 OR M>1E+11 THEN J=12 ELSE M=M*10:E=E-1
  313. 1570  NEXT J
  314. 1575  M=INT(M+0.5):PRINT USING"###,###,###,###";M;:PRINT" E+";RIGHT$(STR$(E),LEN(STR$(E))-1)
  315. 1580  PRINT USING"LOG(E)  =#####.######";INT((LOG(M)+LN*E)*1E+06+0.5)/1E+06
  316. 1585  PRINT USING"LOG(10) =#####.######";INT((E+LG*LOG(M))*1E+06+0.5)/1E+06
  317. 1590  PRINT:GOTO 1535
  318. 1595  '
  319. 1600  ' -- ERRORS & END
  320. 1605  '
  321. 1610  IF ERR=0 THEN 1630
  322. 1615  IF ERR=11 THEN Y=100:RESUME NEXT
  323. 1620  BEEP:BEEP:IF ERR=2 AND ERL=15 THEN PRINT"  ERROR IN YOUR FORMULA:  "Z2$:PRINT"  Please re-enter it in correct BASIC (starting with `Y=' )":RESUME 255
  324. 1625  ON ERROR GOTO 0:GOTO 1635
  325. 1630  RUN "menu.bas"
  326. 1635  END
  327.