atn.c 1.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172
  1. /*
  2. * (c) copyright 1988 by the Vrije Universiteit, Amsterdam, The Netherlands.
  3. * See the copyright notice in the ACK home directory, in the file "Copyright".
  4. *
  5. * Author: Ceriel J.H. Jacobs
  6. */
  7. /* $Id$ */
  8. #define __NO_DEFS
  9. #include <math.h>
  10. #if __STDC__
  11. #include <pc_math.h>
  12. #endif
  13. double
  14. _atn(x)
  15. double x;
  16. {
  17. /* Algorithm and coefficients from:
  18. "Software manual for the elementary functions"
  19. by W.J. Cody and W. Waite, Prentice-Hall, 1980
  20. */
  21. static double p[] = {
  22. -0.13688768894191926929e+2,
  23. -0.20505855195861651981e+2,
  24. -0.84946240351320683534e+1,
  25. -0.83758299368150059274e+0
  26. };
  27. static double q[] = {
  28. 0.41066306682575781263e+2,
  29. 0.86157349597130242515e+2,
  30. 0.59578436142597344465e+2,
  31. 0.15024001160028576121e+2,
  32. 1.0
  33. };
  34. static double a[] = {
  35. 0.0,
  36. 0.52359877559829887307710723554658381, /* pi/6 */
  37. M_PI_2,
  38. 1.04719755119659774615421446109316763 /* pi/3 */
  39. };
  40. int neg = x < 0;
  41. int n;
  42. double g;
  43. if (neg) {
  44. x = -x;
  45. }
  46. if (x > 1.0) {
  47. x = 1.0/x;
  48. n = 2;
  49. }
  50. else n = 0;
  51. if (x > 0.26794919243112270647) { /* 2-sqtr(3) */
  52. n = n + 1;
  53. x = (((0.73205080756887729353*x-0.5)-0.5)+x)/
  54. (1.73205080756887729353+x);
  55. }
  56. /* ??? avoid underflow ??? */
  57. g = x * x;
  58. x += x * g * POLYNOM3(g, p) / POLYNOM4(g, q);
  59. if (n > 1) x = -x;
  60. x += a[n];
  61. return neg ? -x : x;
  62. }