1 /* sin (sine) function with 'long double' argument. 2 3 Copyright (C) 2003-2021 Free Software Foundation, Inc. 4 5 This file is free software: you can redistribute it and/or modify 6 it under the terms of the GNU Lesser General Public License as 7 published by the Free Software Foundation; either version 3 of the 8 License, or (at your option) any later version. 9 10 This file is distributed in the hope that it will be useful, 11 but WITHOUT ANY WARRANTY; without even the implied warranty of 12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 13 GNU Lesser General Public License for more details. 14 15 You should have received a copy of the GNU Lesser General Public License 16 along with this program. If not, see <https://www.gnu.org/licenses/>. */ 17 18 /* s_sinl.c -- long double version of s_sin.c. 19 * Conversion to long double by Jakub Jelinek, jj@ultra.linux.cz. 20 */ 21 22 /* 23 * ==================================================== 24 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. 25 * 26 * Developed at SunPro, a Sun Microsystems, Inc. business. 27 * Permission to use, copy, modify, and distribute this 28 * software is freely granted, provided that this notice 29 * is preserved. 30 * ==================================================== 31 */ 32 33 #include <config.h> 34 35 /* Specification. */ 36 #include <math.h> 37 38 #if HAVE_SAME_LONG_DOUBLE_AS_DOUBLE 39 40 long double 41 sinl (long double x) /* */ 42 { 43 return sin (x); 44 } 45 46 #else 47 48 /* Code based on glibc/sysdeps/ieee754/ldbl-128/s_sinl.c. */ 49 50 /* sinl(x) 51 * Return sine function of x. 52 * 53 * kernel function: 54 * __kernel_sinl ... sine function on [-pi/4,pi/4] 55 * __kernel_cosl ... cosine function on [-pi/4,pi/4] 56 * __ieee754_rem_pio2l ... argument reduction routine 57 * 58 * Method. 59 * Let S,C and T denote the sin, cos and tan respectively on 60 * [-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2 61 * in [-pi/4 , +pi/4], and let n = k mod 4. 62 * We have 63 * 64 * n sin(x) cos(x) tan(x) 65 * ---------------------------------------------------------- 66 * 0 S C T 67 * 1 C -S -1/T 68 * 2 -S -C T 69 * 3 -C S -1/T 70 * ---------------------------------------------------------- 71 * 72 * Special cases: 73 * Let trig be any of sin, cos, or tan. 74 * trig(+-INF) is NaN, with signals; 75 * trig(NaN) is that NaN; 76 * 77 * Accuracy: 78 * TRIG(x) returns trig(x) nearly rounded 79 */ 80 81 # include "trigl.h" 82 83 long double 84 sinl (long double x) /* */ 85 { 86 long double y[2], z = 0.0L; 87 int n; 88 89 /* sinl(NaN) is NaN */ 90 if (isnanl (x)) 91 return x; 92 93 /* |x| ~< pi/4 */ 94 if (x >= -0.7853981633974483096156608458198757210492 95 && x <= 0.7853981633974483096156608458198757210492) 96 return kernel_sinl (x, z, 0); 97 98 /* sinl(Inf) is NaN, sinl(0) is 0 */ 99 else if (x + x == x) 100 return x - x; /* NaN */ 101 102 /* argument reduction needed */ 103 else 104 { 105 n = ieee754_rem_pio2l (x, y); 106 switch (n & 3) 107 { 108 case 0: 109 return kernel_sinl (y[0], y[1], 1); 110 case 1: 111 return kernel_cosl (y[0], y[1]); 112 case 2: 113 return -kernel_sinl (y[0], y[1], 1); 114 default: 115 return -kernel_cosl (y[0], y[1]); 116 } 117 } 118 } 119 120 #endif 121 122 #if 0 123 int 124 main (void) /* */ 125 { 126 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492)); 127 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *29)); 128 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *2)); 129 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *30)); 130 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *4)); 131 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *32)); 132 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *2/3)); 133 printf ("%.16Lg\n", sinl (0.7853981633974483096156608458198757210492 *4/3)); 134 } 135 #endif