cabs.c revision 25c28e83beb90e7c80452a7c818c5e6f73a07dc8
/*
* 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 2005 Sun Microsystems, Inc. All rights reserved.
* Use is subject to license terms.
*/
#include "libm_synonyms.h"
#include <math.h>
#include "complex_wrapper.h"
/*
* If C were the only standard we cared about, cabs could just call
* hypot. Unfortunately, various other standards say that hypot must
* Since cabs should do neither of these things, we have to either
* make hypot a wrapper on another internal function or duplicate
* the hypot implementation here. I've chosen to do the latter.
*/
static const double
zero = 0.0,
double
{
x = D_RE(z);
y = D_IM(z);
/* force ax = |x| ~>~ ay = |y| */
i = ix;
iy = i;
i = lx;
ly = i;
} else {
}
if (nx >= 0x5f3) {
/* x >= 2^500 (x*x or y*y may overflow) */
if (nx == 0x7ff) {
/* inf or NaN, signal of sNaN */
else
} else if (j > 32) {
/* x >> y */
if (j <= 53)
return (ax);
}
iscale = 2;
} else if (ny < 0x23d) {
/* y < 2^-450 (x*x or y*y may underflow) */
return (ay);
return (ax);
if (j > 53) /* x >> y */
iscale = 1;
if (nx == 0) {
return (ax);
} else {
}
if (ny == 0) {
} else {
}
if (j > 32) {
/* x >> y */
if (j <= 53)
}
} else if (j > 32) {
/* x >> y */
if (j <= 53)
}
/*
* First check rounding mode by comparing onep1u*onep1u with onep1u
* + twom53. Make sure the computation is done at run-time.
*/
/* round-to-zero, positive, negative mode */
/* magic formula with less than an ulp error */
} else {
/* round-to-nearest mode */
if (w > ay) {
} else {
}
}
if (iscale > 0) {
if (iscale == 1)
else
}
return (ax);
}