1/* Copyright JS Foundation and other contributors, http://js.foundation
2 *
3 * Licensed under the Apache License, Version 2.0 (the "License");
4 * you may not use this file except in compliance with the License.
5 * You may obtain a copy of the License at
6 *
7 *     http://www.apache.org/licenses/LICENSE-2.0
8 *
9 * Unless required by applicable law or agreed to in writing, software
10 * distributed under the License is distributed on an "AS IS" BASIS
11 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12 * See the License for the specific language governing permissions and
13 * limitations under the License.
14 *
15 * This file is based on work under the following copyright and permission
16 * notice:
17 *
18 *     Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
19 *
20 *     Developed at SunSoft, a Sun Microsystems, Inc. business.
21 *     Permission to use, copy, modify, and distribute this
22 *     software is freely granted, provided that this notice
23 *     is preserved.
24 *
25 *     @(#)e_atan2.c 1.3 95/01/18
26 */
27
28#include "jerry-libm-internal.h"
29
30/* atan2(y,x)
31 *
32 * Method:
33 *      1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
34 *      2. Reduce x to positive by (if x and y are unexceptional):
35 *              ARG (x+iy) = arctan(y/x)           ... if x > 0,
36 *              ARG (x+iy) = pi - arctan[y/(-x)]   ... if x < 0,
37 *
38 * Special cases:
39 *      ATAN2((anything), NaN ) is NaN;
40 *      ATAN2(NAN , (anything) ) is NaN;
41 *      ATAN2(+-0, +(anything but NaN)) is +-0  ;
42 *      ATAN2(+-0, -(anything but NaN)) is +-pi ;
43 *      ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2;
44 *      ATAN2(+-(anything but INF and NaN), +INF) is +-0 ;
45 *      ATAN2(+-(anything but INF and NaN), -INF) is +-pi;
46 *      ATAN2(+-INF,+INF ) is +-pi/4 ;
47 *      ATAN2(+-INF,-INF ) is +-3pi/4;
48 *      ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2;
49 *
50 * Constants:
51 * The hexadecimal values are the intended ones for the following
52 * constants. The decimal values may be used, provided that the
53 * compiler will convert from decimal to binary accurately enough
54 * to produce the hexadecimal values shown.
55 */
56
57#define tiny   1.0e-300
58#define zero   0.0
59#define pi_o_4 7.8539816339744827900E-01 /* 0x3FE921FB, 0x54442D18 */
60#define pi_o_2 1.5707963267948965580E+00 /* 0x3FF921FB, 0x54442D18 */
61#define pi     3.1415926535897931160E+00 /* 0x400921FB, 0x54442D18 */
62#define pi_lo  1.2246467991473531772E-16 /* 0x3CA1A626, 0x33145C07 */
63
64double
65atan2 (double y, double x)
66{
67  double_accessor z;
68  int k, m, hx, hy, ix, iy;
69  unsigned lx, ly;
70
71  hx = __HI (x);
72  ix = hx & 0x7fffffff;
73  lx = __LO (x);
74  hy = __HI (y);
75  iy = hy & 0x7fffffff;
76  ly = __LO (y);
77  if (((ix | ((lx | -lx) >> 31)) > 0x7ff00000) || ((iy | ((ly | -ly) >> 31)) > 0x7ff00000)) /* x or y is NaN */
78  {
79    return x + y;
80  }
81  if (((hx - 0x3ff00000) | lx) == 0) /* x = 1.0 */
82  {
83    return atan (y);
84  }
85  m = ((hy < 0) ? 1 : 0) + ((hx < 0) ? 2 : 0); /* 2 * sign(x) + sign(y) */
86
87  /* when y = 0 */
88  if ((iy | ly) == 0)
89  {
90    switch (m)
91    {
92      case 0:
93      case 1:
94      {
95        return y; /* atan(+-0,+anything) = +-0 */
96      }
97      case 2:
98      {
99        return pi + tiny; /* atan(+0,-anything) = pi */
100      }
101      case 3:
102      {
103        return -pi - tiny; /* atan(-0,-anything) = -pi */
104      }
105    }
106  }
107  /* when x = 0 */
108  if ((ix | lx) == 0)
109  {
110    return (hy < 0) ? -pi_o_2 - tiny : pi_o_2 + tiny;
111  }
112
113  /* when x is INF */
114  if (ix == 0x7ff00000)
115  {
116    if (iy == 0x7ff00000)
117    {
118      switch (m)
119      {
120        case 0: /* atan(+INF,+INF) */
121        {
122          return pi_o_4 + tiny;
123        }
124        case 1: /* atan(-INF,+INF) */
125        {
126          return -pi_o_4 - tiny;
127        }
128        case 2: /* atan(+INF,-INF) */
129        {
130          return 3.0 * pi_o_4 + tiny;
131        }
132        case 3: /* atan(-INF,-INF) */
133        {
134          return -3.0 * pi_o_4 - tiny;
135        }
136      }
137    }
138    else
139    {
140      switch (m)
141      {
142        case 0: /* atan(+...,+INF) */
143        {
144          return zero;
145        }
146        case 1: /* atan(-...,+INF) */
147        {
148          return -zero;
149        }
150        case 2: /* atan(+...,-INF) */
151        {
152          return pi + tiny;
153        }
154        case 3: /* atan(-...,-INF) */
155        {
156          return -pi - tiny;
157        }
158      }
159    }
160  }
161  /* when y is INF */
162  if (iy == 0x7ff00000)
163  {
164    return (hy < 0) ? -pi_o_2 - tiny : pi_o_2 + tiny;
165  }
166
167  /* compute y / x */
168  k = (iy - ix) >> 20;
169  if (k > 60) /* |y / x| > 2**60 */
170  {
171    z.dbl = pi_o_2 + 0.5 * pi_lo;
172  }
173  else if (hx < 0 && k < -60) /* |y| / x < -2**60 */
174  {
175    z.dbl = 0.0;
176  }
177  else /* safe to do y / x */
178  {
179    z.dbl = atan (fabs (y / x));
180  }
181  switch (m)
182  {
183    case 0: /* atan(+,+) */
184    {
185      return z.dbl;
186    }
187    case 1: /* atan(-,+) */
188    {
189      z.as_int.hi ^= 0x80000000;
190      return z.dbl;
191    }
192    case 2: /* atan(+,-) */
193    {
194      return pi - (z.dbl - pi_lo);
195    }
196    /* case 3: */
197    default: /* atan(-,-) */
198    {
199      return (z.dbl - pi_lo) - pi;
200    }
201  }
202} /* atan2 */
203
204#undef tiny
205#undef zero
206#undef pi_o_4
207#undef pi_o_2
208#undef pi
209#undef pi_lo
210