home *** CD-ROM | disk | FTP | other *** search
/ Frozen Fish 1: Amiga / FrozenFish-Apr94.iso / bbs / alib / d3xx / d386 / xlispstat.lha / XLispStat / src2.lzh / XLisp-Stat / commonmath.c < prev    next >
Encoding:
C/C++ Source or Header  |  1990-10-09  |  22.7 KB  |  898 lines
  1. /* commonmath - xlisp math functions modified and augmented to         */
  2. /* correspond more closely to Common Lisp standard                     */
  3. /* XLISP-STAT 2.1 Copyright (c) 1990, by Luke Tierney                  */
  4. /* Additions to Xlisp 2.1, Copyright (c) 1989 by David Michael Betz    */
  5. /* You may give out copies of this software; for conditions see the    */
  6. /* file COPYING included with this distribution.                       */
  7.  
  8. #include "xlisp.h"
  9. #include "osdef.h"
  10. #ifdef ANSI
  11. #include "xlproto.h"
  12. #include "xlsproto.h"
  13. #include "osproto.h"
  14. #else
  15. #include "xlfun.h"
  16. #include "xlsfun.h"
  17. #include "osfun.h"
  18. #endif ANSI
  19. #include "xlvar.h"
  20.  
  21. typedef struct {
  22.   int mode;
  23.   FIXTYPE val, crval, cival;
  24.   FLOTYPE fval, cfrval, cfival;
  25. } Number; /* (Note: definition differs from statistics.c JKL */
  26.  
  27. /* forward declarations */
  28. #ifdef ANSI
  29. LVAL readnumber(Number *),lispnumber(Number *),unary(int),predicate(int),
  30.      compare(int);
  31. void setmode(Number *,int),matchmodes(Number *,Number *),
  32.      get_rem_arg(FLOTYPE *,int *),badiop(void),badfop(void),badcop(void);
  33. #else
  34. LVAL readnumber(),lispnumber(),unary(),predicate(),
  35.      compare();
  36. void setmode(),matchmodes(),
  37.      get_rem_arg(),badiop(),badfop(),badcop();
  38. #endif ANSI
  39.  
  40. #define IN 0
  41. #define FL 1
  42. #define CI 2
  43. #define CF 3
  44.  
  45. /* added JKL */
  46. #define badarg xlbadtype
  47.  
  48. /* binary functions */
  49. LVAL xadd()    { return (binary('+')); } /* + */
  50. LVAL xsub()    { return (binary('-')); } /* - */
  51. LVAL xmul()    { return (binary('*')); } /* * */
  52. LVAL xdiv()    { return (binary('/')); } /* / */
  53. LVAL xmin()    { return (binary('m')); } /* min */
  54. LVAL xmax()    { return (binary('M')); } /* max */
  55.  
  56. LVAL xlogand() { return (logbinary('&')); } /* logand */
  57. LVAL xlogior() { return (logbinary('|')); } /* logior */
  58. LVAL xlogxor() { return (logbinary('^')); } /* logxor */
  59.  
  60. static LVAL readnumber(num)
  61.     Number *num;
  62. {
  63.   LVAL arg = xlgetarg(), real, imag;
  64.   
  65.   if (fixp(arg)) {
  66.     num->mode = IN;
  67.     num->val = getfixnum(arg);
  68.   }
  69.   else if (floatp(arg)) {
  70.     num->mode = FL;
  71.     num->fval = getflonum(arg);
  72.   }
  73.   else if (complexp(arg)) {
  74.     real = realpart(arg);
  75.     imag = imagpart(arg);
  76.     if (fixp(real)) {
  77.       num->mode = CI;
  78.       num->crval = getfixnum(real);
  79.       num->cival = getfixnum(imag);
  80.     }
  81.     else {
  82.       num->mode = CF;
  83.       num->cfrval = makedouble(real);
  84.       num->cfival = makedouble(imag);
  85.     }
  86.   }
  87.   else xlerror("not a number", arg);
  88.   return(arg);
  89. }
  90.  
  91. static void setmode(x, mode)
  92.     Number *x;
  93.     int mode;
  94. {
  95.   switch (mode) {
  96.   case FL:
  97.     if (x->mode != IN) return;
  98.     x->mode = mode;
  99.     x->fval = x->val;
  100.     break;
  101.   case CI:
  102.     if (x->mode != IN) return;
  103.     x->mode = mode;
  104.     x->crval = x->val;
  105.     x->cival = 0;
  106.     break;
  107.   case CF:
  108.     switch (x->mode) {
  109.     case IN:
  110.       x->mode = mode;
  111.       x->cfrval = x->val;
  112.       x->cfival = 0.0;
  113.       break;
  114.     case FL:
  115.       x->mode = mode;
  116.       x->cfrval = x->fval;
  117.       x->cfival = 0.0;
  118.       break;
  119.     case CI:
  120.       x->mode = mode;
  121.       x->cfrval = x->crval;
  122.       x->cfival = x->cival;
  123.       break;
  124.     }
  125.     break;
  126.   }
  127. }
  128.  
  129. static void matchmodes(x, y)
  130.     Number *x, *y;
  131. {
  132.   int mode = x->mode;
  133.   switch (mode) {
  134.   case IN: mode = y->mode; break;
  135.   case FL: if (y->mode == CI || y->mode == CF) mode = CF; break;
  136.   case CI: if (y->mode == FL || y->mode == CF) mode = CF; break;
  137.   case CF: break;
  138.   }
  139.   if (x->mode != mode) setmode(x, mode);
  140.   if (y->mode != mode) setmode(y, mode);
  141. }
  142.  
  143. static LVAL lispnumber(x)
  144.     Number *x;
  145. {
  146.   switch (x->mode) {
  147.   case IN: return(cvfixnum(x->val));
  148.   case FL: return(cvflonum(x->fval));
  149.   case CI: return(newicomplex(x->crval, x->cival));
  150.   case CF: return(newdcomplex(x->cfrval, x->cfival));
  151.   }
  152. }
  153. /* see define above JKL
  154. static void badarg(arg)
  155.     LVAL arg;
  156. {
  157.   xlerror("bad argument type", arg);
  158. }
  159. */
  160. LVAL binary(which)
  161.     int which;
  162. {
  163.   LVAL larg;
  164.   Number val, arg;
  165.   FIXTYPE rtemp, itemp;
  166.   FLOTYPE frtemp, fitemp, magn;
  167.  
  168.   if (xlargc == 1 && (which == '-' || which == '/')) {
  169.     val.mode = IN;
  170.     switch (which) {
  171.     case '-': val.val = 0; break;
  172.     case '/': val.val = 1; break;
  173.     }
  174.   }
  175.   else /*larg =*/ readnumber(&val); /* not used JKL */
  176.   while (moreargs()) {
  177.     larg = readnumber(&arg);
  178.     matchmodes(&val, &arg);
  179.     switch (which) {
  180.     case '+':
  181.       switch (val.mode) {
  182.       case IN: val.val   += arg.val;  break;
  183.       case FL: val.fval  += arg.fval; break;
  184.       case CI: val.crval += arg.crval;   val.cival += arg.cival;   break;
  185.       case CF: val.cfrval += arg.cfrval; val.cfival += arg.cfival; break;
  186.       }
  187.       break;
  188.     case '-':
  189.       switch (val.mode) {
  190.       case IN: val.val   -= arg.val;  break;
  191.       case FL: val.fval  -= arg.fval; break;
  192.       case CI: val.crval -= arg.crval;   val.cival -= arg.cival;   break;
  193.       case CF: val.cfrval -= arg.cfrval; val.cfival -= arg.cfival; break;
  194.       }
  195.       break;
  196.     case '*':
  197.       switch (val.mode) {
  198.       case IN: val.val   *= arg.val;  break;
  199.       case FL: val.fval  *= arg.fval; break;
  200.       case CI:
  201.         rtemp = val.crval * arg.crval - val.cival * arg.cival;
  202.         itemp = val.cival * arg.crval + val.crval * arg.cival;
  203.         val.crval = rtemp; val.cival = itemp;
  204.         break;
  205.       case CF:
  206.         frtemp = val.cfrval * arg.cfrval - val.cfival * arg.cfival;
  207.         fitemp = val.cfival * arg.cfrval + val.cfrval * arg.cfival;
  208.         val.cfrval = frtemp; val.cfival = fitemp;
  209.         break;
  210.       }
  211.       break;
  212.     case '/':
  213.       switch (val.mode) {
  214.       case IN:
  215.         checkizero(arg.val);
  216.         if (val.val % arg.val == 0) {
  217.           val.val /= arg.val;
  218.           break;
  219.         }
  220.         else {
  221.           setmode(&val, FL);
  222.           setmode(&arg, FL);
  223.         }
  224.         /* drop through */
  225.       case FL:
  226.         checkfzero(arg.fval);
  227.         val.fval /= arg.fval;
  228.         break;
  229.       case CI:
  230.         setmode(&val, CF);
  231.         setmode(&arg, CF);
  232.         /* drop through */
  233.       case CF:
  234.         magn = arg.cfrval * arg.cfrval + arg.cfival * arg.cfival;
  235.         checkfzero(magn);
  236.         frtemp = (val.cfrval * arg.cfrval + val.cfival * arg.cfival) / magn;
  237.         fitemp = (val.cfival * arg.cfrval - val.cfrval * arg.cfival) / magn;
  238.         val.cfrval = frtemp; val.cfival = fitemp;
  239.         break;
  240.       }
  241.       break;
  242.     case 'M':
  243.       switch (val.mode) {
  244.       case IN: val.val  = (val.val > arg.val)   ? val.val  : arg.val;  break;
  245.       case FL: val.fval = (val.fval > arg.fval) ? val.fval : arg.fval; break;
  246.       default: badarg(larg);
  247.       }
  248.       break;
  249.     case 'm':
  250.       switch (val.mode) {
  251.       case IN: val.val  = (val.val < arg.val)   ? val.val  : arg.val;  break;
  252.       case FL: val.fval = (val.fval < arg.fval) ? val.fval : arg.fval; break;
  253.       default: badarg(larg);
  254.       }
  255.       break;
  256.     }
  257.   }
  258.   return(lispnumber(&val));
  259. }
  260.  
  261. static void get_rem_arg(fval,mode)
  262.     FLOTYPE *fval;
  263.     int *mode;
  264. {
  265.   LVAL arg;
  266.   
  267.   arg = xlgetarg();
  268.   if (fixp(arg)) {
  269.     *fval = getfixnum(arg);
  270.     *mode = IN;
  271.   }
  272.   else if (floatp(arg)) {
  273.     *fval = getflonum(arg);
  274.     *mode = FL;
  275.   }
  276.   else badarg(arg);
  277. }
  278.  
  279. LVAL xrem()
  280. {
  281.   int mode1, mode2;
  282.   FLOTYPE fval1, fval2, fres;
  283.   
  284.   get_rem_arg(&fval1, &mode1);
  285.   get_rem_arg(&fval2, &mode2);
  286.   xllastarg();
  287.  
  288.   fres = fval1 - fval2 * floor(fval1 / fval2);
  289.   return((mode1 == IN && mode2 == IN) ? cvfixnum((FIXTYPE) fres)
  290.                                       : cvflonum((FLOTYPE) fres));
  291. }
  292.  
  293. LVAL logbinary(which)
  294.     int which;
  295. {
  296.   int val, arg;
  297.   
  298.   switch (which) {
  299.   case '&': val = -1; break;
  300.   case '|': val =  0; break;
  301.   case '^': val =  0; break;
  302.   }
  303.   while (moreargs()) {
  304.     arg = getfixnum(xlgafixnum());
  305.     switch (which) {
  306.     case '&': val &= arg; break;
  307.     case '|': val |= arg; break;
  308.     case '^': val ^= arg; break;
  309.     }
  310.   }
  311.   return(cvfixnum((FIXTYPE) val));
  312. }
  313.  
  314. LVAL xexpt()
  315. {
  316.   LVAL base, power;
  317.   int bsign, psign;
  318.   FIXTYPE b, p, val;
  319.   FLOTYPE fb, fp, fval;
  320.   
  321.   base = xlgetarg();
  322.   power = xlgetarg();
  323.   xllastarg();
  324.   
  325.   if (fixp(base) && fixp(power)) {
  326.     b = getfixnum(base);
  327.     p = getfixnum(power);
  328.     if (p == 0) return(cvfixnum((FIXTYPE) 1));
  329.     if (b == 0 && p > 0) return(cvfixnum((FIXTYPE) 0));
  330.     checkizero(b);
  331.     bsign = (b > 0) ? 1 : -1;
  332.     psign = (p > 0) ? 1 : -1;
  333.     b = (b > 0) ? b : -b;
  334.     p = (p > 0) ? p : -p;
  335.     fval = floor(f_expt((double) b, (double) p) + 0.1); /* to get integer right */
  336.     if (bsign == -1 && p % 2 == 1) fval = -fval;
  337.     if (psign == 1) {
  338.       val = fval;
  339.       if (val == fval) return(cvfixnum((FIXTYPE) val));
  340.       else return(cvflonum((FLOTYPE) fval));    /* to handle precision for large results */
  341.     }
  342.     else {
  343.       checkfzero(fval);
  344.       return(cvflonum((FLOTYPE) 1.0 / fval));
  345.     }
  346.   } 
  347.   else if (floatp(base) && fixp(power)) {
  348.     fb = getflonum(base);
  349.     p = getfixnum(power);
  350.     if (p == 0) return(cvfixnum((FIXTYPE) 1));
  351.     if (fb == 0.0 && p > 0) return(cvflonum((FLOTYPE) 0.0));
  352.     checkfzero(fb);
  353.     bsign = (fb > 0) ? 1 : -1;
  354.     psign = (p > 0) ? 1 : -1;
  355.     fb = (fb > 0) ? fb : -fb;
  356.     p = (p > 0) ? p : -p;
  357.     fval = f_expt((double) fb, (double) p);
  358.     if (bsign == -1 && p % 2 == 1) fval = -fval;
  359.     if (psign == 1) return(cvflonum((FLOTYPE) fval));
  360.     else {
  361.       checkfzero(fval);
  362.       return(cvflonum((FLOTYPE) 1.0 / fval));
  363.     }
  364.   }
  365.   else if ((fixp(base) || floatp(base)) && floatp(power)) {
  366.     fb = makedouble(base);
  367.     fp = getflonum(power);
  368.     if (fp == 0.0) return(cvflonum((FLOTYPE) 1.0));
  369.     if (fb == 0.0 && fp > 0.0) return(cvflonum((FLOTYPE) 0.0));
  370.     if (fb < 0.0)
  371.       return(cvcomplex(cexpt(makecomplex(base), makecomplex(power))));
  372.     psign = (fp > 0) ? 1 : -1;
  373.     fb = (fb > 0) ? fb : -fb;
  374.     fp = (fp > 0) ? fp : -fp;
  375.     fval = f_expt((double) fb, (double) fp);
  376.     if (psign == 1) return(cvflonum((FLOTYPE) fval));
  377.     else {
  378.       checkfzero(fval);
  379.       return(cvflonum((FLOTYPE) 1.0 / fval));
  380.     }    
  381.   }
  382.   else if (complexp(base) || complexp(power))
  383.     return(cvcomplex(cexpt(makecomplex(base), makecomplex(power))));
  384.   else xlfail("bad argument type(s)");
  385. }
  386.  
  387. LVAL xlog()
  388. {
  389.   LVAL arg, base;
  390.   int base_supplied = FALSE;
  391.   double fx, fb;
  392.   
  393.   arg = xlgetarg();
  394.   if (moreargs()) {
  395.     base_supplied = TRUE;
  396.     base = xlgetarg();
  397.   }
  398.   if (base_supplied) {
  399.     if (numberp(arg) && numberp(base)) {
  400.       fx = makedouble(arg);
  401.       fb = makedouble(base);
  402.       if (fx <= 0.0 || fb <= 0.0)
  403.         return(cvcomplex(cdiv(clog(makecomplex(arg)), clog(makecomplex(base)))));
  404.       else return(cvflonum((FLOTYPE) logarithm(fx, fb, TRUE)));
  405.     }
  406.     else if ((numberp(arg) && complexp(base))
  407.              || (complexp(arg) && numberp(base))
  408.              || (complexp(arg) && complexp(base)))
  409.       return(cvcomplex(cdiv(clog(makecomplex(arg)), clog(makecomplex(base)))));
  410.     else xlfail("bad argument type(s)");
  411.   }
  412.   else {
  413.     if (numberp(arg)) {
  414.       fx = makedouble(arg);
  415.       if (fx <= 0.0) return(cvcomplex(clog(makecomplex(arg))));
  416.       else return(cvflonum((FLOTYPE) logarithm(fx, 0.0, FALSE)));
  417.     }
  418.     else if (complexp(arg)) 
  419.       return(cvcomplex(clog(makecomplex(arg))));
  420.     else xlfail("bad argument type(s)");
  421.   }
  422. }
  423.  
  424. /* xgcd - greatest common divisor */
  425. LVAL xgcd()
  426. {
  427.   FIXTYPE m,n,r;
  428.   LVAL arg;
  429.  
  430.   if (!moreargs())                  /* check for identity case */
  431.     return (cvfixnum((FIXTYPE)0));
  432.   arg = xlgafixnum();
  433.   n = getfixnum(arg);
  434.   if (n < (FIXTYPE)0) n = -n;        /* absolute value */
  435.   while (moreargs()) {
  436.     arg = xlgafixnum();
  437.     m = getfixnum(arg);
  438.     if (m == 0 || n == 0) xlfail("zero argument");
  439.     if (m < (FIXTYPE)0) m = -m;        /* absolute value */
  440.     for (;;) {                      /* euclid's algorithm */
  441.       r = m % n;
  442.       if (r == (FIXTYPE) 0)
  443.         break;
  444.       m = n;
  445.       n = r;
  446.     }
  447.   }
  448.   return (cvfixnum(n));
  449. }
  450.  
  451. /* checkizero - check for integer division by zero */
  452. void checkizero(iarg)
  453.   FIXTYPE iarg;
  454. {
  455.   if (iarg == 0)
  456.   xlfail("illegal zero argument");
  457. }
  458.  
  459. /* checkfzero - check for floating point division by zero or log of zero */
  460. void checkfzero(farg)
  461.   FLOTYPE farg;
  462. {
  463.   if (farg == 0.0)
  464.   xlfail("illegal zero argument");
  465. }
  466.  
  467. /* unary functions */
  468. LVAL xlognot() { return (unary('~')); } /* lognot */
  469. LVAL xabs()    { return (unary('A')); } /* abs */
  470. LVAL xadd1()   { return (unary('+')); } /* 1+ */
  471. LVAL xsub1()   { return (unary('-')); } /* 1- */
  472. LVAL xsin()    { return (unary('S')); } /* sin */
  473. LVAL xcos()    { return (unary('C')); } /* cos */
  474. LVAL xtan()    { return (unary('T')); } /* tan */
  475. LVAL xexp()    { return (unary('E')); } /* exp */
  476. LVAL xsqrt()   { return (unary('R')); } /* sqrt */
  477. LVAL xfix()    { return (unary('I')); } /* truncate */
  478. LVAL xfloat()  { return (unary('F')); } /* float */
  479. LVAL xrand()   { return (unary('?')); } /* random */
  480. LVAL xfloor()  { return (unary('_')); } /* floor */
  481. LVAL xceil()   { return (unary('^')); } /* ceiling */
  482. LVAL xround()  { return (unary('r')); } /* round */
  483. LVAL xasin()   { return (unary('s')); } /* asin */
  484. LVAL xacos()   { return (unary('c')); } /* acos */
  485. LVAL xatan()   { return (unary('t')); } /* atan */
  486. LVAL xphase()  { return (unary('P')); } /* phase */
  487.  
  488. /* unary - handle unary operations */
  489. LOCAL LVAL unary(which)
  490.     int which;
  491. {
  492.   FLOTYPE fval;
  493.   FIXTYPE ival;
  494.   Complex cval;
  495.   LVAL arg, real, imag;
  496.   int mode;
  497.   
  498.   /* get the argument */
  499.   arg = xlgetarg();
  500.   if (which == 'F' && moreargs()) xlgaflonum();
  501.   xllastarg();
  502.  
  503.   /* check its type */
  504.   if (fixp(arg)) {
  505.     ival = getfixnum(arg);
  506.     mode = IN;
  507.   }
  508.   else if (floatp(arg)) {
  509.     fval = getflonum(arg);
  510.     mode = FL;
  511.   }
  512.   else if (complexp(arg)) {
  513.     cval = makecomplex(arg);
  514.     real = realpart(arg);
  515.     imag = imagpart(arg);
  516.     if (fixp(realpart(arg))) mode = CI;
  517.     else mode = CF;
  518.   }
  519.   else xlerror("not a number", arg);
  520.  
  521.   switch (which) {
  522.   case '~':
  523.     if (mode == IN) return(cvfixnum((FIXTYPE) ~ival));    
  524.     else badiop();
  525.     break;
  526.   case 'A':
  527.     switch (mode) {
  528.     case IN: return(cvfixnum((FIXTYPE) (ival < 0   ? -ival : ival)));
  529.     case FL: return(cvflonum((FLOTYPE) (fval < 0.0 ? -fval : fval)));
  530.     case CI:
  531.     case CF: return(cvflonum((FLOTYPE) modulus(cval)));
  532.     }
  533.     break;
  534.   case '+':
  535.     switch (mode) {
  536.     case IN: return(cvfixnum((FIXTYPE) ival + 1));
  537.     case FL: return(cvflonum((FLOTYPE) fval + 1.0));
  538.     case CI: return(newicomplex(getfixnum(real) + 1, getfixnum(imag)));
  539.     case CF: return(newdcomplex(getflonum(real) + 1.0, getflonum(imag)));
  540.     }
  541.     break;
  542.   case '-':
  543.     switch (mode) {
  544.     case IN: return(cvfixnum((FIXTYPE) ival - 1));
  545.     case FL: return(cvflonum((FLOTYPE) fval - 1.0));
  546.     case CI: return(newicomplex(getfixnum(real) - 1, getfixnum(imag)));
  547.     case CF: return(newdcomplex(getflonum(real) - 1.0, getflonum(imag)));
  548.     }
  549.     break;
  550.   case 'S':
  551.     switch (mode) {
  552.     case IN: return(cvflonum((FLOTYPE) sin((double) ival)));
  553.     case FL: return(cvflonum((FLOTYPE) sin((double) fval)));
  554.     case CI:
  555.     case CF: return(cvcomplex(csin(cval)));
  556.     }
  557.   case 'C':
  558.     switch (mode) {
  559.     case IN: return(cvflonum((FLOTYPE) cos((double) ival)));
  560.     case FL: return(cvflonum((FLOTYPE) cos((double) fval)));
  561.     case CI:
  562.     case CF: return(cvcomplex(ccos(cval)));
  563.     }
  564.   case 'T': 
  565.     switch (mode) {
  566.     case IN: return(cvflonum((FLOTYPE) tan((double) ival)));
  567.     case FL: return(cvflonum((FLOTYPE) tan((double) fval)));
  568.     case CI:
  569.     case CF: return(cvcomplex(ctan(cval)));
  570.     }
  571.   case 'E':
  572.     switch (mode) {
  573.     case IN: return(cvflonum((FLOTYPE) f_exp((double) ival)));
  574.     case FL: return(cvflonum((FLOTYPE) f_exp((double) fval)));
  575.     case CI: 
  576.     case CF: return(cvcomplex(cexp(cval)));
  577.     }
  578.     break;
  579.   case 'R':
  580.     switch (mode) {
  581.     case IN:
  582.       if (ival < 0) return(cvcomplex(csqrt(makecomplex(arg))));
  583.       else return(cvflonum((FLOTYPE) f_sqrt((double) ival)));
  584.     case FL:
  585.       if (fval < 0) return(cvcomplex(csqrt(makecomplex(arg))));
  586.       else return(cvflonum((FLOTYPE) f_sqrt(fval)));
  587.     case CI:
  588.     case CF: return(cvcomplex(csqrt(cval)));
  589.     }
  590.     break;
  591.   case 'I':
  592.     switch (mode) {
  593.     case IN: return (cvfixnum((FIXTYPE) ival));
  594.     case FL: return (cvfixnum((FIXTYPE) fval));
  595.     default: badcop();
  596.     }
  597.     break;
  598.   case 'F': 
  599.     switch (mode) {
  600.     case IN: return (cvflonum((FLOTYPE) ival));
  601.     case FL: return (cvflonum((FLOTYPE) fval));
  602.     default: badcop();
  603.     }
  604.     break;
  605.   case '?':
  606.     switch (mode) {
  607.     case IN: 
  608.       if (ival <= 0) badiop();
  609.       return (cvfixnum((FIXTYPE) osrand((int) ival)));
  610.     case FL: 
  611.       if (fval <= 0.0) badfop();
  612.       return (cvflonum((FLOTYPE) fval * uni()));
  613.     default: badcop();
  614.     }
  615.     break;
  616.   case '_':
  617.     switch (mode) {
  618.     case IN: return (cvfixnum((FIXTYPE) ival));
  619.     case FL: return (cvfixnum((FIXTYPE) floor(fval)));
  620.     default: badcop();
  621.     }
  622.     break;
  623.   case '^':
  624.     switch (mode) {
  625.     case IN: return (cvfixnum((FIXTYPE) ival));
  626.     case FL: return (cvfixnum((FIXTYPE) ceil(fval)));
  627.     default: badcop();
  628.     }
  629.     break;
  630.   case 'r':
  631.     switch (mode) {
  632.     case IN: return (cvfixnum((FIXTYPE) ival));
  633.     case FL: return (cvfixnum((FIXTYPE) floor(fval + 0.5)));
  634.     default: badcop();
  635.     }
  636.     break;
  637.   case 's':
  638.     switch (mode) {
  639.     case IN:
  640.       fval = ival;
  641.       /* drop through */
  642.     case FL:
  643.       if (fval > 1.0 || fval < -1.0) 
  644.         return(cvcomplex(casin(makecomplex(arg))));
  645.       else return(cvflonum((FLOTYPE) asin(fval)));
  646.     case CI:
  647.     case CF: return(cvcomplex(casin(cval)));
  648.     }
  649.     break;
  650.   case 'c':
  651.     switch (mode) {
  652.     case IN:
  653.       fval = ival;
  654.       /* drop through */
  655.     case FL:
  656.       if (fval > 1.0 || fval < -1.0) 
  657.         return(cvcomplex(cacos(makecomplex(arg))));
  658.       else return(cvflonum((FLOTYPE) acos(fval)));
  659.     case CI:
  660.     case CF: return(cvcomplex(cacos(cval)));
  661.     }
  662.     break;
  663.   case 't':
  664.     switch (mode) {
  665.     case IN: fval = ival; /* drop through */
  666.     case FL: return(cvflonum((FLOTYPE) atan(fval)));
  667.     case CI:
  668.     case CF: return(cvcomplex(catan(cval)));
  669.     }
  670.     break;
  671.   case 'P':
  672.     switch (mode) {
  673.     case IN: return(cvflonum((FLOTYPE) (ival >= 0) ? 0.0 : PI));
  674.     case FL: return(cvflonum((FLOTYPE) (fval >= 0.0) ? 0.0 : PI));
  675.     case CI:
  676.     case CF: return(cvflonum((FLOTYPE) phase(cval)));
  677.     }
  678.     break;
  679.   default: xlfail("unsupported operation");
  680.   }
  681. }
  682.  
  683. /* unary predicates */
  684. LVAL xminusp() { return (predicate('-')); } /* minusp */
  685. LVAL xzerop()  { return (predicate('Z')); } /* zerop */
  686. LVAL xplusp()  { return (predicate('+')); } /* plusp */
  687. LVAL xevenp()  { return (predicate('E')); } /* evenp */
  688. LVAL xoddp()   { return (predicate('O')); } /* oddp */
  689.  
  690.  
  691. /* predicate - handle a predicate function */
  692. LOCAL LVAL predicate(fcn)
  693.   int fcn;
  694. {
  695.   FLOTYPE fval;
  696.   FIXTYPE ival;
  697.   LVAL arg;
  698.  
  699.   /* get the argument */
  700.   arg = xlgetarg();
  701.   xllastarg();
  702.  
  703.   /* check the argument type */
  704.   if (fixp(arg)) {
  705.     ival = getfixnum(arg);
  706.     switch (fcn) {
  707.     case '-': ival = (ival < 0); break;
  708.     case 'Z': ival = (ival == 0); break;
  709.     case '+': ival = (ival > 0); break;
  710.     case 'E': ival = ((ival & 1) == 0); break;
  711.     case 'O': ival = ((ival & 1) != 0); break;
  712.     default:  badiop();
  713.     }
  714.   }
  715.   else if (floatp(arg)) {
  716.     fval = getflonum(arg);
  717.     switch (fcn) {
  718.     case '-': ival = (fval < 0); break;
  719.     case 'Z': ival = (fval == 0); break;
  720.     case '+': ival = (fval > 0); break;
  721.     default:  badfop();
  722.     }
  723.   }
  724.   else
  725.     badarg(arg);
  726.  
  727.   /* return the result value */
  728.   return (ival ? true : NIL);
  729. }
  730.  
  731. /* comparison functions */
  732. LVAL xlss() { return (compare('<'));  } /* < */
  733. LVAL xleq() { return (compare('L'));  } /* <= */
  734. LVAL xequ() { return (ccompare('=')); } /* = */
  735. LVAL xneq() { return (ccompare('#')); } /* /= */
  736. LVAL xgeq() { return (compare('G'));  } /* >= */
  737. LVAL xgtr() { return (compare('>'));  } /* > */
  738.  
  739. /* compare - common compare function */
  740. LOCAL LVAL compare(fcn)
  741.   int fcn;
  742. {
  743.     FIXTYPE icmp,ival,iarg;
  744.     FLOTYPE fcmp,fval,farg;
  745.     LVAL arg;
  746.     int mode;
  747.  
  748.     /* get the first argument */
  749.     arg = xlgetarg();
  750.  
  751.     /* set the type of the first argument */
  752.     if (fixp(arg)) {
  753.     ival = getfixnum(arg);
  754.     mode = 'I';
  755.     }
  756.     else if (floatp(arg)) {
  757.     fval = getflonum(arg);
  758.     mode = 'F';
  759.     }
  760.     else
  761.       badarg(arg);
  762.  
  763.     /* handle each remaining argument */
  764.     for (icmp = TRUE; icmp && moreargs(); ival = iarg, fval = farg) {
  765.  
  766.     /* get the next argument */
  767.     arg = xlgetarg();
  768.  
  769.     /* check its type */
  770.     if (fixp(arg)) {
  771.         switch (mode) {
  772.         case 'I':
  773.             iarg = getfixnum(arg);
  774.             break;
  775.         case 'F':
  776.             farg = (FLOTYPE)getfixnum(arg);
  777.         break;
  778.         }
  779.     }
  780.     else if (floatp(arg)) {
  781.         switch (mode) {
  782.         case 'I':
  783.             fval = (FLOTYPE)ival;
  784.         farg = getflonum(arg);
  785.         mode = 'F';
  786.         break;
  787.         case 'F':
  788.             farg = getflonum(arg);
  789.         break;
  790.         }
  791.     }
  792.     else
  793.       badarg(arg);
  794.  
  795.     /* compute result of the compare */
  796.     switch (mode) {
  797.     case 'I':
  798.         icmp = ival - iarg;
  799.         switch (fcn) {
  800.         case '<':    icmp = (icmp < 0); break;
  801.         case 'L':    icmp = (icmp <= 0); break;
  802.         case '=':    icmp = (icmp == 0); break;
  803.         case '#':    icmp = (icmp != 0); break;
  804.         case 'G':    icmp = (icmp >= 0); break;
  805.         case '>':    icmp = (icmp > 0); break;
  806.         }
  807.         break;
  808.     case 'F':
  809.         fcmp = fval - farg;
  810.         switch (fcn) {
  811.         case '<':    icmp = (fcmp < 0.0); break;
  812.         case 'L':    icmp = (fcmp <= 0.0); break;
  813.         case '=':    icmp = (fcmp == 0.0); break;
  814.         case '#':    icmp = (fcmp != 0.0); break;
  815.         case 'G':    icmp = (fcmp >= 0.0); break;
  816.         case '>':    icmp = (fcmp > 0.0); break;
  817.         }
  818.         break;
  819.     }
  820.     }
  821.  
  822.     /* return the result */
  823.     return (icmp ? true : NIL);
  824. }
  825.  
  826. LVAL ccompare(which)
  827.     int which;
  828. {
  829. /*LVAL larg; */
  830.   Number val, arg;
  831.   int icmp;
  832.   
  833.   switch (which) {
  834.   case '=': icmp = TRUE;  break;
  835.   case '#': icmp = FALSE; break;
  836.   }
  837. /*larg = */ readnumber(&val); /* not used JKL */
  838.   while (moreargs()) {
  839. /*  larg = */ readnumber(&arg); /* not used JKL */
  840.     matchmodes(&val, &arg);
  841.     switch (which) {
  842.     case '=':
  843.       switch (val.mode) {
  844.       case IN: icmp = icmp && val.val  == arg.val;  break;
  845.       case FL: icmp = icmp && val.fval == arg.fval; break;
  846.       case CI: icmp = icmp && val.crval == arg.crval && val.cival == arg.cival; break;
  847.       case CF: icmp = icmp && val.cfrval == arg.cfrval && val.cfival == arg.cfival; break;
  848.       }
  849.       break;
  850.     case '#':
  851.       switch (val.mode) {
  852.       case IN: icmp = icmp || val.val  != arg.val;  break;
  853.       case FL: icmp = icmp || val.fval != arg.fval; break;
  854.       case CI: icmp = icmp || val.crval != arg.crval || val.cival != arg.cival; break;
  855.       case CF: icmp = icmp || val.cfrval != arg.cfrval || val.cfival != arg.cfival; break;
  856.       }
  857.       break;
  858.     }
  859.   }
  860.   return((icmp) ? true : NIL);
  861. }
  862.  
  863. /* badiop - bad integer operation */
  864. LOCAL void badiop()
  865. {
  866.   xlfail("bad integer operation");
  867. }
  868.  
  869. /* badfop - bad floating point operation */
  870. LOCAL void badfop()
  871. {
  872.   xlfail("bad floating point operation");
  873. }
  874.  
  875. /* badcop - bad complex number operation */
  876. LOCAL void badcop()
  877. {
  878.   xlfail("bad complex number operation");
  879. }
  880.  
  881. /* two argument logarithm */
  882. double logarithm(x, base, base_supplied)
  883.      FLOTYPE x, base;
  884.      int base_supplied;
  885. {
  886.   double lbase;
  887.   if (x <= 0.0) xlfail("logarithm of a nonpositive number");
  888.   if (base_supplied) {
  889.     if (base <= 0.0) xlfail("logarithm to a nonpositive base");
  890.     else {
  891.       lbase = f_log(base);
  892.       if (lbase == 0.0) xlfail("logarith to a unit base");
  893.       else return((f_log(x)/lbase));
  894.     }
  895.   }
  896.   else return (f_log(x));
  897. }
  898.