home *** CD-ROM | disk | FTP | other *** search
/ Atari FTP / ATARI_FTP_0693.zip / ATARI_FTP_0693 / Mint / mntlib32.zoo / atof.c < prev    next >
C/C++ Source or Header  |  1993-02-28  |  17KB  |  793 lines

  1. /*
  2.  *     double strtod (str, endptr);
  3.  *     const char *str;
  4.  *     char **endptr;
  5.  *        if !NULL, on return, points to char in str where conv. stopped
  6.  *
  7.  *     double atof (str)
  8.  *     const char *str;
  9.  *
  10.  * recognizes:
  11.  [spaces] [sign] digits [ [.] [ [digits] [ [e|E|d|D] [space|sign] [int][F|L]]]]
  12.  *
  13.  * returns:
  14.  *    the number
  15.  *        on overflow: HUGE_VAL and errno = ERANGE
  16.  *        on underflow: -HUGE_VAL and errno = ERANGE
  17.  *
  18.  * bugs:
  19.  *   naive over/underflow detection
  20.  *
  21.  *    ++jrb bammi@dsrgsun.ces.cwru.edu
  22.  *
  23.  * 07/29/89:
  24.  *    ok, before you beat me over the head, here is a hopefully
  25.  *    much improved one, with proper regard for rounding/precision etc
  26.  *    (had to read up on everything i did'nt want to know in K. V2)
  27.  *    the old naive coding is at the end bracketed by ifdef __OLD__.
  28.  *    modeled after peter housels posting on comp.os.minix.
  29.  *    thanks peter!
  30.  *
  31.  *    mjr: 68881 version added
  32.  *
  33.  *    schwab: 68881 version replaced with new implementation that
  34.  *        uses fmovep for maximum precision, no bits lost anymore!
  35.  */
  36.  
  37. #if !defined (__M68881__) && !defined (sfp004)
  38.  
  39. #if !(defined(unix) || defined(minix))
  40. #include <stddef.h>
  41. #include <stdlib.h>
  42. #include <float.h>
  43. #endif
  44. #include <errno.h>
  45. #include <assert.h>
  46. #include <math.h>
  47.  
  48. #ifdef minix
  49. #include "lib.h"
  50. #endif
  51.  
  52. #ifndef _COMPILER_H
  53. #include <compiler.h>
  54. #endif
  55.  
  56. extern int errno;
  57. #if (defined(unix) || defined(minix))
  58. #define NULL     ((void *)0)
  59. #endif
  60.  
  61. #define Ise(c)        ((c == 'e') || (c == 'E') || (c == 'd') || (c == 'D'))
  62. #define Isdigit(c)    ((c <= '9') && (c >= '0'))
  63. #define Isspace(c)    ((c == ' ') || (c == '\t'))
  64. #define Issign(c)    ((c == '-') || (c == '+'))
  65. #define IsValidTrail(c) ((c == 'F') || (c == 'L'))
  66. #define Val(c)        ((c - '0'))
  67.  
  68. #define MAXDOUBLE    DBL_MAX
  69. #define MINDOUBLE    DBL_MIN
  70.  
  71. #define MAXF  1.797693134862316
  72. #define MINF  2.225073858507201
  73. #define MAXE  308
  74. #define MINE  (-308)
  75.  
  76. /* another alias for ieee double */
  77. struct ldouble {
  78.     unsigned long hi, lo;
  79. };
  80.  
  81. static int __ten_mul __PROTO((double *acc, int digit));
  82. static double __adjust __PROTO((double *acc, int dexp, int sign));
  83.  
  84. #ifdef __OLD__
  85. static double __ten_pow __PROTO((double r, int e));
  86. #endif
  87.  
  88. /*
  89.  * mul 64 bit accumulator by 10 and add digit
  90.  * algorithm:
  91.  *    10x = 2( 4x + x ) == ( x<<2 + x) << 1
  92.  *    result = 10x + digit
  93.  */
  94. static int __ten_mul(acc, digit)
  95. double *acc;
  96. int digit;
  97. {
  98.     register unsigned long d0, d1, d2, d3;
  99.     register          short i;
  100.     
  101.     d2 = d0 = ((struct ldouble *)acc)->hi;
  102.     d3 = d1 = ((struct ldouble *)acc)->lo;
  103.  
  104.     /* check possibility of overflow */
  105. /*    if( (d0 & 0x0FFF0000L) >= 0x0ccc0000L ) */
  106. /*    if( (d0 & 0x70000000L) != 0 ) */
  107.     if( (d0 & 0xF0000000L) != 0 )
  108.     /* report overflow somehow */
  109.     return 1;
  110.     
  111.     /* 10acc == 2(4acc + acc) */
  112.     for(i = 0; i < 2; i++)
  113.     {  /* 4acc == ((acc) << 2) */
  114.     asm volatile("    lsll    #1,%1;
  115.              roxll    #1,%0"    /* shift L 64 bit acc 1bit */
  116.         : "=d" (d0), "=d" (d1)
  117.         : "0"  (d0), "1"  (d1) );
  118.     }
  119.  
  120.     /* 4acc + acc */
  121.     asm volatile(" addl    %2,%0" : "=d" (d1) : "0" (d1), "d" (d3));
  122.     asm volatile(" addxl   %2,%0" : "=d" (d0) : "0" (d0), "d" (d2));
  123.  
  124.     /* (4acc + acc) << 1 */
  125.     asm volatile("  lsll    #1,%1;
  126.              roxll   #1,%0"    /* shift L 64 bit acc 1bit */
  127.         : "=d" (d0), "=d" (d1)
  128.         : "0"  (d0), "1"  (d1) );
  129.  
  130.     /* add in digit */
  131.     d2 = 0;
  132.     d3 = digit;
  133.     asm volatile(" addl    %2,%0" : "=d" (d1) : "0" (d1), "d" (d3));
  134.     asm volatile(" addxl   %2,%0" : "=d" (d0) : "0" (d0), "d" (d2));
  135.  
  136.  
  137.     /* stuff result back into acc */
  138.     ((struct ldouble *)acc)->hi = d0;
  139.     ((struct ldouble *)acc)->lo = d1;
  140.  
  141.     return 0;    /* no overflow */
  142. }
  143.  
  144. #include "flonum.h"
  145.  
  146. static double __adjust(acc, dexp, sign)
  147. double *acc;    /* the 64 bit accumulator */
  148. int     dexp;    /* decimal exponent       */
  149. int    sign;    /* sign flag          */
  150. {
  151.     register unsigned long  d0, d1, d2, d3;
  152.     register          short i;
  153.     register           short bexp = 0; /* binary exponent */
  154.     unsigned short tmp[4];
  155.     double r;
  156.  
  157. #ifdef __STDC__
  158.     double __normdf( double d, int exp, int sign, int rbits );
  159.     double ldexp(double d, int exp);
  160. #else
  161.     extern double __normdf();
  162.     extern double ldexp();
  163. #endif    
  164.     d0 = ((struct ldouble *)acc)->hi;
  165.     d1 = ((struct ldouble *)acc)->lo;
  166.     while(dexp != 0)
  167.     {    /* something to do */
  168.     if(dexp > 0)
  169.     { /* need to scale up by mul */
  170.         while(d0 > 429496729 ) /* 2**31 / 5 */
  171.         {    /* possibility of overflow, div/2 */
  172.         asm volatile(" lsrl    #1,%1;
  173.                     roxrl    #1,%0"
  174.             : "=d" (d1), "=d" (d0)
  175.             : "0"  (d1), "1"  (d0));
  176.         bexp++;
  177.         }
  178.         /* acc * 10 == 2(4acc + acc) */
  179.         d2 = d0;
  180.         d3 = d1;
  181.         for(i = 0; i < 2; i++)
  182.         {  /* 4acc == ((acc) << 2) */
  183.         asm volatile("    lsll    #1,%1;
  184.                  roxll    #1,%0"    /* shift L 64 bit acc 1bit */
  185.                  : "=d" (d0), "=d" (d1)
  186.                  : "0"  (d0), "1"  (d1) );
  187.         }
  188.  
  189.         /* 4acc + acc */
  190.         asm volatile(" addl    %2,%0" : "=d" (d1) : "0" (d1), "d" (d3));
  191.         asm volatile(" addxl   %2,%0" : "=d" (d0) : "0" (d0), "d" (d2));
  192.  
  193.         /* (4acc + acc) << 1 */
  194.         bexp++;    /* increment exponent to effectively acc * 10 */
  195.         dexp--;
  196.     }
  197.     else /* (dexp < 0) */
  198.     { /* scale down by 10 */
  199.         while((d0 & 0xC0000000L) == 0)
  200.         {    /* preserve prec by ensuring upper bits are set before div */
  201.         asm volatile(" lsll  #1,%1;
  202.                     roxll #1,%0" /* shift L to move bits up */
  203.             : "=d" (d0), "=d" (d1)
  204.             : "0"  (d0), "1"  (d1) );
  205.         bexp--;    /* compensate for the shift */
  206.         }
  207.         /* acc/10 = acc/5/2 */
  208.         *((unsigned long *)&tmp[0]) = d0;
  209.         *((unsigned long *)&tmp[2]) = d1;
  210.         d2 = (unsigned long)tmp[0];
  211.         asm volatile(" divu #5,%0" : "=d" (d2) : "0" (d2));
  212.         tmp[0] = (unsigned short)d2;    /* the quotient only */
  213.         for(i = 1; i < 4; i++)
  214.         {
  215.         d2 = (d2 & 0xFFFF0000L) | (unsigned long)tmp[i]; /* REM|next */
  216.         asm volatile(" divu #5,%0" : "=d" (d2) : "0" (d2));
  217.         tmp[i] = (unsigned short)d2;
  218.         }
  219.         d0 = *((unsigned long *)&tmp[0]);
  220.         d1 = *((unsigned long *)&tmp[2]);
  221.         /* acc/2 */
  222.         bexp--;    /* div/2 taken care of by decrementing binary exp */
  223.         dexp++;
  224.     }
  225.     }
  226.     
  227.     /* stuff the result back into acc */
  228.     ((struct ldouble *)acc)->hi = d0;
  229.     ((struct ldouble *)acc)->lo = d1;
  230.  
  231.     /* normalize it */
  232.     r = __normdf( *acc, ((0x03ff - 1) + 53), (sign)? -1 : 0, 0 );
  233.     /* now shove in the binary exponent */
  234.     return ldexp(r, bexp);
  235. }
  236.  
  237. /* flags */
  238. #define SIGN    0x01
  239. #define ESIGN    0x02
  240. #define DECP    0x04
  241. #define CONVF    0x08
  242.  
  243. double strtod (s, endptr)
  244. const register char *s;
  245. char **endptr;
  246. {
  247.     double         accum = 0.0;
  248.     register short flags = 0;
  249.     register short texp  = 0;
  250.     register short e     = 0;
  251.     double       zero = 0.0;
  252.     const char        *tmp;
  253.  
  254.     assert ((s != NULL));
  255.  
  256.     if(endptr != NULL) *endptr = (char *)s;
  257.     while(Isspace(*s)) s++;
  258.     if(*s == '\0')
  259.     {    /* just leading spaces */
  260.     return zero;
  261.     }
  262.  
  263.     if(Issign(*s))
  264.     {
  265.     if(*s == '-') flags = SIGN;
  266.     if(*++s == '\0')
  267.     {   /* "+|-" : should be an error ? */
  268.         return zero;
  269.     }
  270.     }
  271.  
  272.     do {
  273.     if (Isdigit(*s))
  274.     {
  275.         flags |= CONVF;
  276.         if( __ten_mul(&accum, Val(*s)) ) texp++;
  277.         if(flags & DECP) texp--;
  278.     }
  279.     else if(*s == '.')
  280.     {
  281.         if (flags & DECP)  /* second decimal point */
  282.         break;
  283.         flags |= DECP;
  284.     }
  285.     else
  286.         break;
  287.     s++;
  288.     } while (1);
  289.  
  290.     if(Ise(*s))
  291.     {
  292.     if(*++s != '\0') /* skip e|E|d|D */
  293.     {  /* ! ([s]xxx[.[yyy]]e)  */
  294.          tmp = s;
  295.         while(Isspace(*s)) s++; /* Ansi allows spaces after e */
  296.         if(*s != '\0')
  297.         { /*  ! ([s]xxx[.[yyy]]e[space])  */
  298.  
  299.         if(Issign(*s))
  300.             if(*s++ == '-') flags |= ESIGN;
  301.  
  302.         if(*s != '\0')
  303.         { /*  ! ([s]xxx[.[yyy]]e[s])  -- error?? */
  304.  
  305.             for(; Isdigit(*s); s++)
  306.             e = (((e << 2) + e) << 1) + Val(*s);
  307.  
  308.             if(IsValidTrail(*s)) s++;
  309.         
  310.             /* dont care what comes after this */
  311.             if(flags & ESIGN)
  312.             texp -= e;
  313.             else
  314.             texp += e;
  315.         }
  316.         }
  317.         if(s == tmp) s++;    /* back up pointer for a trailing e|E|d|D */
  318.     }
  319.     }
  320.  
  321.     if((endptr != NULL) && (flags & CONVF))
  322.     *endptr = (char *) s;
  323.     if(accum == zero) return zero;
  324.  
  325.     return __adjust(&accum, (int)texp, (int)(flags & SIGN));
  326. }
  327.  
  328. double atof(s)
  329. const char *s;
  330. {
  331. #ifdef __OLD__
  332.     extern double strtod();
  333. #endif
  334.     return strtod(s, (char **)NULL);
  335. }
  336.  
  337.  
  338. /* old naive coding */
  339. #ifdef __OLD__
  340. static double __ten_pow(r, e)
  341. double r;
  342. register int e;
  343. {
  344.     if(e < 0)
  345.     for(; e < 0; e++) r /= 10.0;
  346.     else
  347.     for(; e > 0; --e) r *= 10.0;
  348.     return r;
  349. }
  350.  
  351. #define RET(X)     {goto ret;}
  352.  
  353. double strtod (s, endptr)
  354. const register char *s;
  355. char **endptr;
  356. {
  357.     double f = 0.1, r = 0.0, accum = 0.0;
  358.     register int  e = 0, esign = 0, sign = 0;
  359.     register int texp = 0, countexp;
  360.     
  361.     assert ((s != NULL));
  362.     
  363.     while(Isspace(*s)) s++;
  364.     if(*s == '\0') RET(r);    /* just spaces */
  365.     
  366.     if(Issign(*s))
  367.     {
  368.     if(*s == '-') sign = 1;
  369.     s++;
  370.     if(*s == '\0') RET(r); /* "+|-" : should be an error ? */
  371.     }
  372.     countexp = 0;
  373.     while(Isdigit(*s))
  374.     {
  375.     if(!countexp && (*s != '0')) countexp = 1;
  376.     accum = (accum * 10.0) + Val(*s);
  377.     /* should check for overflow here somehow */
  378.     s++; 
  379.     if(countexp) texp++;
  380.     }
  381.     r = (sign ? (-accum) : accum);
  382.     if(*s == '\0') RET(r);  /* [s]xxx */
  383.     
  384.     countexp = (texp == 0);
  385.     if(*s == '.')
  386.     {
  387.     s++;
  388.     if(*s == '\0') RET(r); /* [s]xxx. */
  389.     if(!Ise(*s))
  390.     {
  391.         while(Isdigit(*s))
  392.         {
  393.         if(countexp && (*s == '0')) --texp;
  394.         if(countexp && (*s != '0')) countexp = 0;
  395.         accum = accum + (Val(*s) * f);
  396.         f = f / 10.0;
  397.         /* overflow (w + f) ? */
  398.         s++;
  399.         }
  400.         r = (sign ? (-accum) : accum);
  401.         if(*s == '\0') RET(r); /* [s]xxx.yyy  */
  402.     }
  403.     }
  404.     if(!Ise(*s)) RET(r);    /* [s]xxx[.[yyy]]  */
  405.     
  406.     s++; /* skip e */
  407.     if(*s == '\0') RET(r); /* [s]xxx[.[yyy]]e  */
  408.     
  409.     while(Isspace(*s)) s++;
  410.     if(*s == '\0') RET(r); /* [s]xxx[.[yyy]]e[space]  */
  411.     
  412.     if(Issign(*s))
  413.     {
  414.     if(*s == '-') esign = 1;
  415.     s++;
  416.     if(*s == '\0') RET(r); /* [s]xxx[.[yyy]]e[s]  */
  417.     }
  418.     
  419.     while(Isdigit(*s))
  420.     {
  421.     e = (e * 10) + Val(*s);
  422.     s++;
  423.     }
  424.     /* dont care what comes after this */
  425.     if(esign) e = (-e);
  426.     texp += e;
  427.     
  428.     /* overflow ? */ /* FIXME */
  429.     if( texp > MAXE)
  430.     {
  431.     if( ! ((texp == (MAXE+1)) && (accum <= MAXF)))
  432.     {
  433.         errno = ERANGE;
  434.         r = ((sign) ? -HUGE_VAL : HUGE_VAL);
  435.         RET(r);
  436.     }
  437.     }
  438.     
  439.     /* underflow -- in reality occurs before this */ /* FIXME */
  440.     if(texp < MINE)
  441.     {
  442.     errno = ERANGE;
  443.     r = ((sign) ? -HUGE_VAL : HUGE_VAL);
  444.     RET(r);
  445.     }
  446.     r = __ten_pow(r, e);
  447.     /* fall through */
  448.     
  449.     /* all returns end up here, with return value in r */
  450.   ret:
  451.     if(endptr != NULL)
  452.     *endptr = s;
  453.     return r;
  454. }
  455. #endif /* __OLD__ */
  456.  
  457. #else /* __M68881__ || sfp004 */
  458.  
  459. #include <string.h>
  460. #include <ctype.h>
  461. #include <math.h>
  462. #include <float.h>
  463. #include <errno.h>
  464. #include "flonum.h"
  465.  
  466. #ifdef sfp004
  467. /* Note: not tested -- schwab */
  468.  
  469. /* static */ double Float __PROTO ((char *));
  470.  
  471. asm ("
  472. comm = -6;
  473. resp = -16;
  474. .text
  475. .even
  476. _Float:
  477.     lea    0xfffffa50:w,a0            | fpu address
  478.     movel    sp@(4),a1
  479.     movew    #0x4c00,a0@(comm)        | fmovep -> fp0
  480. 1:    cmpw    #0x8900,a0@(resp)
  481.     beq    1b
  482.     movel    a1@+,a0@            | load data
  483.     movel    a1@+,a0@
  484.     movel    a1@+,a0@
  485.  
  486.     movew    #0x7400,a0@(comm)        | fmoved fp0 -> d0/d1
  487. 1:    cmpw    #0x8900,a0@(resp)
  488.     beq    1b
  489.     movel    a0@,d0
  490.     movel    a0@,d1
  491.     rts
  492. ");
  493. #endif
  494.  
  495. /* Format of packed decimal (from left to right):
  496.  
  497.     1 Bit: sign of mantissa
  498.     1 Bit: sign of exponent
  499.     2 Bits zero
  500.    12 Bits: three digits exponent
  501.     4 Bits unused, fourth (higher order) digit of exponent
  502.     8 Bits zero
  503.    68 Bits: 17 digits of mantissa, decimal point after first digit
  504.   --------
  505.    96 Bits == 12 Bytes
  506.  
  507.    All digits in BCD format.  */
  508.  
  509. double
  510. strtod (const char *str, const char **endptr)
  511. {
  512.   char packed_buf[12];
  513.   char *p;
  514.   int exponent, i, digits_seen;
  515.   union double_di di;
  516.  
  517.   if (endptr)
  518.     *endptr = str;
  519.   while (isspace (*str))
  520.     str++;
  521.   p = packed_buf;
  522.   for (i = 0; i < sizeof (packed_buf); i++)
  523.     *p++ = 0;
  524.   if (*str == '-')
  525.     {
  526.       packed_buf[0] = 0x80;
  527.       str++;
  528.     }
  529.   else if (*str == '+')
  530.     str++;
  531.   else if (*str == '\0')
  532.     return 0.0;
  533.   i = 0;
  534.   exponent = -1;
  535.   digits_seen = 0;
  536.   p = packed_buf + 3;
  537.   /* Ignore leading 0's. */
  538.   while (*str == '0')
  539.     {
  540.       digits_seen++;
  541.       str++;
  542.     }
  543.   while (isdigit (*str))
  544.     {
  545.       digits_seen++;
  546.       if (i < 17)
  547.     {
  548.       if (i & 1)
  549.         *p = (*str - '0') << 4;
  550.       else
  551.         *p++ |= *str - '0';
  552.       i++;
  553.     }
  554.       exponent++;
  555.       str++;
  556.     }
  557.   if (*str == '.')
  558.     {
  559.       str++;
  560.       if (i == 0)
  561.     {
  562.       /* Skip leading 0's.  */
  563.       while (*str == '0')
  564.         {
  565.           digits_seen++;
  566.           exponent--;
  567.           str++;
  568.         }
  569.     }
  570.       while (isdigit (*str))
  571.     {
  572.       digits_seen++;
  573.       if (i < 17)
  574.         {
  575.           if (i++ & 1)
  576.         *p = (*str - '0') << 4;
  577.           else
  578.         *p++ |= *str - '0';
  579.         }
  580.       str++;
  581.     }
  582.     }
  583.   /* Check that there were any digits.  */
  584.   if (!digits_seen)
  585.     return 0.0;
  586.  
  587.   if (*str == 'e' || *str == 'E' || *str == 'd' || *str == 'D')
  588.     {
  589.       const char *eptr = str;
  590.       int exp_given, exp_neg;
  591.  
  592.       str++;
  593.       while (isspace (*str))
  594.     str++;
  595.       exp_neg = 0;
  596.       if (*str == '-')
  597.     {
  598.       exp_neg = 1;
  599.       str++;
  600.     }
  601.       else if (*str == '+')
  602.     str++;
  603.       if (!isdigit (*str))
  604.     {
  605.       str = eptr;
  606.       goto convert;
  607.     }
  608.       while (*str == '0')
  609.     str++;
  610.       exp_given = 0;
  611.       while (isdigit (*str) && exp_given < 900)
  612.     {
  613.       exp_given = exp_given * 10 + *str - '0';
  614.       str++;
  615.     }
  616.       while (isdigit (*str))
  617.     str++;
  618.       if (exp_given >= 900)
  619.     {
  620.       /* Must be overflow/underflow.  */
  621.       if (endptr)
  622.         *endptr = str;
  623.       if (exp_neg)
  624.         return 0.0;
  625.       goto overflow;
  626.     }
  627.       if (exp_neg)
  628.     exponent -= exp_given;
  629.       else
  630.     exponent += exp_given;
  631.     }
  632.  convert:
  633.   if (endptr)
  634.     *endptr = str;
  635.   if (exponent < 0)
  636.     {
  637.       packed_buf[0] |= 0x40;
  638.       exponent = -exponent;
  639.     }
  640.   packed_buf[1] = exponent % 10;
  641.   packed_buf[1] |= ((exponent /= 10) % 10) << 4;
  642.   packed_buf[0] |= exponent / 10 % 10;
  643. #ifdef sfp004
  644.   di.d = Float (packed_buf);
  645. #else
  646.   __asm ("fmovep %1,%0" : "=f" (di.d) : "m" (packed_buf[0]));
  647. #endif
  648.   /* Check for overflow.  */
  649.   if ((di.i[0] & 0x7fffffff) >= 0x7ff00000)
  650.     {
  651.     overflow:
  652.       errno = ERANGE;
  653.       if (packed_buf[0] & 0x80)
  654.     return -DBL_MAX;
  655.       else
  656.     return DBL_MAX;
  657.     }
  658.   return di.d;
  659. }
  660.  
  661. double
  662. atof (const char *str)
  663. {
  664.   return strtod (str, NULL);
  665. }
  666.  
  667. #endif /* __M68881__ || sfp004 */
  668.  
  669. #ifdef TEST
  670. #if 0
  671. #ifdef __MSHORT__
  672. #error "please run this test in 32 bit int mode"
  673. #endif
  674. #endif
  675.  
  676. #define NTEST 10000L
  677.  
  678. #ifdef __MSHORT__
  679. # define    RAND_MAX    (0x7FFF)    /* maximum value from rand() */
  680. #else
  681. # define    RAND_MAX    (0x7FFFFFFFL)    /* maximum value from rand() */
  682. #endif
  683.  
  684. main()
  685. {
  686.     
  687.     double expected, result, e, max_abs_err;
  688.     char buf[128];
  689.     register long i, errs;
  690.     register long s;
  691. #ifdef __STDC__
  692.     double atof(const char *);
  693.     int rand(void);
  694. #else
  695.     extern double atof();
  696.     extern int rand();
  697. #endif
  698.  
  699. #if 0
  700.     expected = atof("3.14159265358979e23");
  701.     expected = atof("3.141");
  702.     expected = atof(".31415"); 
  703.     printf("%f\n\n", expected);
  704.     expected = atof("3.1415"); 
  705.     printf("%f\n\n", expected);
  706.     expected = atof("31.415"); 
  707.     printf("%f\n\n", expected);
  708.     expected = atof("314.15"); 
  709.     printf("%f\n\n", expected);
  710.  
  711.     expected = atof(".31415"); 
  712.     printf("%f\n\n", expected);
  713.     expected = atof(".031415"); 
  714.     printf("%f\n\n", expected);
  715.     expected = atof(".0031415"); 
  716.     printf("%f\n\n", expected);
  717.     expected = atof(".00031415"); 
  718.     printf("%f\n\n", expected);
  719.     expected = atof(".000031415"); 
  720.     printf("%f\n\n", expected);
  721.  
  722.     expected = atof("-3.1415e-9"); 
  723.     printf("%20.15e\n\n", expected);
  724.  
  725.     expected = atof("+3.1415e+009"); 
  726.     printf("%20.15e\n\n", expected);
  727. #endif
  728.  
  729.     expected = atof("+3.123456789123456789"); 
  730.     printf("%30.25e\n\n", expected);
  731.  
  732.     expected = atof(".000003123456789123456789"); 
  733.     printf("%30.25e\n\n", expected);
  734.  
  735.     expected = atof("3.1234567891234567890000000000"); 
  736.     printf("%30.25e\n\n", expected);
  737.  
  738.     expected = atof("9.22337999999999999999999999999999999999999999"); 
  739.     printf("%47.45e\n\n", expected);
  740.  
  741.     expected = atof("1.0000000000000000000"); 
  742.     printf("%25.19e\n\n", expected);
  743.  
  744.     expected = atof("1.00000000000000000000"); 
  745.     printf("%26.20e\n\n", expected);
  746.  
  747.     expected = atof("1.000000000000000000000"); 
  748.     printf("%27.21e\n\n", expected);
  749.  
  750.     expected = atof("1.000000000000000000000000"); 
  751.     printf("%30.24e\n\n", expected);
  752.  
  753.  
  754. #if 0
  755.     expected = atof("1.7e+308");
  756.     if(errno != 0)
  757.     {
  758.     printf("%d\n", errno);
  759.     }
  760.     else    printf("1.7e308 OK %g\n", expected);
  761.     expected = atof("1.797693e308");    /* anything gt looses */
  762.     if(errno != 0)
  763.     {
  764.     printf("%d\n", errno);
  765.     }
  766.     else    printf("Max OK %g\n", expected);
  767.     expected = atof("2.225073858507201E-307");
  768.     if(errno != 0)
  769.     {
  770.     printf("%d\n", errno, expected);
  771.     }
  772.     else    printf("Min OK %g\n", expected);
  773. #endif
  774.     
  775.     max_abs_err = 0.0;
  776.     for(errs = 0, i = 0; i < NTEST; i++)
  777.     {
  778.     expected = (double)(s = rand()) / (double)rand();
  779.     if(s > (RAND_MAX >> 1)) expected = -expected;
  780.     sprintf(buf, "%.14e", expected);
  781.     result = atof(buf);
  782.     e = (expected == 0.0) ? result : (result - expected)/expected;
  783.     if(e < 0) e = (-e);
  784.     if(e > 1.0e-6) 
  785.     {
  786.         errs++; printf("%.14e %s %.14e (%.14e)\n", expected, buf, result, e);
  787.     }
  788.     if (e > max_abs_err) max_abs_err = e;
  789.     }
  790.     printf("%ld Error(s), Max abs err %.14e\n", errs, max_abs_err);
  791. }
  792. #endif /* TEST */
  793.