Subversion Repositories Kolibri OS

Rev

Go to most recent revision | Blame | Last modification | View Log | Download | RSS feed

  1. /*                                                      sinhl.c
  2.  *
  3.  *      Hyperbolic sine, long double precision
  4.  *
  5.  *
  6.  *
  7.  * SYNOPSIS:
  8.  *
  9.  * long double x, y, sinhl();
  10.  *
  11.  * y = sinhl( x );
  12.  *
  13.  *
  14.  *
  15.  * DESCRIPTION:
  16.  *
  17.  * Returns hyperbolic sine of argument in the range MINLOGL to
  18.  * MAXLOGL.
  19.  *
  20.  * The range is partitioned into two segments.  If |x| <= 1, a
  21.  * rational function of the form x + x**3 P(x)/Q(x) is employed.
  22.  * Otherwise the calculation is sinh(x) = ( exp(x) - exp(-x) )/2.
  23.  *
  24.  *
  25.  *
  26.  * ACCURACY:
  27.  *
  28.  *                      Relative error:
  29.  * arithmetic   domain     # trials      peak         rms
  30.  *    IEEE       -2,2       10000       1.5e-19     3.9e-20
  31.  *    IEEE     +-10000      30000       1.1e-19     2.8e-20
  32.  *
  33.  */
  35. /*
  36. Cephes Math Library Release 2.7:  January, 1998
  37. Copyright 1984, 1991, 1998 by Stephen L. Moshier
  38. */
  39.  
  40. /*
  41. Modified for mingw
  42. 2002-07-22 Danny Smith <dannysmith@users.sourceforge.net>
  43. */
  44.  
  45. #ifdef __MINGW32__
  46. #include "cephes_mconf.h"
  47. #else
  48. #include "mconf.h"
  49. #endif
  50.  
  51. #ifndef _SET_ERRNO
  52. #define _SET_ERRNO(x)
  53. #endif
  54.  
  55. #ifdef UNK
  56. static long double P[] = {
  57.  1.7550769032975377032681E-6L,
  58.  4.1680702175874268714539E-4L,
  59.  3.0993532520425419002409E-2L,
  60.  9.9999999999999999998002E-1L,
  61. };
  62. static long double Q[] = {
  63.  1.7453965448620151484660E-8L,
  64. -5.9116673682651952419571E-6L,
  65.  1.0599252315677389339530E-3L,
  66. -1.1403880487744749056675E-1L,
  67.  6.0000000000000000000200E0L,
  68. };
  69. #endif
  70.  
  71. #ifdef IBMPC
  72. static const unsigned short P[] = {
  73. 0xec6a,0xd942,0xfbb3,0xeb8f,0x3feb, XPD
  74. 0x365e,0xb30a,0xe437,0xda86,0x3ff3, XPD
  75. 0x8890,0x01f6,0x2612,0xfde6,0x3ff9, XPD
  76. 0x0000,0x0000,0x0000,0x8000,0x3fff, XPD
  77. };
  78. static const unsigned short Q[] = {
  79. 0x4edd,0x4c21,0xad09,0x95ed,0x3fe5, XPD
  80. 0x4376,0x9b70,0xd605,0xc65c,0xbfed, XPD
  81. 0xc8ad,0x5d21,0x3069,0x8aed,0x3ff5, XPD
  82. 0x9c32,0x6374,0x2d4b,0xe98d,0xbffb, XPD
  83. 0x0000,0x0000,0x0000,0xc000,0x4001, XPD
  84. };
  85. #endif
  86.  
  87. #ifdef MIEEE
  88. static long P[] = {
  89. 0x3feb0000,0xeb8ffbb3,0xd942ec6a,
  90. 0x3ff30000,0xda86e437,0xb30a365e,
  91. 0x3ff90000,0xfde62612,0x01f68890,
  92. 0x3fff0000,0x80000000,0x00000000,
  93. };
  94. static long Q[] = {
  95. 0x3fe50000,0x95edad09,0x4c214edd,
  96. 0xbfed0000,0xc65cd605,0x9b704376,
  97. 0x3ff50000,0x8aed3069,0x5d21c8ad,
  98. 0xbffb0000,0xe98d2d4b,0x63749c32,
  99. 0x40010000,0xc0000000,0x00000000,
  100. };
  101. #endif
  102.  
  103. #ifndef __MINGW32__
  104. extern long double MAXNUML, MAXLOGL, MINLOGL, LOGE2L;
  105. #ifdef ANSIPROT
  106. extern long double fabsl ( long double );
  107. extern long double expl ( long double );
  108. extern long double polevll ( long double, void *, int );
  109. extern long double p1evll ( long double, void *, int );
  110. #else
  111. long double fabsl(), expl(), polevll(), p1evll();
  112. #endif
  113. #ifdef INFINITIES
  114. extern long double INFINITYL;
  115. #endif
  116. #ifdef NANS
  117. extern long double NANL;
  118. #endif
  119. #endif /* __MINGW32__ */
  120.  
  121. long double sinhl(x)
  122. long double x;
  123. {
  124. long double a;
  125.  
  126. #ifdef MINUSZERO
  127. if( x == 0.0 )
  128.         return(x);
  129. #endif
  130. #ifdef NANS
  131. if (isnanl(x))
  132.         {
  133.         _SET_ERRNO(EDOM);
  134.         }
  135. #endif
  136. a = fabsl(x);
  137. if( (x > (MAXLOGL + LOGE2L)) || (x > -(MINLOGL-LOGE2L) ) )
  138.         {
  139.         mtherr( "sinhl", DOMAIN );
  140.         _SET_ERRNO(ERANGE);
  141. #ifdef INFINITIES
  142.         if( x > 0.0L )
  143.                 return( INFINITYL );
  144.         else
  145.                 return( -INFINITYL );
  146. #else
  147.         if( x > 0.0L )
  148.                 return( MAXNUML );
  149.         else
  150.                 return( -MAXNUML );
  151. #endif
  152.         }
  153. if( a > 1.0L )
  154.         {
  155.         if( a >= (MAXLOGL - LOGE2L) )
  156.                 {
  157.                 a = expl(0.5L*a);
  158.                 a = (0.5L * a) * a;
  159.                 if( x < 0.0L )
  160.                         a = -a;
  161.                 return(a);
  162.                 }
  163.         a = expl(a);
  164.         a = 0.5L*a - (0.5L/a);
  165.         if( x < 0.0L )
  166.                 a = -a;
  167.         return(a);
  168.         }
  169.  
  170. a *= a;
  171. return( x + x * a * (polevll(a,P,3)/polevll(a,Q,4)) );
  172. }
  173.