home *** CD-ROM | disk | FTP | other *** search
/ OS/2 Shareware BBS: Science / Science.zip / gmt_os2.zip / src / math / readme < prev    next >
Encoding:
Text File  |  1996-09-03  |  8.5 KB  |  262 lines
  1.  
  2.       ********************************
  3.        * Announcing FDLIBM Version 5  *
  4.       ********************************
  5. ============================================================
  6.             FDLIBM
  7. ============================================================
  8.   (developed at SunSoft, a Sun Microsystems, Inc. business.)
  9.  
  10. What's new in FDLIBM 5.2?
  11. BUGS FIXED
  12.     1. Little endian bug in frexp (affect only little endian machine):
  13.        in file s_frexp.c, last line of program frexp before exit
  14.                 *(int*)&x = hx;
  15.        should read
  16.                 *(n0+(int*)&x) = hx;
  17.  
  18.     2. jn(-1,x) is three times larger that the actual answer:
  19.        in file e_jn.c, the line 
  20.          sign = 1 - ((n&1)<<2);
  21.        should read
  22.          sign = 1 - ((n&1)<<1);
  23.  
  24.     3. Compiler failure on non-standard code
  25.        J.T. Conklin found that gcc optimizing out the manipulation of doubles
  26.        via pointer bashing of the form
  27.         double x = 0;
  28.         *(((int*)&x)+n0)=0x7fff0000;
  29.         foo(x);
  30.        C experts confirmed that the behavior of *(((int*)&x)+n0)=0x7fff0000
  31.        is undefined. By replacing n0 with a constant 0 or 1, the GCC "knows" 
  32.        that the assignment is modifying the double, and "does the right thing."
  33.        Thus, in FDLIBM 5.2, we replace n0 with a constant and use a macro 
  34.        __HI() and __LO() with #ifdef __LITTLE_ENDIAN to avoid the above problem.
  35.  
  36.     4. Performance issue on rem_pio2
  37.        An attempt to speed up the argument reduction in the trig function is the
  38.        consider pi/4 < x < 3pi/4 a special case. This was done in the file
  39.        e_rem_pio2.c
  40.     
  41.     
  42. FDLIBM (Freely Distributable LIBM) is a C math library 
  43. for machines that support IEEE 754 floating-point arithmetic. 
  44. In this release, only double precision is supported.
  45.  
  46. FDLIBM is intended to provide a reasonably portable (see 
  47. assumptions below), reference quality (below one ulp for
  48. major functions like sin,cos,exp,log) math library 
  49. (libm.a).  For a copy of FDLIBM, please send a message "send index from fdlibm"
  50. to netlib@research.bell-labs.com.
  51.  
  52. --------------
  53. 1. ASSUMPTIONS
  54. --------------
  55. FDLIBM (double precision version) assumes:
  56.  a.  IEEE 754 style (if not precise compliance) arithmetic;
  57.  b.  32 bit 2's complement integer arithmetic;
  58.  c.  Each double precision floating-point number must be in IEEE 754 
  59.      double format, and that each number can be retrieved as two 32-bit 
  60.      integers through the using of pointer bashing as in the example 
  61.      below:
  62.  
  63.      Example: let y = 2.0
  64.     double fp number y:     2.0
  65.     IEEE double format:    0x4000000000000000
  66.  
  67.     Referencing y as two integers:
  68.     *(int*)&y,*(1+(int*)&y) =    {0x40000000,0x0} (on sparc)
  69.                     {0x0,0x40000000} (on 386)
  70.  
  71.     Note: Four macros are defined in fdlibm.h to handle this kind of
  72.     retrieving:
  73.  
  74.     __HI(x)        the high part of a double x 
  75.             (sign,exponent,the first 21 significant bits)
  76.     __LO(x)        the least 32 significant bits of x
  77.     __HIp(x)    same as __HI except that the argument is a pointer
  78.             to a double
  79.     __LOp(x)    same as __LO except that the argument is a pointer
  80.             to a double
  81.     
  82.     To ensure obtaining correct ordering, one must define  __LITTLE_ENDIAN
  83.     during compilation for little endian machine (like 386,486). The 
  84.     default is big endian.
  85.  
  86.     If the behavior of pointer bashing is undefined, one may hack on the 
  87.     macro in fdlibm.h.
  88.     
  89.   d. IEEE exceptions may trigger "signals" as is common in Unix
  90.      implementations. 
  91.  
  92. -------------------
  93. 2. EXCEPTION CASES
  94. -------------------
  95. All exception cases in the FDLIBM functions will be mapped
  96. to one of the following four exceptions:
  97.  
  98.    +-huge*huge, +-tiny*tiny,    +-1.0/0.0,    +-0.0/0.0
  99.     (overflow)    (underflow)  (divided-by-zero)     (invalid)
  100.  
  101. For example, log(0) is a singularity and is thus mapped to 
  102.     -1.0/0.0 = -infinity.
  103. That is, FDLIBM's log will compute -one/zero and return the
  104. computed value.  On an IEEE machine, this will trigger the 
  105. divided-by-zero exception and a negative infinity is returned by 
  106. default.
  107.  
  108. Similarly, exp(-huge) will be mapped to tiny*tiny to generate
  109. an underflow signal. 
  110.  
  111.  
  112. --------------------------------
  113. 3. STANDARD CONFORMANCE WRAPPER 
  114. --------------------------------
  115. The default FDLIBM functions (compiled with -D_IEEE_LIBM flag)  
  116. are in "IEEE spirit" (i.e., return the most reasonable result in 
  117. floating-point arithmetic). If one wants FDLIBM to comply with
  118. standards like SVID, X/OPEN, or POSIX/ANSI, then one can 
  119. create a multi-standard compliant FDLIBM. In this case, each
  120. function in FDLIBM is actually a standard compliant wrapper
  121. function.  
  122.  
  123. File organization:
  124.     1. For FDLIBM's kernel (internal) function,
  125.         File name    Entry point
  126.         ---------------------------
  127.         k_sin.c        __kernel_sin
  128.         k_tan.c        __kernel_tan
  129.         ---------------------------
  130.     2. For functions that have no standards conflict 
  131.         File name    Entry point
  132.         ---------------------------
  133.         s_sin.c        sin
  134.         s_erf.c        erf
  135.         ---------------------------
  136.     3. Ieee754 core functions
  137.         File name    Entry point
  138.         ---------------------------
  139.         e_exp.c        __ieee754_exp
  140.         e_sinh.c    __ieee754_sinh
  141.         ---------------------------
  142.     4. Wrapper functions
  143.         File name    Entry point
  144.         ---------------------------
  145.         w_exp.c        exp
  146.         w_sinh.c    sinh
  147.         ---------------------------
  148.  
  149. Wrapper functions will twist the result of the ieee754 
  150. function to comply to the standard specified by the value 
  151. of _LIB_VERSION 
  152.     if _LIB_VERSION = _IEEE_, return the ieee754 result;
  153.     if _LIB_VERSION = _SVID_, return SVID result;
  154.     if _LIB_VERSION = _XOPEN_, return XOPEN result;
  155.     if _LIB_VERSION = _POSIX_, return POSIX/ANSI result.
  156. (These are macros, see fdlibm.h for their definition.)
  157.  
  158.  
  159. --------------------------------
  160. 4. HOW TO CREATE FDLIBM's libm.a
  161. --------------------------------
  162. There are two types of libm.a. One is IEEE only, and the other is
  163. multi-standard compliant (supports IEEE,XOPEN,POSIX/ANSI,SVID).
  164.  
  165. To create the IEEE only libm.a, use 
  166.         make "CFLAGS = -D_IEEE_LIBM"     
  167. This will create an IEEE libm.a, which is smaller in size, and 
  168. somewhat faster.
  169.  
  170. To create a multi-standard compliant libm, use
  171.     make "CFLAGS = -D_IEEE_MODE"   --- multi-standard fdlibm: default
  172.                      to IEEE
  173.     make "CFLAGS = -D_XOPEN_MODE"  --- multi-standard fdlibm: default
  174.                      to X/OPEN
  175.     make "CFLAGS = -D_POSIX_MODE"  --- multi-standard fdlibm: default
  176.                      to POSIX/ANSI
  177.     make "CFLAGS = -D_SVID3_MODE"  --- multi-standard fdlibm: default
  178.                      to SVID
  179.  
  180.  
  181. Here is how one makes a SVID compliant libm.
  182.     Make the library by
  183.         make "CFLAGS = -D_SVID3_MODE".
  184.     The libm.a of FDLIBM will be multi-standard compliant and 
  185.     _LIB_VERSION is initialized to the value _SVID_ . 
  186.  
  187.     example1:
  188.     ---------
  189.         main()
  190.         {
  191.         double y0();
  192.         printf("y0(1e300) = %1.20e\n",y0(1e300));
  193.         exit(0);
  194.         }
  195.  
  196.     % cc example1.c libm.a
  197.     % a.out
  198.     y0: TLOSS error
  199.     y0(1e300) = 0.00000000000000000000e+00
  200.  
  201.  
  202. It is possible to change the default standard in multi-standard 
  203. fdlibm. Here is an example of how to do it:
  204.     example2:
  205.     ---------
  206.     #include "fdlibm.h"    /* must include FDLIBM's fdlibm.h */
  207.     main()
  208.     {
  209.         double y0();
  210.         _LIB_VERSION =  _IEEE_;
  211.         printf("IEEE: y0(1e300) = %1.20e\n",y0(1e300));
  212.         _LIB_VERSION = _XOPEN_;
  213.         printf("XOPEN y0(1e300) = %1.20e\n",y0(1e300));
  214.         _LIB_VERSION = _POSIX_;
  215.         printf("POSIX y0(1e300) = %1.20e\n",y0(1e300));
  216.         _LIB_VERSION = _SVID_;
  217.         printf("SVID  y0(1e300) = %1.20e\n",y0(1e300));
  218.         exit(0);
  219.     }
  220.  
  221.     % cc example2.c libm.a
  222.     % a.out
  223.       IEEE: y0(1e300) = -1.36813604503424810557e-151
  224.       XOPEN y0(1e300) = 0.00000000000000000000e+00
  225.       POSIX y0(1e300) = 0.00000000000000000000e+00
  226.       y0: TLOSS error
  227.       SVID  y0(1e300) = 0.00000000000000000000e+00
  228.  
  229. Note:    Here _LIB_VERSION is a global variable. If global variables 
  230.     are forbidden, then one should modify fdlibm.h to change
  231.     _LIB_VERSION to be a global constant. In this case, one
  232.     may not change the value of _LIB_VERSION as in example2.
  233.  
  234. ---------------------------
  235. 5. NOTES ON PORTING FDLIBM
  236. ---------------------------
  237.     Care must be taken when installing FDLIBM over existing
  238.     libm.a.
  239.     All co-existing function prototypes must agree, otherwise
  240.     users will encounter mysterious failures.
  241.  
  242.     So far, the only known likely conflict is the declaration 
  243.     of the IEEE recommended function scalb:
  244.  
  245.         double scalb(double,double)    (1)    SVID3 defined
  246.         double scalb(double,int)    (2)    IBM,DEC,...
  247.  
  248.     FDLIBM follows Sun definition and use (1) as default. 
  249.     If one's existing libm.a uses (2), then one may raise
  250.     the flags _SCALB_INT during the compilation of FDLIBM
  251.     to get the correct function prototype.
  252.     (E.g., make "CFLAGS = -D_IEEE_LIBM -D_SCALB_INT".)
  253.     NOTE that if -D_SCALB_INT is raised, it won't be SVID3
  254.     conformant.
  255.  
  256. --------------
  257. 6. PROBLEMS ?
  258. --------------
  259. Please send comments and bug report to: 
  260.         fdlibm-comments@sunpro.eng.sun.com
  261.  
  262.