home *** CD-ROM | disk | FTP | other *** search
/ Usenet 1994 January / usenetsourcesnewsgroupsinfomagicjanuary1994.iso / sources / unix / volume10 / ifp / part03 / interp / F_arith.c next >
Encoding:
C/C++ Source or Header  |  1987-07-05  |  9.0 KB  |  378 lines
  1.  
  2. /****** F_arith.c *****************************************************/
  3. /**                                                                  **/
  4. /**                    University of Illinois                        **/
  5. /**                                                                  **/
  6. /**                Department of Computer Science                    **/
  7. /**                                                                  **/
  8. /**   Tool: IFP                         Version: 0.5                 **/
  9. /**                                                                  **/
  10. /**   Author:  Arch D. Robison          Date:   May 1, 1985          **/
  11. /**                                                                  **/
  12. /**   Revised by: Arch D. Robison       Date:  June 4, 1986          **/
  13. /**                                                                  **/
  14. /**   Principal Investigators: Prof. R. H. Campbell                  **/
  15. /**                            Prof. W. J. Kubitz                    **/
  16. /**                                                                  **/
  17. /**                                                                  **/
  18. /**------------------------------------------------------------------**/
  19. /**   (C) Copyright 1987  University of Illinois Board of Trustees   **/
  20. /**                       All Rights Reserved.                       **/
  21. /**********************************************************************/
  22.  
  23. #include <stdio.h>
  24. #include <math.h>
  25. #include <errno.h>
  26. #include "struct.h"
  27. #include "node.h"
  28.  
  29. #if OPSYS!=CTSS
  30. extern int errno;      /* exists somewhere in UNIX */
  31. #endif
  32.  
  33. /* NOTE - function Dyadic assumes integers are in two's complement form! */
  34.  
  35. private F_Minus (), F_AddN (), Monadic (), Dyadic (), F_Sum ();
  36.  
  37. private OpDef OpArith [] = {
  38. #if OPSYS!=CTSS
  39.    {"ln",       0,      Monadic},
  40.    {"exp",      1,      Monadic},
  41.    {"sqrt",     2,      Monadic},
  42.    {"sin",      3,      Monadic},
  43.    {"cos",      4,      Monadic},
  44.    {"tan",      5,      Monadic},
  45.    {"arcsin",   6,      Monadic},
  46.    {"arccos",   7,      Monadic},
  47.    {"arctan",   8,      Monadic},
  48. #endif
  49.    {"minus",    -1,     F_Minus},
  50.    {"add1",     1,      F_AddN},
  51.    {"sub1",     -1,     F_AddN},
  52.    {"+",        0,      Dyadic},
  53.    {"-",        1,      Dyadic},
  54.    {"*",        2,      Dyadic},
  55.    {"%",        3,      Dyadic},
  56. #if OPSYS!=CTSS
  57.    {"mod",      4,      Dyadic},
  58.    {"div",      5,      Dyadic},
  59. #endif
  60.    {"min",      6,      Dyadic},
  61.    {"max",      7,      Dyadic},
  62. #if OPSYS!=CTSS
  63.    {"power",    8,      Dyadic},
  64. #endif
  65.    {"sum",      -1,     F_Sum}
  66. };
  67.  
  68.  
  69. /*
  70.  * Monadic
  71.  *
  72.  * Evaluate a monadic function
  73.  *
  74.  * Input
  75.  *      InOut = argument to apply function
  76.  *      Op = operation - see array F_Name in code for values
  77.  *
  78.  * Output
  79.  *      InOut = result of applying function
  80.  */
  81. private Monadic (InOut,Op)
  82.    ObjectPtr InOut;
  83.    int Op;
  84.    {
  85.       double X,Z;
  86.       register int E;
  87.  
  88.       if (GetDouble (InOut,&X)) FunError ("not numeric",InOut);
  89.       else {
  90.      E = 0;
  91.      switch (Op) {
  92. #if OPSYS!=CTSS
  93.         case 0:                     /* base e log */
  94.            if (X <= 0) E = EDOM;
  95.            else Z = log (X);
  96.            break;
  97.         case 1:                     /* base e power */
  98.            if (X >= LNMAXFLOAT) E = ERANGE;
  99.            else Z = exp (X);
  100.            break;
  101.         case 2:                     /* square root */
  102.            if (X < 0) E = EDOM;
  103.            else Z = sqrt (X);
  104.            break;
  105.         case 3:                     /* sin */
  106.            Z = sin (X);
  107.            break;
  108.         case 4:                     /* cos */
  109.            Z = cos (X);
  110.            break;
  111.         case 5:                     /* tan */
  112.            Z = tan (X);
  113.            break;
  114.         case 6:                     /* arcsin */
  115.            Z = asin (X);
  116.            E = errno;
  117.            break;
  118.         case 7:                     /* arccos */
  119.            Z = acos (X);
  120.            E = errno;
  121.            break;
  122.         case 8:                     /* arctan */
  123.            Z = atan (X);
  124.            E = errno;
  125.            break;
  126. #endif /* OPSYS!=CTSS */
  127.         case 9:                     /* minus */
  128.            Z = -X;
  129.            E = 0;
  130.            break;
  131.      }
  132.      switch (E) {
  133. #if OPSYS!=CTSS
  134.         case EDOM:
  135.            FunError ("domain error",InOut);
  136.            break;
  137.         case ERANGE:
  138.            FunError ("range error",InOut);
  139.            break;
  140. #endif
  141.         default:
  142.            InOut->Tag = FLOAT;
  143.            InOut->Float = Z;
  144.            break;
  145.      }
  146.       }
  147.    }
  148.  
  149.  
  150. private F_Minus (InOut)
  151.    register ObjectPtr InOut;
  152.    {
  153.       if (InOut->Tag == INT && InOut->Int != FPMaxInt+1)
  154.      InOut->Int = - InOut->Int;
  155.       else Monadic (InOut,9);
  156.    }
  157.  
  158.  
  159. /*
  160.  * F_Sum
  161.  */
  162. private F_Sum (InOut)
  163.    ObjectPtr InOut;
  164.    {
  165.       Object S;
  166.       register ListPtr P;
  167.  
  168.       switch (InOut->Tag) {
  169.      default:
  170.         FunError (ArgNotSeq,InOut);
  171.         return;
  172.      case LIST:
  173.         S.Tag = INT;
  174.         S.Int = 0;
  175.         for (P=InOut->List; P!=NULL; P=P->Next) {
  176.            if (P->Val.Tag != INT && P->Val.Tag != FLOAT) {
  177.           FunError ("non-numeric sequence",InOut); 
  178.           return;
  179.            }
  180.            if (S.Tag == INT) {
  181.           if (P->Val.Tag == INT) {
  182.  
  183.              /* Both arguments are integers. See if we can avoid    */
  184.              /* floating arithmetic.                                */
  185.  
  186.              FPint Zi = S.Int + P->Val.Int;
  187.              if ((S.Int ^ P->Val.Int) < 0 || (S.Int^Zi)) 
  188.              /* arithmetic overflow occured - float result */;
  189.              else {
  190.             S.Int = Zi; 
  191.             continue;
  192.              }
  193.           }
  194.           S.Float = S.Int; 
  195.           S.Tag = FLOAT;
  196.            }
  197.            S.Float += P->Val.Tag==INT ? P->Val.Int : P->Val.Float;
  198.         }
  199.         break;
  200.       }
  201.       RepObject (InOut,&S);
  202.    }
  203.  
  204. /*
  205.  * Dyadic
  206.  *
  207.  * Evaluate a dyadic function
  208.  *
  209.  * Input
  210.  *      InOut = argument to apply function
  211.  *      Op = operation - see case statement in code for possibilities
  212.  *
  213.  * Output
  214.  *      InOut = result of applying function
  215.  *
  216.  * The author sold his anti-GOTO morals for speed.
  217.  */
  218. private Dyadic (InOut,Op)
  219.    register ObjectPtr InOut;
  220.    register int Op;
  221.    {
  222.       double X,Y,Z;
  223.       register FPint Xi,Yi,Zi;
  224.       register ListPtr P,Q;
  225.       static char *DivZero = "division by zero";
  226.  
  227.       if (InOut->Tag != LIST ||
  228.       NULL == (P=InOut->List) ||
  229.       NULL == (Q=P->Next) ||
  230.       Q->Next != NULL ||
  231.       NotNumPair (P->Val.Tag,Q->Val.Tag)) {
  232.  
  233.      FunError ("not a numeric pair",InOut);
  234.      return;
  235.       }
  236.  
  237.       if (IntPair (P->Val.Tag,Q->Val.Tag)) {
  238.  
  239.      /* Both arguments are integers. See if we can avoid floating point */
  240.      /* arithmetic.                                                     */
  241.  
  242.      Xi = P->Val.Int;
  243.      Yi = Q->Val.Int;
  244.  
  245.      switch (Op) {
  246.  
  247.         case 0:
  248.            /* assume two's complement arithmetic */
  249.            Zi = Xi+Yi;
  250.            if (((Xi ^ Yi) | ~(Xi ^ Zi)) < 0) goto RetInt;
  251.            break;
  252.            /* else arithmetic overflow occured */
  253.  
  254.         case 1:
  255.            /* assume two's complement arithmetic */
  256.            Zi = Xi - Yi;
  257.            if (((Xi ^ Yi) & (Xi ^ Zi)) >= 0) goto RetInt;
  258.            /* else arithmetic overflow occured */
  259.            break;
  260.  
  261.         case 2:
  262.            Zi = Xi * Yi;
  263.            if (Yi==0 || Zi/Yi == Xi) goto RetInt;
  264.            /* else arithmetic overflow occured */
  265.            break;
  266.  
  267.      /* case 3: division  result always FLOAT */
  268.  
  269. #if OPSYS!=CTSS
  270.         case 4:                     /* mod */
  271.            if (Xi >= 0 && Yi > 0) {
  272.           Zi = Xi % Yi;
  273.           goto RetInt;
  274.            }
  275.            break;
  276.  
  277.         case 5:                     /* div */
  278.            if (Xi >= 0 && Yi > 0) {
  279.           Zi = Xi / Yi;
  280.           goto RetInt;
  281.            }
  282.            break;
  283. #endif /* OPSYS!=CTSS */
  284.  
  285.         case 6:
  286.            Zi = Xi > Yi ? Yi : Xi;
  287.            goto RetInt;
  288.  
  289.         case 7:
  290.            Zi = Xi < Yi ? Yi : Xi;
  291.            goto RetInt;
  292.  
  293.      /* case 8: power result always FLOAT */
  294.      }
  295.       }
  296.  
  297.       X = P->Val.Tag==INT ? P->Val.Int : P->Val.Float;
  298.       Y = Q->Val.Tag==INT ? Q->Val.Int : Q->Val.Float;
  299.  
  300.       switch (Op) {
  301.      case 0: Z = X + Y; break;
  302.      case 1: Z = X - Y; break;
  303.      case 2: Z = X * Y; break;
  304.      case 3: 
  305.         if (Y==0.0) {
  306.            FunError (DivZero,InOut);
  307.            return;
  308.         }
  309.         Z = X / Y; 
  310.         break;
  311. #if OPSYS!=CTSS
  312.      case 4:
  313.         Z = Y==0.0 ? 0.0 : X - floor (X / Y) * Y;   /* mod */
  314.         break;
  315.      case 5:
  316.         if (Y==0.0) {                               /* div */
  317.            FunError (DivZero,InOut);
  318.            return;
  319.         }
  320.         Z = floor (X / Y);
  321.         break;
  322. #endif
  323.      case 6: Z = X > Y ? Y:X; break;
  324.      case 7: Z = X > Y ? X:Y; break;
  325. #if OPSYS!=CTSS
  326.      case 8: Z = pow (X,Y);   break;
  327. #endif
  328.       }
  329.       InOut->Tag = FLOAT;
  330.       InOut->Float = Z;
  331.  
  332.    Return:
  333.       DelLPtr (P);
  334.       return;
  335.  
  336.    RetInt: 
  337.       InOut->Tag = INT;
  338.       InOut->Int = Zi;
  339.       goto Return;
  340.    }
  341.  
  342.  
  343. /*
  344.  * F_Add1
  345.  */
  346. private F_AddN (InOut,N)
  347.    register ObjectPtr InOut;
  348.    int N;
  349.    {
  350.       register FPint K;
  351.  
  352.       switch (InOut->Tag) {
  353.      case INT:
  354.         K = InOut->Int + N;
  355.         if (N >= 0 ? InOut->Int <= K : InOut->Int >  K) {
  356.            InOut->Int = K;
  357.            return;
  358.         }
  359.         /* else integer overflow - convert and drop down */
  360.         InOut->Float = ((FPfloat) InOut->Int);
  361.         InOut->Tag = FLOAT;
  362.      case FLOAT:
  363.         InOut->Float = InOut->Float + N;
  364.         break;
  365.      default:
  366.         FunError ("not a number",InOut);
  367.         break;
  368.       }
  369.    }
  370.  
  371. void D_arith ()
  372.    {
  373.       GroupDef (OpArith,OpCount (OpArith), ArithNode);
  374.    }
  375.  
  376. /************************** end of F_arith.c **************************/
  377.  
  378.