test-invtrig.c revision 181206
1181076Sdas/*- 2181076Sdas * Copyright (c) 2008 David Schultz <das@FreeBSD.org> 3181076Sdas * All rights reserved. 4181076Sdas * 5181076Sdas * Redistribution and use in source and binary forms, with or without 6181076Sdas * modification, are permitted provided that the following conditions 7181076Sdas * are met: 8181076Sdas * 1. Redistributions of source code must retain the above copyright 9181076Sdas * notice, this list of conditions and the following disclaimer. 10181076Sdas * 2. Redistributions in binary form must reproduce the above copyright 11181076Sdas * notice, this list of conditions and the following disclaimer in the 12181076Sdas * documentation and/or other materials provided with the distribution. 13181076Sdas * 14181076Sdas * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND 15181076Sdas * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE 16181076Sdas * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE 17181076Sdas * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE 18181076Sdas * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 19181076Sdas * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS 20181076Sdas * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 21181076Sdas * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 22181076Sdas * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY 23181076Sdas * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF 24181076Sdas * SUCH DAMAGE. 25181076Sdas */ 26181076Sdas 27181076Sdas/* 28181076Sdas * Tests for corner cases in the inverse trigonometric functions. Some 29181076Sdas * accuracy tests are included as well, but these are very basic 30181076Sdas * sanity checks, not intended to be comprehensive. 31181076Sdas */ 32181076Sdas 33181076Sdas#include <sys/cdefs.h> 34181076Sdas__FBSDID("$FreeBSD: head/tools/regression/lib/msun/test-invtrig.c 181206 2008-08-02 19:21:33Z das $"); 35181076Sdas 36181076Sdas#include <assert.h> 37181076Sdas#include <fenv.h> 38181076Sdas#include <float.h> 39181076Sdas#include <math.h> 40181076Sdas#include <stdio.h> 41181076Sdas 42181076Sdas#define ALL_STD_EXCEPT (FE_DIVBYZERO | FE_INEXACT | FE_INVALID | \ 43181076Sdas FE_OVERFLOW | FE_UNDERFLOW) 44181076Sdas 45181076Sdas#define LEN(a) (sizeof(a) / sizeof((a)[0])) 46181076Sdas 47181076Sdas#pragma STDC FENV_ACCESS ON 48181076Sdas 49181076Sdas/* 50181076Sdas * Test that a function returns the correct value and sets the 51181076Sdas * exception flags correctly. A tolerance specifying the maximum 52181076Sdas * relative error allowed may be specified. For the 'testall' 53181076Sdas * functions, the tolerance is specified in ulps. 54181076Sdas * 55181076Sdas * These are macros instead of functions so that assert provides more 56181076Sdas * meaningful error messages. 57181076Sdas */ 58181076Sdas#define test_tol(func, x, result, tol, excepts) do { \ 59181076Sdas volatile long double _in = (x), _out = (result); \ 60181076Sdas assert(feclearexcept(FE_ALL_EXCEPT) == 0); \ 61181076Sdas assert(fpequal(func(_in), _out, (tol))); \ 62181076Sdas assert((func, fetestexcept(ALL_STD_EXCEPT) == (excepts))); \ 63181076Sdas} while (0) 64181076Sdas#define test(func, x, result, excepts) \ 65181076Sdas test_tol(func, (x), (result), 0, (excepts)) 66181076Sdas 67181076Sdas#define testall_tol(prefix, x, result, tol, excepts) do { \ 68181076Sdas test_tol(prefix, (double)(x), (double)(result), \ 69181076Sdas (tol) * ldexp(1.0, 1 - DBL_MANT_DIG), (excepts)); \ 70181076Sdas test_tol(prefix##f, (float)(x), (float)(result), \ 71181076Sdas (tol) * ldexpf(1.0, 1 - FLT_MANT_DIG), (excepts)); \ 72181076Sdas test_tol(prefix##l, (x), (result), \ 73181076Sdas (tol) * ldexpl(1.0, 1 - LDBL_MANT_DIG), (excepts)); \ 74181076Sdas} while (0) 75181076Sdas#define testall(prefix, x, result, excepts) \ 76181076Sdas testall_tol(prefix, (x), (result), 0, (excepts)) 77181076Sdas 78181076Sdas#define test2_tol(func, y, x, result, tol, excepts) do { \ 79181076Sdas volatile long double _iny = (y), _inx = (x), _out = (result); \ 80181076Sdas assert(feclearexcept(FE_ALL_EXCEPT) == 0); \ 81181076Sdas assert(fpequal(func(_iny, _inx), _out, (tol))); \ 82181076Sdas assert((func, fetestexcept(ALL_STD_EXCEPT) == (excepts))); \ 83181076Sdas} while (0) 84181076Sdas#define test2(func, y, x, result, excepts) \ 85181076Sdas test2_tol(func, (y), (x), (result), 0, (excepts)) 86181076Sdas 87181076Sdas#define testall2_tol(prefix, y, x, result, tol, excepts) do { \ 88181076Sdas test2_tol(prefix, (double)(y), (double)(x), (double)(result), \ 89181076Sdas (tol) * ldexp(1.0, 1 - DBL_MANT_DIG), (excepts)); \ 90181076Sdas test2_tol(prefix##f, (float)(y), (float)(x), (float)(result), \ 91181076Sdas (tol) * ldexpf(1.0, 1 - FLT_MANT_DIG), (excepts)); \ 92181076Sdas test2_tol(prefix##l, (y), (x), (result), \ 93181076Sdas (tol) * ldexpl(1.0, 1 - LDBL_MANT_DIG), (excepts)); \ 94181076Sdas} while (0) 95181076Sdas#define testall2(prefix, y, x, result, excepts) \ 96181076Sdas testall2_tol(prefix, (y), (x), (result), 0, (excepts)) 97181076Sdas 98181076Sdaslong double 99181076Sdaspi = 3.14159265358979323846264338327950280e+00L, 100181076Sdaspio3 = 1.04719755119659774615421446109316766e+00L, 101181076Sdasc3pi = 9.42477796076937971538793014983850839e+00L, 102181076Sdasc5pi = 1.57079632679489661923132169163975140e+01L, 103181076Sdasc7pi = 2.19911485751285526692385036829565196e+01L, 104181076Sdasc5pio3 = 5.23598775598298873077107230546583851e+00L, 105181076Sdassqrt2m1 = 4.14213562373095048801688724209698081e-01L; 106181076Sdas 107181076Sdas/* 108181076Sdas * Determine whether x and y are equal to within a relative error of tol, 109181076Sdas * with two special rules: 110181076Sdas * +0.0 != -0.0 111181076Sdas * NaN == NaN 112181076Sdas */ 113181076Sdasint 114181076Sdasfpequal(long double x, long double y, long double tol) 115181076Sdas{ 116181076Sdas fenv_t env; 117181076Sdas int ret; 118181076Sdas 119181076Sdas if (isnan(x) && isnan(y)) 120181076Sdas return (1); 121181076Sdas if (signbit(x) != signbit(y)) 122181076Sdas return (0); 123181076Sdas if (x == y) 124181076Sdas return (1); 125181076Sdas if (tol == 0) 126181076Sdas return (0); 127181076Sdas 128181076Sdas /* Hard case: need to check the tolerance. */ 129181076Sdas feholdexcept(&env); 130181076Sdas ret = fabsl(x - y) <= fabsl(y * tol); 131181076Sdas fesetenv(&env); 132181076Sdas return (ret); 133181076Sdas} 134181076Sdas 135181076Sdas/* 136181076Sdas * Test special case inputs in asin(), acos() and atan(): signed 137181076Sdas * zeroes, infinities, and NaNs. 138181076Sdas */ 139181076Sdasstatic void 140181076Sdastest_special(void) 141181076Sdas{ 142181076Sdas 143181076Sdas testall(asin, 0.0, 0.0, 0); 144181076Sdas testall(acos, 0.0, pi / 2, FE_INEXACT); 145181076Sdas testall(atan, 0.0, 0.0, 0); 146181076Sdas testall(asin, -0.0, -0.0, 0); 147181076Sdas testall(acos, -0.0, pi / 2, FE_INEXACT); 148181076Sdas testall(atan, -0.0, -0.0, 0); 149181076Sdas 150181076Sdas testall(asin, INFINITY, NAN, FE_INVALID); 151181076Sdas testall(acos, INFINITY, NAN, FE_INVALID); 152181076Sdas testall(atan, INFINITY, pi / 2, FE_INEXACT); 153181076Sdas testall(asin, -INFINITY, NAN, FE_INVALID); 154181076Sdas testall(acos, -INFINITY, NAN, FE_INVALID); 155181076Sdas testall(atan, -INFINITY, -pi / 2, FE_INEXACT); 156181076Sdas 157181076Sdas testall(asin, NAN, NAN, 0); 158181076Sdas testall(acos, NAN, NAN, 0); 159181076Sdas testall(atan, NAN, NAN, 0); 160181076Sdas} 161181076Sdas 162181076Sdas/* 163181076Sdas * Test special case inputs in atan2(), where the exact value of y/x is 164181076Sdas * zero or non-finite. 165181076Sdas */ 166181076Sdasstatic void 167181076Sdastest_special_atan2(void) 168181076Sdas{ 169181076Sdas long double z; 170181076Sdas int e; 171181076Sdas 172181076Sdas testall2(atan2, 0.0, -0.0, pi, FE_INEXACT); 173181076Sdas testall2(atan2, -0.0, -0.0, -pi, FE_INEXACT); 174181076Sdas testall2(atan2, 0.0, 0.0, 0.0, 0); 175181076Sdas testall2(atan2, -0.0, 0.0, -0.0, 0); 176181076Sdas 177181076Sdas testall2(atan2, INFINITY, -INFINITY, c3pi / 4, FE_INEXACT); 178181076Sdas testall2(atan2, -INFINITY, -INFINITY, -c3pi / 4, FE_INEXACT); 179181076Sdas testall2(atan2, INFINITY, INFINITY, pi / 4, FE_INEXACT); 180181076Sdas testall2(atan2, -INFINITY, INFINITY, -pi / 4, FE_INEXACT); 181181076Sdas 182181076Sdas /* Tests with one input in the range (0, Inf]. */ 183181076Sdas z = 1.23456789L; 184181076Sdas for (e = FLT_MIN_EXP - FLT_MANT_DIG; e <= FLT_MAX_EXP; e++) { 185181076Sdas test2(atan2f, 0.0, ldexpf(z, e), 0.0, 0); 186181076Sdas test2(atan2f, -0.0, ldexpf(z, e), -0.0, 0); 187181076Sdas test2(atan2f, 0.0, ldexpf(-z, e), (float)pi, FE_INEXACT); 188181076Sdas test2(atan2f, -0.0, ldexpf(-z, e), (float)-pi, FE_INEXACT); 189181076Sdas test2(atan2f, ldexpf(z, e), 0.0, (float)pi / 2, FE_INEXACT); 190181076Sdas test2(atan2f, ldexpf(z, e), -0.0, (float)pi / 2, FE_INEXACT); 191181076Sdas test2(atan2f, ldexpf(-z, e), 0.0, (float)-pi / 2, FE_INEXACT); 192181076Sdas test2(atan2f, ldexpf(-z, e), -0.0, (float)-pi / 2, FE_INEXACT); 193181076Sdas } 194181076Sdas for (e = DBL_MIN_EXP - DBL_MANT_DIG; e <= DBL_MAX_EXP; e++) { 195181076Sdas test2(atan2, 0.0, ldexp(z, e), 0.0, 0); 196181076Sdas test2(atan2, -0.0, ldexp(z, e), -0.0, 0); 197181076Sdas test2(atan2, 0.0, ldexp(-z, e), (double)pi, FE_INEXACT); 198181076Sdas test2(atan2, -0.0, ldexp(-z, e), (double)-pi, FE_INEXACT); 199181076Sdas test2(atan2, ldexp(z, e), 0.0, (double)pi / 2, FE_INEXACT); 200181076Sdas test2(atan2, ldexp(z, e), -0.0, (double)pi / 2, FE_INEXACT); 201181076Sdas test2(atan2, ldexp(-z, e), 0.0, (double)-pi / 2, FE_INEXACT); 202181076Sdas test2(atan2, ldexp(-z, e), -0.0, (double)-pi / 2, FE_INEXACT); 203181076Sdas } 204181076Sdas for (e = LDBL_MIN_EXP - LDBL_MANT_DIG; e <= LDBL_MAX_EXP; e++) { 205181076Sdas test2(atan2l, 0.0, ldexpl(z, e), 0.0, 0); 206181076Sdas test2(atan2l, -0.0, ldexpl(z, e), -0.0, 0); 207181076Sdas test2(atan2l, 0.0, ldexpl(-z, e), pi, FE_INEXACT); 208181076Sdas test2(atan2l, -0.0, ldexpl(-z, e), -pi, FE_INEXACT); 209181076Sdas test2(atan2l, ldexpl(z, e), 0.0, pi / 2, FE_INEXACT); 210181076Sdas test2(atan2l, ldexpl(z, e), -0.0, pi / 2, FE_INEXACT); 211181076Sdas test2(atan2l, ldexpl(-z, e), 0.0, -pi / 2, FE_INEXACT); 212181076Sdas test2(atan2l, ldexpl(-z, e), -0.0, -pi / 2, FE_INEXACT); 213181076Sdas } 214181076Sdas 215181076Sdas /* Tests with one input in the range (0, Inf). */ 216181076Sdas for (e = FLT_MIN_EXP - FLT_MANT_DIG; e <= FLT_MAX_EXP - 1; e++) { 217181076Sdas test2(atan2f, ldexpf(z, e), INFINITY, 0.0, 0); 218181076Sdas test2(atan2f, ldexpf(-z,e), INFINITY, -0.0, 0); 219181076Sdas test2(atan2f, ldexpf(z, e), -INFINITY, (float)pi, FE_INEXACT); 220181076Sdas test2(atan2f, ldexpf(-z,e), -INFINITY, (float)-pi, FE_INEXACT); 221181076Sdas test2(atan2f, INFINITY, ldexpf(z,e), (float)pi/2, FE_INEXACT); 222181076Sdas test2(atan2f, INFINITY, ldexpf(-z,e), (float)pi/2, FE_INEXACT); 223181076Sdas test2(atan2f, -INFINITY, ldexpf(z,e), (float)-pi/2,FE_INEXACT); 224181076Sdas test2(atan2f, -INFINITY, ldexpf(-z,e),(float)-pi/2,FE_INEXACT); 225181076Sdas } 226181076Sdas for (e = DBL_MIN_EXP - DBL_MANT_DIG; e <= DBL_MAX_EXP - 1; e++) { 227181076Sdas test2(atan2, ldexp(z, e), INFINITY, 0.0, 0); 228181076Sdas test2(atan2, ldexp(-z,e), INFINITY, -0.0, 0); 229181076Sdas test2(atan2, ldexp(z, e), -INFINITY, (double)pi, FE_INEXACT); 230181076Sdas test2(atan2, ldexp(-z,e), -INFINITY, (double)-pi, FE_INEXACT); 231181076Sdas test2(atan2, INFINITY, ldexp(z,e), (double)pi/2, FE_INEXACT); 232181076Sdas test2(atan2, INFINITY, ldexp(-z,e), (double)pi/2, FE_INEXACT); 233181076Sdas test2(atan2, -INFINITY, ldexp(z,e), (double)-pi/2,FE_INEXACT); 234181076Sdas test2(atan2, -INFINITY, ldexp(-z,e),(double)-pi/2,FE_INEXACT); 235181076Sdas } 236181076Sdas for (e = LDBL_MIN_EXP - LDBL_MANT_DIG; e <= LDBL_MAX_EXP - 1; e++) { 237181076Sdas test2(atan2l, ldexpl(z, e), INFINITY, 0.0, 0); 238181076Sdas test2(atan2l, ldexpl(-z,e), INFINITY, -0.0, 0); 239181076Sdas test2(atan2l, ldexpl(z, e), -INFINITY, pi, FE_INEXACT); 240181076Sdas test2(atan2l, ldexpl(-z,e), -INFINITY, -pi, FE_INEXACT); 241181076Sdas test2(atan2l, INFINITY, ldexpl(z, e), pi / 2, FE_INEXACT); 242181076Sdas test2(atan2l, INFINITY, ldexpl(-z, e), pi / 2, FE_INEXACT); 243181076Sdas test2(atan2l, -INFINITY, ldexpl(z, e), -pi / 2, FE_INEXACT); 244181076Sdas test2(atan2l, -INFINITY, ldexpl(-z, e), -pi / 2, FE_INEXACT); 245181076Sdas } 246181076Sdas} 247181076Sdas 248181076Sdas/* 249181076Sdas * Test various inputs to asin(), acos() and atan() and verify that the 250181076Sdas * results are accurate to within 1 ulp. 251181076Sdas */ 252181076Sdasstatic void 253181076Sdastest_accuracy(void) 254181076Sdas{ 255181076Sdas 256181076Sdas /* We expect correctly rounded results for these basic cases. */ 257181076Sdas testall(asin, 1.0, pi / 2, FE_INEXACT); 258181076Sdas testall(acos, 1.0, 0, 0); 259181076Sdas testall(atan, 1.0, pi / 4, FE_INEXACT); 260181076Sdas testall(asin, -1.0, -pi / 2, FE_INEXACT); 261181076Sdas testall(acos, -1.0, pi, FE_INEXACT); 262181076Sdas testall(atan, -1.0, -pi / 4, FE_INEXACT); 263181076Sdas 264181076Sdas /* 265181076Sdas * Here we expect answers to be within 1 ulp, although inexactness 266181076Sdas * in the input, combined with double rounding, could cause larger 267181076Sdas * errors. 268181076Sdas */ 269181076Sdas 270181076Sdas testall_tol(asin, sqrtl(2) / 2, pi / 4, 1, FE_INEXACT); 271181076Sdas testall_tol(acos, sqrtl(2) / 2, pi / 4, 1, FE_INEXACT); 272181076Sdas testall_tol(asin, -sqrtl(2) / 2, -pi / 4, 1, FE_INEXACT); 273181076Sdas testall_tol(acos, -sqrtl(2) / 2, c3pi / 4, 1, FE_INEXACT); 274181076Sdas 275181076Sdas testall_tol(asin, sqrtl(3) / 2, pio3, 1, FE_INEXACT); 276181076Sdas testall_tol(acos, sqrtl(3) / 2, pio3 / 2, 1, FE_INEXACT); 277181076Sdas testall_tol(atan, sqrtl(3), pio3, 1, FE_INEXACT); 278181076Sdas testall_tol(asin, -sqrtl(3) / 2, -pio3, 1, FE_INEXACT); 279181076Sdas testall_tol(acos, -sqrtl(3) / 2, c5pio3 / 2, 1, FE_INEXACT); 280181076Sdas testall_tol(atan, -sqrtl(3), -pio3, 1, FE_INEXACT); 281181076Sdas 282181076Sdas testall_tol(atan, sqrt2m1, pi / 8, 1, FE_INEXACT); 283181076Sdas testall_tol(atan, -sqrt2m1, -pi / 8, 1, FE_INEXACT); 284181076Sdas} 285181076Sdas 286181076Sdas/* 287181076Sdas * Test inputs to atan2() where x is a power of 2. These are easy cases 288181076Sdas * because y/x is exact. 289181076Sdas */ 290181076Sdasstatic void 291181076Sdastest_p2x_atan2(void) 292181076Sdas{ 293181076Sdas 294181076Sdas testall2(atan2, 1.0, 1.0, pi / 4, FE_INEXACT); 295181076Sdas testall2(atan2, 1.0, -1.0, c3pi / 4, FE_INEXACT); 296181076Sdas testall2(atan2, -1.0, 1.0, -pi / 4, FE_INEXACT); 297181076Sdas testall2(atan2, -1.0, -1.0, -c3pi / 4, FE_INEXACT); 298181076Sdas 299181076Sdas testall2_tol(atan2, sqrt2m1 * 2, 2.0, pi / 8, 1, FE_INEXACT); 300181076Sdas testall2_tol(atan2, sqrt2m1 * 2, -2.0, c7pi / 8, 1, FE_INEXACT); 301181076Sdas testall2_tol(atan2, -sqrt2m1 * 2, 2.0, -pi / 8, 1, FE_INEXACT); 302181076Sdas testall2_tol(atan2, -sqrt2m1 * 2, -2.0, -c7pi / 8, 1, FE_INEXACT); 303181076Sdas 304181076Sdas testall2_tol(atan2, sqrtl(3) * 0.5, 0.5, pio3, 1, FE_INEXACT); 305181076Sdas testall2_tol(atan2, sqrtl(3) * 0.5, -0.5, pio3 * 2, 1, FE_INEXACT); 306181076Sdas testall2_tol(atan2, -sqrtl(3) * 0.5, 0.5, -pio3, 1, FE_INEXACT); 307181076Sdas testall2_tol(atan2, -sqrtl(3) * 0.5, -0.5, -pio3 * 2, 1, FE_INEXACT); 308181076Sdas} 309181076Sdas 310181076Sdas/* 311181076Sdas * Test inputs very close to 0. 312181076Sdas */ 313181076Sdasstatic void 314181076Sdastest_tiny(void) 315181076Sdas{ 316181076Sdas float tiny = 0x1.23456p-120f; 317181076Sdas 318181076Sdas testall(asin, tiny, tiny, FE_INEXACT); 319181076Sdas testall(acos, tiny, pi / 2, FE_INEXACT); 320181076Sdas testall(atan, tiny, tiny, FE_INEXACT); 321181076Sdas 322181076Sdas testall(asin, -tiny, -tiny, FE_INEXACT); 323181076Sdas testall(acos, -tiny, pi / 2, FE_INEXACT); 324181076Sdas testall(atan, -tiny, -tiny, FE_INEXACT); 325181076Sdas 326181076Sdas /* Test inputs to atan2() that would cause y/x to underflow. */ 327181076Sdas test2(atan2f, 0x1.0p-100, 0x1.0p100, 0.0, FE_INEXACT | FE_UNDERFLOW); 328181076Sdas test2(atan2, 0x1.0p-1000, 0x1.0p1000, 0.0, FE_INEXACT | FE_UNDERFLOW); 329181076Sdas test2(atan2l, ldexpl(1.0, 100 - LDBL_MAX_EXP), 330181076Sdas ldexpl(1.0, LDBL_MAX_EXP - 100), 0.0, FE_INEXACT | FE_UNDERFLOW); 331181076Sdas test2(atan2f, -0x1.0p-100, 0x1.0p100, -0.0, FE_INEXACT | FE_UNDERFLOW); 332181076Sdas test2(atan2, -0x1.0p-1000, 0x1.0p1000, -0.0, FE_INEXACT | FE_UNDERFLOW); 333181076Sdas test2(atan2l, -ldexpl(1.0, 100 - LDBL_MAX_EXP), 334181076Sdas ldexpl(1.0, LDBL_MAX_EXP - 100), -0.0, FE_INEXACT | FE_UNDERFLOW); 335181076Sdas test2(atan2f, 0x1.0p-100, -0x1.0p100, (float)pi, FE_INEXACT); 336181076Sdas test2(atan2, 0x1.0p-1000, -0x1.0p1000, (double)pi, FE_INEXACT); 337181076Sdas test2(atan2l, ldexpl(1.0, 100 - LDBL_MAX_EXP), 338181076Sdas -ldexpl(1.0, LDBL_MAX_EXP - 100), pi, FE_INEXACT); 339181076Sdas test2(atan2f, -0x1.0p-100, -0x1.0p100, (float)-pi, FE_INEXACT); 340181076Sdas test2(atan2, -0x1.0p-1000, -0x1.0p1000, (double)-pi, FE_INEXACT); 341181076Sdas test2(atan2l, -ldexpl(1.0, 100 - LDBL_MAX_EXP), 342181076Sdas -ldexpl(1.0, LDBL_MAX_EXP - 100), -pi, FE_INEXACT); 343181076Sdas} 344181076Sdas 345181076Sdas/* 346181076Sdas * Test very large inputs to atan(). 347181076Sdas */ 348181076Sdasstatic void 349181076Sdastest_atan_huge(void) 350181076Sdas{ 351181076Sdas float huge = 0x1.23456p120; 352181076Sdas 353181076Sdas testall(atan, huge, pi / 2, FE_INEXACT); 354181076Sdas testall(atan, -huge, -pi / 2, FE_INEXACT); 355181076Sdas 356181076Sdas /* Test inputs to atan2() that would cause y/x to overflow. */ 357181076Sdas test2(atan2f, 0x1.0p100, 0x1.0p-100, (float)pi / 2, FE_INEXACT); 358181076Sdas test2(atan2, 0x1.0p1000, 0x1.0p-1000, (double)pi / 2, FE_INEXACT); 359181076Sdas test2(atan2l, ldexpl(1.0, LDBL_MAX_EXP - 100), 360181076Sdas ldexpl(1.0, 100 - LDBL_MAX_EXP), pi / 2, FE_INEXACT); 361181076Sdas test2(atan2f, -0x1.0p100, 0x1.0p-100, (float)-pi / 2, FE_INEXACT); 362181076Sdas test2(atan2, -0x1.0p1000, 0x1.0p-1000, (double)-pi / 2, FE_INEXACT); 363181076Sdas test2(atan2l, -ldexpl(1.0, LDBL_MAX_EXP - 100), 364181076Sdas ldexpl(1.0, 100 - LDBL_MAX_EXP), -pi / 2, FE_INEXACT); 365181206Sdas 366181206Sdas test2(atan2f, 0x1.0p100, -0x1.0p-100, (float)pi / 2, FE_INEXACT); 367181206Sdas test2(atan2, 0x1.0p1000, -0x1.0p-1000, (double)pi / 2, FE_INEXACT); 368181206Sdas test2(atan2l, ldexpl(1.0, LDBL_MAX_EXP - 100), 369181206Sdas -ldexpl(1.0, 100 - LDBL_MAX_EXP), pi / 2, FE_INEXACT); 370181206Sdas test2(atan2f, -0x1.0p100, -0x1.0p-100, (float)-pi / 2, FE_INEXACT); 371181206Sdas test2(atan2, -0x1.0p1000, -0x1.0p-1000, (double)-pi / 2, FE_INEXACT); 372181206Sdas test2(atan2l, -ldexpl(1.0, LDBL_MAX_EXP - 100), 373181206Sdas -ldexpl(1.0, 100 - LDBL_MAX_EXP), -pi / 2, FE_INEXACT); 374181076Sdas} 375181076Sdas 376181076Sdas/* 377181076Sdas * Test that sin(asin(x)) == x, and similarly for acos() and atan(). 378181076Sdas * You need to have a working sinl(), cosl(), and tanl() for these 379181076Sdas * tests to pass. 380181076Sdas */ 381181076Sdasstatic long double 382181076Sdassinasinf(float x) 383181076Sdas{ 384181076Sdas 385181076Sdas return (sinl(asinf(x))); 386181076Sdas} 387181076Sdas 388181076Sdasstatic long double 389181076Sdassinasin(double x) 390181076Sdas{ 391181076Sdas 392181076Sdas return (sinl(asin(x))); 393181076Sdas} 394181076Sdas 395181076Sdasstatic long double 396181076Sdassinasinl(long double x) 397181076Sdas{ 398181076Sdas 399181076Sdas return (sinl(asinl(x))); 400181076Sdas} 401181076Sdas 402181076Sdasstatic long double 403181076Sdascosacosf(float x) 404181076Sdas{ 405181076Sdas 406181076Sdas return (cosl(acosf(x))); 407181076Sdas} 408181076Sdas 409181076Sdasstatic long double 410181076Sdascosacos(double x) 411181076Sdas{ 412181076Sdas 413181076Sdas return (cosl(acos(x))); 414181076Sdas} 415181076Sdas 416181076Sdasstatic long double 417181076Sdascosacosl(long double x) 418181076Sdas{ 419181076Sdas 420181076Sdas return (cosl(acosl(x))); 421181076Sdas} 422181076Sdas 423181076Sdasstatic long double 424181076Sdastanatanf(float x) 425181076Sdas{ 426181076Sdas 427181076Sdas return (tanl(atanf(x))); 428181076Sdas} 429181076Sdas 430181076Sdasstatic long double 431181076Sdastanatan(double x) 432181076Sdas{ 433181076Sdas 434181076Sdas return (tanl(atan(x))); 435181076Sdas} 436181076Sdas 437181076Sdasstatic long double 438181076Sdastanatanl(long double x) 439181076Sdas{ 440181076Sdas 441181076Sdas return (tanl(atanl(x))); 442181076Sdas} 443181076Sdas 444181076Sdasstatic void 445181076Sdastest_inverse(void) 446181076Sdas{ 447181076Sdas float i; 448181076Sdas 449181076Sdas for (i = -1; i <= 1; i += 0x1.0p-12f) { 450181076Sdas testall_tol(sinasin, i, i, 2, i == 0 ? 0 : FE_INEXACT); 451181076Sdas /* The relative error for cosacos is very large near x=0. */ 452181076Sdas if (fabsf(i) > 0x1.0p-4f) 453181076Sdas testall_tol(cosacos, i, i, 16, i == 1 ? 0 : FE_INEXACT); 454181076Sdas testall_tol(tanatan, i, i, 2, i == 0 ? 0 : FE_INEXACT); 455181076Sdas } 456181076Sdas} 457181076Sdas 458181076Sdasint 459181076Sdasmain(int argc, char *argv[]) 460181076Sdas{ 461181076Sdas 462181076Sdas printf("1..7\n"); 463181076Sdas 464181076Sdas test_special(); 465181076Sdas printf("ok 1 - special\n"); 466181076Sdas 467181076Sdas test_special_atan2(); 468181076Sdas printf("ok 2 - atan2 special\n"); 469181076Sdas 470181076Sdas test_accuracy(); 471181076Sdas printf("ok 3 - accuracy\n"); 472181076Sdas 473181076Sdas test_p2x_atan2(); 474181076Sdas printf("ok 4 - atan2 p2x\n"); 475181076Sdas 476181076Sdas test_tiny(); 477181076Sdas printf("ok 5 - tiny inputs\n"); 478181076Sdas 479181076Sdas test_atan_huge(); 480181076Sdas printf("ok 6 - atan huge inputs\n"); 481181076Sdas 482181076Sdas test_inverse(); 483181076Sdas printf("ok 7 - inverse\n"); 484181076Sdas 485181076Sdas return (0); 486181076Sdas} 487