/*
* CDDL HEADER START
*
* The contents of this file are subject to the terms of the
* Common Development and Distribution License (the "License").
* You may not use this file except in compliance with the License.
*
* You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
* See the License for the specific language governing permissions
* and limitations under the License.
*
* When distributing Covered Code, include this CDDL HEADER in each
* file and include the License file at usr/src/OPENSOLARIS.LICENSE.
* If applicable, add the following below this CDDL HEADER, with the
* fields enclosed by brackets "[]" replaced with your own identifying
* information: Portions Copyright [yyyy] [name of copyright owner]
*
* CDDL HEADER END
*/
/*
* Copyright 2011 Nexenta Systems, Inc. All rights reserved.
*/
/*
* Copyright 2006 Sun Microsystems, Inc. All rights reserved.
* Use is subject to license terms.
*/
#include "libm.h" /* __k_clog_rl/__k_atan2l */
#include "complex_wrapper.h"
#include "longdouble.h"
#if defined(__sparc)
#define HALF(x) ((int *) &x)[0] = 0
#define LAST(x) ((int *) &x)[0]
#endif
/* INDENT OFF */
static const long double
#if defined(__x86)
/* 43 significant bits, 21 trailing zeros */
#else /* sparc */
/* 0x3FF962E4 2FEFA39E F35793C7 00000000 */
ln2hil = 0.693147180559945309417231592858066493070671489074L,
ln2lol = 5.28600110075004828645286235820646730106802446566153e-25L,
#endif
/* INDENT ON */
/*
* Assuming |t[0]| > |t[1]| and |t[2]| > |t[3]|, sum4fpl subroutine
* compute t[0] + t[1] + t[2] + t[3] into two long double fp numbers.
*/
{
/*
* Rearrange ti so that |t1| >= |t2| >= |t3| >= |t4|
*/
} else {
}
} else
t3 = t;
}
/* summing r = t1 + t2 + t3 + t4 to w1 + w2 */
return (w1);
}
long double x, y, u, v, t, c, s, r;
x = LD_RE(z);
y = LD_IM(z);
u = LD_RE(w);
v = LD_IM(w);
j = 0;
if (v == zero) { /* z**(real) */
if (u == one) { /* (anything) ** 1 is itself */
} else if (u == zero) { /* (anything) ** 0 is 1 */
} else if (y == zero) { /* real ** real */
/* -x ** u is exp(i*pi*u)*pow(x,u) */
r = powl(-x, u);
sincospil(u, &s, &c);
} else
else {
if (x == zero)
r = fabsl(y);
else
t = atan2pil(y, x);
sincospil(t * u, &s, &c);
}
if (hx >= 0) {
sincospil(t * u, &s, &c);
} else if ((LAST(u) & 3) == 0) {
sincospil(t * u, &s, &c);
} else {
r = (hy >= 0)? u : -u;
t = -0.25L * r;
w1 = r + t;
}
/* |x| >= 1 or |u| < 16383 */
else /* special treatment */
j = 2;
if (j == 0) {
}
} else
j = 1;
if (j == 0)
return (ans);
}
/*
* non-zero imag part(s) with inf component(s) yields NaN
*/
} else {
k = 0; /* no scaling */
u *= 1.52587890625000000000e-05L;
v *= 1.52587890625000000000e-05L;
k = 1; /* scale u and v by 2**-16 */
}
/*
* Use similated higher precision arithmetic to compute:
* r = u * log(hypot(x, y)) - v * atan2(y, x)
* q = u * atan2(y, x) + v * log(hypot(x, y))
*/
/* compute q = u * atan2(y, x) + v * log(hypot(x, y)) */
if (j != 2) {
if (j == 1) { /* v = 0 */
w1 = b[0] + b[1];
} else {
}
if (k == 1) /* square j times */
for (i = 0; i < 10; i++) {
t1 = s * c;
c = (c + s) * (c - s);
}
}
/* compute r = u * (t1, t2) - v * (t3, t4) */
if (j == 1) { /* v = 0 */
w1 = b[0] + b[1];
} else {
}
/* scale back unless w1 is large enough to cause exception */
}
/* compute exp(w1 + w2) */
k = 0;
r = one;
else { /* compute exp(w1 + w2) */
t1 = (long double) k;
}
if (c != zero) c *= r;
if (s != zero) s *= r;
if (k != 0) {
c = scalbnl(c, k);
s = scalbnl(s, k);
}
}
return (ans);
}