support.c
上传用户:nvosite88
上传日期:2007-01-17
资源大小:4983k
文件大小:9k
源码类别:

VxWorks

开发平台:

C/C++

  1. /* support.c - math routines */
  2. /* Copyright 1991-2001 Wind River Systems, Inc. */
  3. /*
  4. modification history
  5. --------------------
  6. 01j,03oct01,to   use IEEE little endian for ARM
  7. 01i,20apr01,h_k  fixed P_INDEX for ARM big_endian.
  8. 01h,11jan99,dra  removed non-ANSI aliasing code.
  9. 01g,12feb97,jpd  added support for ARM with early Cygnus release.
  10. 01f,23sep92,kdl  removed define of sparcHardSqrt (now in sqrt.c).
  11. 01e,23sep92,kdl  removed sqrt(), placed in sqrt.c.
  12. 01d,20sep92,smb  documentation additions.
  13. 01c,08jul92,smb  merged with ANSI and documentation.
  14.                  set EDOM in sqrt().
  15. 01b,02jul92,jwt  changed sysHardSqrt to sparcHardSqrt; removed warnings.
  16. 01a,25jun91,jwt  created sqrt() function select - hard versus soft.
  17. */
  18. /*
  19. DESCRIPTION
  20. * Copyright (c) 1985 Regents of the University of California.
  21. * All rights reserved.
  22. *
  23. * Redistribution and use in source and binary forms are permitted
  24. * provided that the above copyright notice and this paragraph are
  25. * duplicated in all such forms and that any documentation,
  26. * advertising materials, and other materials related to such
  27. * distribution and use acknowledge that the software was developed
  28. * by the University of California, Berkeley.  The name of the
  29. * University may not be used to endorse or promote products derived
  30. * from this software without specific prior written permission.
  31. * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
  32. * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
  33. * WARRANTIES OF MERCHANTIBILITY AND FITNESS FOR A PARTICULAR PURPOSE.
  34. *
  35. * All recipients should regard themselves as participants in an ongoing
  36. * research project and hence should feel obligated to report their
  37. * experiences (good or bad) with these elementary function codes, using
  38. * the sendbug(8) program, to the authors.
  39. *
  40. * Some IEEE standard 754 recommended functions and remainder and sqrt for
  41. * supporting the C elementary functions.
  42. ******************************************************************************
  43. * WARNING:
  44. *      These codes are developed (in double) to support the C elementary
  45. * functions temporarily. They are not universal, and some of them are very
  46. * slow (in particular, drem and sqrt is extremely inefficient). Each
  47. * computer system should have its implementation of these functions using
  48. * its own assembler.
  49. ******************************************************************************
  50. *
  51. * IEEE 754 required operations:
  52. *     drem(x,p)
  53. *              returns  x REM y  =  x - [x/y]*y , where [x/y] is the integer
  54. *              nearest x/y; in half way case, choose the even one.
  55. *     sqrt(x)
  56. *              returns the square root of x correctly rounded according to
  57. * the rounding mod.
  58. *
  59. * IEEE 754 recommended functions:
  60. * (a) copysign(x,y)
  61. *              returns x with the sign of y.
  62. * (b) scalb(x,N)
  63. *              returns  x * (2**N), for integer values N.
  64. * (c) logb(x)
  65. *              returns the unbiased exponent of x, a signed integer in
  66. *              double precision, except that logb(0) is -INF, logb(INF)
  67. *              is +INF, and logb(NAN) is that NAN.
  68. * (d) finite(x)
  69. *              returns the value TRUE if -INF < x < +INF and returns
  70. *              FALSE otherwise.
  71. *
  72. *
  73. * CODED IN C BY K.C. NG, 11/25/84;
  74. * REVISED BY K.C. NG on 1/22/85, 2/13/85, 3/24/85.
  75. *
  76. * 25jun91,jwt  created sqrt() function select - hard versus soft
  77. SEE ALSO: American National Standard X3.159-1989
  78. NOMANUAL
  79. */
  80. #include "vxWorks.h"
  81. #include "math.h"
  82. #include "private/mathP.h"
  83. #include "errno.h"
  84. #if (_BYTE_ORDER == _LITTLE_ENDIAN)
  85. #define P_INDEX 3
  86. #else /* little_endian */
  87. #define P_INDEX 0
  88. #endif /* little_endian */
  89. typedef union
  90.     {
  91.     unsigned short s[4];
  92.     double d;
  93.     } tdouble;
  94. #if defined(vax)||defined(tahoe)      /* VAX D format */
  95. #include <errno.h>
  96.     static unsigned short msign=0x7fff , mexp =0x7f80 ;
  97.     static short  prep1=57, gap=7, bias=129           ;
  98.     static double novf=1.7E38, nunf=3.0E-39, zero=0.0 ;
  99. #else /* defined(vax)||defined(tahoe) */
  100.     static unsigned short msign=0x7fff, mexp =0x7ff0  ;
  101.     static short prep1=54, gap=4, bias=1023           ;
  102.     static double novf=1.7E308, nunf=3.0E-308,zero=0.0;
  103. #endif /* defined(vax)||defined(tahoe) */
  104. /*****************************************************************************
  105. * scalb -
  106. *
  107. * RETURN:
  108. * NOMANUAL
  109. */
  110. double scalb(x,N)
  111. double x; int N;
  112. {
  113.         int k;
  114.         double scalb();
  115. tdouble y;
  116. y.d = x;
  117.         if( y.d == zero )
  118.     return(y.d);
  119. #if defined(vax)||defined(tahoe)
  120.         if( (k= y.s[P_INDEX] & mexp ) != ~msign )
  121.     {
  122.             if (N < -260)
  123. return(nunf*nunf);
  124.     else if (N > 260) {
  125. extern double infnan(),copysign();
  126. return(copysign(infnan(ERANGE),y.d));
  127.     }
  128. #else /* defined(vax)||defined(tahoe) */
  129.         if( (k= y.s[P_INDEX] & mexp ) != mexp )
  130.     {
  131.             if( N<-2100)
  132. return(nunf*nunf);
  133.     else if(N>2100)
  134. return(novf+novf);
  135.             if( k == 0 )
  136. {
  137. y.d *= scalb(1.0,(int)prep1);
  138. N -= prep1;
  139. return(scalb(y.d,N));
  140. }
  141. #endif /* defined(vax)||defined(tahoe) */
  142.             if((k = (k>>gap)+ N) > 0 )
  143.                 if( k < (mexp>>gap) )
  144.     y.s[P_INDEX] = (y.s[P_INDEX]&~mexp) | (k<<gap);
  145.                 else
  146.     x=novf+novf;               /* overflow */
  147.             else
  148.                 if( k > -prep1 )
  149.                                         /* gradual underflow */
  150.                     {
  151.     y.s[P_INDEX]=(y.s[P_INDEX]&~mexp)|(short)(1<<gap);
  152.     y.d *= scalb(1.0,k-1);
  153.     }
  154.                 else
  155.                 return(nunf*nunf);
  156.             }
  157.         return(y.d);
  158. }
  159. /*****************************************************************************
  160. * copysign - determine the sign
  161. *
  162. * RETURN:
  163. * NOMANUAL
  164. */
  165. double copysign(x,y)
  166. double x,y;
  167. {
  168. tdouble a, b;
  169. a.d = x;
  170. b.d = y;
  171. #if defined(vax)||defined(tahoe)
  172.         if ( (a.s[P_INDEX] & mexp) == 0 ) return(a.d);
  173. #endif /* defined(vax)||defined(tahoe) */
  174.         a.s[P_INDEX] = ( a.s[P_INDEX] & msign ) | ( b.s[P_INDEX] & ~msign );
  175.         return(a.d);
  176. }
  177. /*****************************************************************************
  178. * logb -
  179. *
  180. * RETURN:
  181. * NOMANUAL
  182. */
  183. double logb(x)
  184. double x;
  185. {
  186. tdouble a;
  187. short k;
  188. a.d = x;
  189. #if defined(vax)||defined(tahoe)
  190.         return (int)(((a.s[P_INDEX]&mexp)>>gap)-bias);
  191. #else /* defined(vax)||defined(tahoe) */
  192.         if( (k= a.s[P_INDEX] & mexp ) != mexp )
  193.             if ( k != 0 )
  194.                 return ( (k>>gap) - bias );
  195.             else if( a.d != zero)
  196.                 return ( -1022.0 );
  197.             else
  198.                 return(-(1.0/zero));
  199.         else if(a.d != a.d)
  200.             return(a.d);
  201.         else
  202.             {a.s[P_INDEX] &= msign; return(a.d);}
  203. #endif /* defined(vax)||defined(tahoe) */
  204. }
  205. /*****************************************************************************
  206. * finite - the finite value for this machine 
  207. *
  208. * RETURN:
  209. * NOMANUAL
  210. */
  211. int finite(x)
  212. double x;
  213. {
  214. tdouble a;
  215. a.d = x;
  216. #if defined(vax)||defined(tahoe)
  217.         return(1);
  218. #else /* defined(vax)||defined(tahoe) */
  219.         return( ( a.s[P_INDEX] & mexp ) != mexp );
  220. #endif /* defined(vax)||defined(tahoe) */
  221. }
  222. /*****************************************************************************
  223. * drem - remainder
  224. *
  225. * RETURN:
  226. * NOMANUAL
  227. */
  228. double drem(x,p)
  229. double x,p;
  230. {
  231.         short sign;
  232.         double hp, drem(),scalb();
  233. tdouble tmp, xd, pd, dp;
  234.         unsigned short  k;
  235. xd.d = x;
  236. pd.d = p;
  237.         pd.s[P_INDEX] &= msign;
  238. #if defined(vax)||defined(tahoe)
  239.         if( ( xd.s[P_INDEX] & mexp ) == ~msign ) /* is x a reserved operand? */
  240. #else /* defined(vax)||defined(tahoe) */
  241.         if( ( xd.s[P_INDEX] & mexp ) == mexp )
  242. #endif /* defined(vax)||defined(tahoe) */
  243. return  (xd.d-pd.d)-(xd.d-pd.d); /* create nan if x is inf */
  244. if (pd.d == zero) {
  245. #if defined(vax)||defined(tahoe)
  246. extern double infnan();
  247. return(infnan(EDOM));
  248. #else /* defined(vax)||defined(tahoe) */
  249. return zero/zero;
  250. #endif /* defined(vax)||defined(tahoe) */
  251. }
  252. #if defined(vax)||defined(tahoe)
  253.         if( ( pd.s[P_INDEX] & mexp ) == ~msign ) /* is p a reserved operand? */
  254. #else /* defined(vax)||defined(tahoe) */
  255.         if( ( pd.s[P_INDEX] & mexp ) == mexp )
  256. #endif /* defined(vax)||defined(tahoe) */
  257. { if (pd.d != pd.d) return pd.d; else return xd.d;}
  258.         else  if ( ((pd.s[P_INDEX] & mexp)>>gap) <= 1 )
  259.                 /* subnormal p, or almost subnormal p */
  260.             { double b; b=scalb(1.0,(int)prep1);
  261.               pd.d *= b; xd.d = drem(xd.d,pd.d); xd.d *= b; return(drem(xd.d,pd.d)/b);}
  262.         else  if ( pd.d >= novf/2)
  263.             { pd.d /= 2 ; xd.d /= 2; return(drem(xd.d,pd.d)*2);}
  264.         else
  265.             {
  266.                 dp.d=pd.d+pd.d;
  267. hp=pd.d/2;
  268.                 sign= xd.s[P_INDEX] & ~msign ;
  269.                 xd.s[P_INDEX] &= msign       ;
  270.                 while ( xd.d > dp.d )
  271.                     {
  272.                         k=(xd.s[P_INDEX] & mexp) - (dp.s[P_INDEX] & mexp) ;
  273.                         tmp.d = dp.d ;
  274.                         tmp.s[P_INDEX] += k ;
  275. #if defined(vax)||defined(tahoe)
  276.                         if( xd.d < tmp.d ) tmp.s[P_INDEX] -= 128 ;
  277. #else /* defined(vax)||defined(tahoe) */
  278.                         if( xd.d < tmp.d ) tmp.s[P_INDEX] -= 16 ;
  279. #endif /* defined(vax)||defined(tahoe) */
  280.                         xd.d -= tmp.d ;
  281.                     }
  282.                 if ( xd.d > hp )
  283.                     { xd.d -= pd.d ;  if ( xd.d >= hp ) xd.d -= pd.d ; }
  284. #if defined(vax)||defined(tahoe)
  285. if (xd.d)
  286. #endif /* defined(vax)||defined(tahoe) */
  287. xd.s[P_INDEX] ^= sign;
  288.                 return( xd.d);
  289.             }
  290. }