atan.go 2.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990
  1. package math32
  2. func Atan(x float32) float32 {
  3. if x == 0 {
  4. return x
  5. }
  6. if x > 0 {
  7. return satan(x)
  8. }
  9. return -satan(-x)
  10. }
  11. // The original C code, the long comment, and the constants below were
  12. // from http://netlib.sandia.gov/cephes/cmath/atan.c, available from
  13. // http://www.netlib.org/cephes/cmath.tgz.
  14. // The go code is a version of the original C.
  15. //
  16. // atan.c
  17. // Inverse circular tangent (arctangent)
  18. //
  19. // SYNOPSIS:
  20. // double x, y, atan();
  21. // y = atan( x );
  22. //
  23. // DESCRIPTION:
  24. // Returns radian angle between -pi/2 and +pi/2 whose tangent is x.
  25. //
  26. // Range reduction is from three intervals into the interval from zero to 0.66.
  27. // The approximant uses a rational function of degree 4/5 of the form
  28. // x + x**3 P(x)/Q(x).
  29. //
  30. // ACCURACY:
  31. // Relative error:
  32. // arithmetic domain # trials peak rms
  33. // DEC -10, 10 50000 2.4e-17 8.3e-18
  34. // IEEE -10, 10 10^6 1.8e-16 5.0e-17
  35. //
  36. // Cephes Math Library Release 2.8: June, 2000
  37. // Copyright 1984, 1987, 1989, 1992, 2000 by Stephen L. Moshier
  38. //
  39. // The readme file at http://netlib.sandia.gov/cephes/ says:
  40. // Some software in this archive may be from the book _Methods and
  41. // Programs for Mathematical Functions_ (Prentice-Hall or Simon & Schuster
  42. // International, 1989) or from the Cephes Mathematical Library, a
  43. // commercial product. In either event, it is copyrighted by the author.
  44. // What you see here may be used freely but it comes with no support or
  45. // guarantee.
  46. //
  47. // The two known misprints in the book are repaired here in the
  48. // source listings for the gamma function and the incomplete beta
  49. // integral.
  50. //
  51. // Stephen L. Moshier
  52. // moshier@na-net.ornl.gov
  53. // xatan evaluates a series valid in the range [0, 0.66].
  54. func xatan(x float32) float32 {
  55. const (
  56. P0 = -8.750608600031904122785e-01
  57. P1 = -1.615753718733365076637e+01
  58. P2 = -7.500855792314704667340e+01
  59. P3 = -1.228866684490136173410e+02
  60. P4 = -6.485021904942025371773e+01
  61. Q0 = +2.485846490142306297962e+01
  62. Q1 = +1.650270098316988542046e+02
  63. Q2 = +4.328810604912902668951e+02
  64. Q3 = +4.853903996359136964868e+02
  65. Q4 = +1.945506571482613964425e+02
  66. )
  67. z := x * x
  68. z = z * ((((P0*z+P1)*z+P2)*z+P3)*z + P4) / (((((z+Q0)*z+Q1)*z+Q2)*z+Q3)*z + Q4)
  69. z = x*z + x
  70. return z
  71. }
  72. // satan reduces its argument (known to be positive)
  73. // to the range [0, 0.66] and calls xatan.
  74. func satan(x float32) float32 {
  75. const (
  76. Morebits float32 = 6.123233995736765886130e-17 // pi/2 = PIO2 + Morebits
  77. Tan3pio8 float32 = 2.41421356237309504880 // tan(3*pi/8)
  78. )
  79. if x <= 0.66 {
  80. return xatan(x)
  81. }
  82. if x > Tan3pio8 {
  83. return Pi/2 - xatan(1/x) + Morebits
  84. }
  85. return Pi/4 + xatan((x-1)/(x+1)) + 0.5*Morebits
  86. }