tan.go 3.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126
  1. package math32
  2. /*
  3. Floating-point tangent.
  4. */
  5. // The original C code, the long comment, and the constants
  6. // below were from http://netlib.sandia.gov/cephes/cmath/sin.c,
  7. // available from http://www.netlib.org/cephes/cmath.tgz.
  8. // The go code is a simplified version of the original C.
  9. //
  10. // tan.c
  11. //
  12. // Circular tangent
  13. //
  14. // SYNOPSIS:
  15. //
  16. // double x, y, tan();
  17. // y = tan( x );
  18. //
  19. // DESCRIPTION:
  20. //
  21. // Returns the circular tangent of the radian argument x.
  22. //
  23. // Range reduction is modulo pi/4. A rational function
  24. // x + x**3 P(x**2)/Q(x**2)
  25. // is employed in the basic interval [0, pi/4].
  26. //
  27. // ACCURACY:
  28. // Relative error:
  29. // arithmetic domain # trials peak rms
  30. // DEC +-1.07e9 44000 4.1e-17 1.0e-17
  31. // IEEE +-1.07e9 30000 2.9e-16 8.1e-17
  32. //
  33. // Partial loss of accuracy begins to occur at x = 2**30 = 1.074e9. The loss
  34. // is not gradual, but jumps suddenly to about 1 part in 10e7. Results may
  35. // be meaningless for x > 2**49 = 5.6e14.
  36. // [Accuracy loss statement from sin.go comments.]
  37. //
  38. // Cephes Math Library Release 2.8: June, 2000
  39. // Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
  40. //
  41. // The readme file at http://netlib.sandia.gov/cephes/ says:
  42. // Some software in this archive may be from the book _Methods and
  43. // Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
  44. // International, 1989) or from the Cephes Mathematical Library, a
  45. // commercial product. In either event, it is copyrighted by the author.
  46. // What you see here may be used freely but it comes with no support or
  47. // guarantee.
  48. //
  49. // The two known misprints in the book are repaired here in the
  50. // source listings for the gamma function and the incomplete beta
  51. // integral.
  52. //
  53. // Stephen L. Moshier
  54. // moshier@na-net.ornl.gov
  55. // tan coefficients
  56. var _tanP = [...]float32{
  57. -1.30936939181383777646e4, // 0xc0c992d8d24f3f38
  58. 1.15351664838587416140e6, // 0x413199eca5fc9ddd
  59. -1.79565251976484877988e7, // 0xc1711fead3299176
  60. }
  61. var _tanQ = [...]float32{
  62. 1.00000000000000000000e0,
  63. 1.36812963470692954678e4, //0x40cab8a5eeb36572
  64. -1.32089234440210967447e6, //0xc13427bc582abc96
  65. 2.50083801823357915839e7, //0x4177d98fc2ead8ef
  66. -5.38695755929454629881e7, //0xc189afe03cbe5a31
  67. }
  68. func Tan(x float32) float32 {
  69. return tan(x)
  70. }
  71. func tan(x float32) float32 {
  72. const (
  73. PI4A = 7.85398125648498535156e-1 // 0x3fe921fb40000000, Pi/4 split into three parts
  74. PI4B = 3.77489470793079817668e-8 // 0x3e64442d00000000,
  75. PI4C = 2.69515142907905952645e-15 // 0x3ce8469898cc5170,
  76. )
  77. // special cases
  78. switch {
  79. case x == 0 || IsNaN(x):
  80. return x // return ±0 || NaN()
  81. case IsInf(x, 0):
  82. return NaN()
  83. }
  84. // make argument positive but save the sign
  85. sign := false
  86. if x < 0 {
  87. x = -x
  88. sign = true
  89. }
  90. var j uint64
  91. var y, z float32
  92. if x >= reduceThreshold {
  93. j, z = trigReduce(x)
  94. } else {
  95. j = uint64(x * (4 / Pi)) // integer part of x/(Pi/4), as integer for tests on the phase angle
  96. y = float32(j) // integer part of x/(Pi/4), as float
  97. /* map zeros and singularities to origin */
  98. if j&1 == 1 {
  99. j++
  100. y++
  101. }
  102. z = ((x - y*PI4A) - y*PI4B) - y*PI4C
  103. }
  104. zz := z * z
  105. if zz > 1e-14 {
  106. y = z + z*(zz*(((_tanP[0]*zz)+_tanP[1])*zz+_tanP[2])/((((zz+_tanQ[1])*zz+_tanQ[2])*zz+_tanQ[3])*zz+_tanQ[4]))
  107. } else {
  108. y = z
  109. }
  110. if j&2 == 2 {
  111. y = -1 / y
  112. }
  113. if sign {
  114. y = -y
  115. }
  116. return y
  117. }