'===========================================================================
' Subject: FUNCTION GRAPHER                  Date: 11-24-96 (00:00)       
' Author:  Toshihiro Horie                   Code: QB, QBasic, PDS        
' Origin:  www.ocf.berkley.edu/~horie/     Packet: ALGOR.ABC
'===========================================================================
' Improved Function Grapher                     '
' Log, scaling, incr, and root bugs fixed.      '
' Copyright (c) 1996 Toshihiro Horie            '

CONST version = 6.72
DEFSNG A-Z
DECLARE FUNCTION atan! (x!)
DECLARE FUNCTION asin! (x!)
DECLARE FUNCTION acos! (x!)
DECLARE FUNCTION csc! (x!)
DECLARE FUNCTION sec! (x!)
DECLARE FUNCTION cot! (x!)
DECLARE FUNCTION log10 (x) : DECLARE FUNCTION ln (x)
DECLARE FUNCTION sinh (x) : DECLARE FUNCTION csch (x)
DECLARE FUNCTION cosh (x) : DECLARE FUNCTION sech (x)
DECLARE FUNCTION tanh (x) : DECLARE FUNCTION coth (x)
DECLARE FUNCTION acosh! (x!) : DECLARE FUNCTION asech! (x!)
DECLARE FUNCTION asinh! (x!) : DECLARE FUNCTION acsch! (x!)
DECLARE FUNCTION acoth! (x!) : DECLARE FUNCTION atanh! (x!)
DECLARE SUB grid (x1!, y1!, XS!, ys!, XC!, YC!, XN!, yn!, xm!, ym!)
DECLARE SUB radianmarks (XC, YC, XS, ys, XN)

CLEAR , , 32767: 'allocate more stack to error handler
ON ERROR GOTO 60

'ln(x) divided by LOG(10)/LOG(e) gives common log
'Common LOG is defined in function LOG10
'LN natural log is defined in func LN
CONST l = 2.3025851#
CONST e = 2.7182818#
CONST pi = 3.14159265357989#
CONST DX = 1 / 1048576

'##########################################################
  XS = 640: ys = 480:        'dimensions of the screen
  xm = 10: ym = 16:          'X AND Y DIMENSIONS OF GRID
'##########################################################


SCREEN 12: CLS : COLOR 15: LOCATE 1, 1
PRINT "FUNCTION GRAPHER Version "; version
PRINT "by Toshihiro Horie Rev 04/06/96"
COLOR 7: LOCATE 24, 4: PRINT "Scale: 1 green unit= 1"
GOSUB VARS
XC = XS \ 2: YC = ys \ 2
XN = XS / xm: yn = INT(XS / ym * 4 / 3)
incr = xm / XC * 4: M1 = 0
CALL grid(x1, y1, XS, ys, XC, YC, XN, yn, xm, ym)
REM =======================START GRAPHING!!!=================================
xi = incr
DEF SEG = 0: ' if scroll lock is on, then double the speed
IF (PEEK(&H417) AND &H10) THEN incr = incr * 2: DEF SEG

'DRAW FIRST POINT.......................
x = -XC / XN - 1
GOSUB equation
XP = x * XN + XC
YP = -Y * yn + YC
PSET (XP, YP), 15

DO
LOCATE 1, 67:  PRINT TIME$
'....................................................KEY CHECK....
in$ = INKEY$: xmode = 1
  IF in$ = "*" THEN xmode = 16
  IF in$ = "+" THEN xmode = 4
  IF in$ = "-" THEN xmode = 1 / 4
  IF in$ = "." THEN xmode = 8
  IF in$ = " " THEN xmode = -4: 'reverse
  IF in$ = "," THEN xmode = 1 / 8
xi = incr * xmode

IF UCASE$(in$) = "Q" THEN GOSUB ASK
IF in$ = CHR$(27) THEN END
'INCREMENT THE X, FIND SLOPE............
'if flat curve then speed up
'IF ABS(M1) < .1 THEN xmode = 1.5
x = x + xi
GOSUB equation
M1 = (Y - yold) / (x - xold)
'.................................................BOUNDARY CHECK..1
IF Y > (YC / yn) + 16 OR Y < (-YC / yn) - 16 THEN
   Errflag2 = 1
   xi = incr * 2
   GOTO SKIP
 ELSE
   xi = incr
END IF
'TRANSLATE TO SCREEN COORDINATES..................................1.5
XP = x * XN + XC
YP = -Y * yn + YC
IF Errflag2 = 1 THEN PSET (XP, YP), 15: Errflag2 = 0
IF Errflag = 1 THEN PSET (XP, YP), 15: Errflag = 0
'DRAW LINE F(X)...................................SLOPE CHECK.....2
IF ABS(M1) > 64 THEN xi = incr / 8
IF ABS(M1) > 256 THEN
        ' vert asymptote - probably
        'LINE (XP, 0)-(XP, 479), 7, , &HCCCC
        xi = incr / 16
        PSET (XP, YP), 15
ELSE
        xi = incr
        LINE -(XP, YP), 15
END IF
'CIRCLE (XP, YP), 3, 12
'PSET (XP - 1, YP), 14: PSET (XP + 1, YP), 14: PSET (XP, YP), 14
'FIND CRITICAL POINTS.............................................3
IF MOLD * M1 < 0 THEN '(Mean Value Theorem)
      CIRCLE (xold * XN + XC, -yold * yn + YC), 4, 11
      COLOR 11
      LOCATE 28, 4: PRINT USING "(+###.###"; xold;
      PRINT USING ",+#####.###) "; yold
      COLOR 15
      PSET (XP, YP), 15
END IF
SKIP:
'SHOW SLOPE AND X,Y COORDINATE....................................4
LOCATE 25, 4: PRINT USING "(+###.##"; x;
PRINT USING ",+#####.##) "; Y
IF SGN(Y) * SGN(yold) <= 0 THEN '...............x-roots (misses mins)
      LOCATE 27, 4: PRINT USING "RootX=###.##"; xold
END IF
MOLD = (Y - yold) / (x - xold)
xold = x: yold = Y
IF x <> 0 THEN LOCATE 26, 4: PRINT USING "slope=+####.##   "; MOLD

skip2:
LOOP WHILE x < (XC / XN)
END
'======================================================================

ASK:
xold = x
LOCATE 1, 1: PRINT STRING$(80, 255);
LOCATE 1, 1: INPUT "X COORDINATE"; x: GOSUB equation: y0 = Y
x = x + DX: GOSUB equation: der1 = (Y - y0) / (DX)
LOCATE 1, 38: PRINT USING "Y COORDINATE IS +######.###"; y0
LOCATE 2, 35: PRINT USING "1st DERIVATIVE  IS +######.####"; der1
Y = 0: y0 = 0: der1 = 0: x = xold
RETURN

60 :
xi = .05
XPE = x * XN + XC: CIRCLE (XPE, YC), 8, 8
LOCATE 25, 5: PRINT USING "Error: ####.##   "; x
Errflag = 1
RESUME skip2

VARS:
'fill in the polynomial's coefficients here
'Don't forget to take out the apostrophe
'before the Y=C4*X^4+C3*X^3... equation
'note:CA/X, CB/X^2, etc. causes overflows near x=0
C4 = .2
C3 = 0
C2 = -1
C1 = 0
C0 = 2
CA = 0
CB = 0
CC = 0
RETURN

equation:
'=================================================================================================
'Use LOG10(X) instead of LOG(X) for common logs!!!!
Y = e ^ x
'  y = C4 * X ^ 4 + C3 * X ^ 3 + C2 * X ^ 2 + C1 * X + C0 + CA / X + CB / X ^ 2 + CC / X ^ 3
'          4th         3rd         2nd         1st       0th      -1th        -2nd       -3rd
'=================================================================================================
RETURN

FUNCTION acos (x)
'0<=y<=pi
IF x < 0 THEN
      acos = ATN(SQR(1 - x * x) / x) + pi
ELSEIF x = 0 THEN
      acos = pi / 2
ELSE
      acos = ATN(SQR(1 - x * x) / x) '(normal)
END IF
END FUNCTION

FUNCTION acosh (x)
'x >= 1
acosh = ln(x + SQR(x ^ 2 - 1))
END FUNCTION

FUNCTION acoth (x)
'³x³ > 1
acoth = .5 * ln((x + 1) / (x - 1))
END FUNCTION

FUNCTION acsch (x)
'x <> 0
acsch = ln(1 / x + SQR(1 + x ^ 2) / ABS(x))
END FUNCTION

FUNCTION asech (x)
'0 < x ó 1
asech = ln((1 + SQR(1 - x ^ 2)) / x)
IF x > 1 THEN END
END FUNCTION

FUNCTION asin (x)
asin = ATN(x / SQR(1 - x * x))
END FUNCTION

FUNCTION asinh (x)
asinh = ln(x + SQR(x ^ 2 + 1))
END FUNCTION

FUNCTION atan (x)
atan = ATN(x)
END FUNCTION

FUNCTION atanh (x)
'³x³ < 1
IF x >= 1 THEN END
atanh = .5 * ln((1 + x) / (1 - x))
END FUNCTION

FUNCTION cosh (x)
cosh = (e ^ x + e ^ -x) / 2
END FUNCTION

FUNCTION cot (x)
cot = 1 / TAN(x)
END FUNCTION

FUNCTION coth (x)
'undefined at x=0
coth = 1 / tanh(x)
END FUNCTION

FUNCTION csc (x)
csc = 1 / SIN(x)
END FUNCTION

FUNCTION csch (x)
'undefined at x=0
csch = 1 / sinh(x)
END FUNCTION

SUB eqs
'  y = (PI / 2 - X) * TAN(X)
'  y = sinh(x)
'  Y = atanh(x)
'  Y = tanh(X)
'  Y = (2 * x) / (SQR(x ^ 2 + x + 1))              'AP CALC AB 1995 #1B
'  Y = (X + 2) / (X ^ 2 + X + 1) ^ 1.5             '1st der #1B
'  Y = e ^ (-X ^ 2)                                      'Bell curve
'  Y = -2 * X * e ^ (-X ^ 2)                             '1st der
'  Y = 2 * X * e ^ (-X ^ 2) * (2 * X ^ 2 - 1)            '2nd der
'  Y = X * ln(X)
' --above are test questions from ch7
'  Y = (1 / X ^ 2) ^ X
'  Y = (2 ^ (COS(X) - 2)) / X
'  Y = (3 ^ (SIN(X) - 1)) / X
'  Y = X ^ (X + 1)
'  Y = X ^ ln(X)
'  Y = X ^ (-SQR(3))
'  Y = X ^ (SQR(2))
'  Y = (1 + 1 / X) ^ X        'e AS X->infinity
'  Y = X ^ (1 / X)
'  Y = X ^ (1 / ln(X))
'  Y = (ln(X)) ^ (-2)
'  Y = 1 / x - 1 / SQR(x)               '\________
'  Y = X * (1 - COS(X)) / (X - SIN(X))  '\/~~~^~~~\/
'  Y = (ln(X)) ^ 2
'  Y = e ^ x - x + 3 * SIN(5 * x) - 2
'  Y = 1 - 2 * COS(x) ^ 3

END SUB

SUB grid (x1, y1, XS, ys, XC, YC, XN, yn, xm, ym)

FOR x1 = XC TO XS STEP XN: LINE (x1, 0)-(x1, ys), 10, , &HAAAA: NEXT x1: REM Vertical down
FOR x1 = XC TO 0 STEP -XN: LINE (x1, 0)-(x1, ys), 10, , &HAAAA: NEXT x1: REM Vertical up
FOR y1 = YC TO ys STEP yn: LINE (0, y1)-(XS, y1), 2, , &HAAAA: NEXT y1:  REM Horizontal right
FOR y1 = YC TO 0 STEP -yn: LINE (0, y1)-(XS, y1), 2, , &HAAAA: NEXT y1:  REM Horizontal left

CN = -1
FOR x1 = XC TO XS STEP XN
        CN = CN + 1: IF CN MOD 5 = 0 THEN CL = 12 ELSE CL = 14
        LINE (x1, YC)-(x1, YC + 3), CL
        LINE (x1 + 1, YC)-(x1 + 1, YC + 3), CL
NEXT x1
CN = -1
FOR x1 = XC TO 0 STEP -XN
        CN = CN + 1: IF CN MOD 5 = 0 THEN CL = 12 ELSE CL = 14
        LINE (x1, YC)-(x1, YC + 3), CL
        LINE (x1 + 1, YC)-(x1 + 1, YC + 3), CL
NEXT x1
CN = -1
FOR y1 = YC TO ys STEP yn
        CN = CN + 1: IF CN MOD 5 = 0 THEN CL = 12 ELSE CL = 14
        LINE (XC, y1)-(XC + 4, y1), CL
NEXT y1
CN = -1
FOR y1 = YC TO 0 STEP -yn
        CN = CN + 1: IF CN MOD 5 = 0 THEN CL = 12 ELSE CL = 14
        LINE (XC, y1)-(XC + 4, y1), CL
NEXT y1
COLOR 15: LOCATE 30, 1
'PRINT "± KEYBOARD GUIDE:  [*] for WARP, [+] for FAST, [-] FOR SLOW, [Escape] to stop. ±";
LOCATE 29, 60: PRINT "[" + LTRIM$(STR$(-xm \ 2)) + "," + LTRIM$(STR$(-ym \ 2)) + "]";
PRINT "x[" + LTRIM$(STR$(xm \ 2)) + "," + LTRIM$(STR$(ym \ 2)) + "]";

LINE (XC, 0)-(XC, ys), 15: LINE (0, YC)-(XS, YC), 15:          REM Center

radianmarks XC, YC, XS, ys, XN

END SUB

FUNCTION ln (x)
'the LOG function in QBASIC returns natural log for some odd reason...
ln = LOG(x)
END FUNCTION

FUNCTION log10 (x)
log10 = LOG(x) / l
END FUNCTION

SUB radianmarks (XC, YC, XS, ys, XN)
CL1 = 9: CL2 = 4
CN = -1: FOR x1 = XC TO XS STEP (XN * pi / 2)
        CN = CN + 1: IF CN MOD 2 = 0 THEN CL = CL1 ELSE CL = CL2
        LINE (x1, YC)-(x1, YC - 3), CL
        DRAW "D2R1L2R1"
        'LINE (X1 + 1, YC)-(X1 + 1, YC - 3), CL
NEXT x1
CN = -1: FOR x1 = XC TO 0 STEP (-XN * pi / 2)
        CN = CN + 1: IF CN MOD 2 = 0 THEN CL = CL1 ELSE CL = CL2
        LINE (x1, YC)-(x1, YC - 3), CL
        DRAW "D2R1L2R1"
        'LINE (X1 - 1, YC)-(X1 + 1, YC - 3), CL
NEXT x1
END SUB

FUNCTION sec (x)
sec = 1 / COS(x)
END FUNCTION

FUNCTION sech (x)
sech = 1 / cosh(x)
END FUNCTION

FUNCTION sinh (x)
sinh = (e ^ x - e ^ -x) / 2
END FUNCTION

FUNCTION tanh (x)
tanh = sinh(x) / cosh(x)
END FUNCTION