1 /* Implementation of the degree trignometric functions COSD, SIND, TAND. 2 Copyright (C) 2020-2021 Free Software Foundation, Inc. 3 Contributed by Steven G. Kargl <kargl@gcc.gnu.org> 4 5 This file is part of the GNU Fortran runtime library (libgfortran). 6 7 Libgfortran is free software; you can redistribute it and/or 8 modify it under the terms of the GNU General Public 9 License as published by the Free Software Foundation; either 10 version 3 of the License, or (at your option) any later version. 11 12 Libgfortran is distributed in the hope that it will be useful, 13 but WITHOUT ANY WARRANTY; without even the implied warranty of 14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 GNU General Public License for more details. 16 17 Under Section 7 of GPL version 3, you are granted additional 18 permissions described in the GCC Runtime Library Exception, version 19 3.1, as published by the Free Software Foundation. 20 21 You should have received a copy of the GNU General Public License and 22 a copy of the GCC Runtime Library Exception along with this program; 23 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see 24 <http://www.gnu.org/licenses/>. */ 25 26 #include "libgfortran.h" 27 28 #include <math.h> 29 30 /* Body of library functions which are cannot be implemented on the current 31 * platform because it lacks a capability, such as an underlying trigonometric 32 * function (sin, cos, tan) or C99 floating-point function (fabs, fmod). */ 33 #define STRINGIFY_EXPAND(x) #x 34 #define ERROR_RETURN(f, k, x) runtime_error (#f " is unavailable for" \ 35 " REAL(KIND=" STRINGIFY_EXPAND(k) ") because the system math library" \ 36 " lacks support for it"); \ 37 RETURN(x) 38 39 /* 40 For real x, let {x}_P or x_P be the closest representible number in the 41 floating point representation which uses P binary bits of fractional 42 precision (with IEEE rounding semantics). 43 44 Similarly, let f_P(x) be shorthand for {f(x)}_P. 45 46 Let ulp_P(x) be the unit of least precision for x: in other words the 47 maximal value of |a_P - b_P| where a_P <= x <= b_P and a_P != b_P. 48 49 Let x ~= y <-> | x - y | < ulp_P(x - y). 50 51 Let deg(x) be the value of x radians in degrees. 52 53 Values for each precision P were selected as follows. 54 55 56 COSD_SMALL = 2**{-N} such that for all x <= COSD_SMALL: 57 58 * cos(deg(x)) ~= 1, or equivalently: 59 60 | 1 - cos(deg(x)) | < ulp_P(1). 61 62 Unfortunately for SIND (and therefore TAND) a similar relation is only 63 possible for REAL(4) and REAL(8). With REAL(10) and REAL(16), enough 64 precision is available such that sin_P(x) != x_P for some x less than any 65 value. (There are values where this equality holds, but the distance has 66 inflection points.) 67 68 For REAL(4) and REAL(8), we can select SIND_SMALL such that: 69 70 * sin(deg(x)) ~= deg(x), or equivalently: 71 72 | deg(x) - sin(deg(x)) | < ulp_P(deg(x)). 73 74 */ 75 76 #ifdef HAVE_GFC_REAL_4 77 78 /* Build _gfortran_sind_r4, _gfortran_cosd_r4, and _gfortran_tand_r4 */ 79 80 #define KIND 4 81 #define TINY 0x1.p-100 /* ~= 7.889e-31 */ 82 #define COSD_SMALL 0x1.p-7 /* = 7.8125e-3 */ 83 #define SIND_SMALL 0x1.p-5 /* = 3.125e-2 */ 84 #define COSD30 8.66025388e-01 85 #define PIO180H 1.74560547e-02 /* high 12 bits. */ 86 #define PIO180L -2.76216747e-06 /* Next 24 bits. */ 87 88 #if defined(HAVE_FABSF) && defined(HAVE_FMODF) && defined(HAVE_COPYSIGNF) 89 90 #ifdef HAVE_SINF 91 #define ENABLE_SIND 92 #endif 93 94 #ifdef HAVE_COSF 95 #define ENABLE_COSD 96 #endif 97 98 #ifdef HAVE_TANF 99 #define ENABLE_TAND 100 #endif 101 102 #endif /* HAVE_FABSF && HAVE_FMODF && HAVE_COPYSIGNF */ 103 104 #ifdef GFC_REAL_4_INFINITY 105 #define HAVE_INFINITY_KIND 106 #endif 107 108 #include "trigd_lib.inc" 109 110 #undef KIND 111 #undef TINY 112 #undef COSD_SMALL 113 #undef SIND_SMALL 114 #undef COSD30 115 #undef PIO180H 116 #undef PIO180L 117 #undef ENABLE_SIND 118 #undef ENABLE_COSD 119 #undef ENABLE_TAND 120 #undef HAVE_INFINITY_KIND 121 122 #endif /* HAVE_GFC_REAL_4... */ 123 124 125 #ifdef HAVE_GFC_REAL_8 126 127 /* Build _gfortran_sind_r8, _gfortran_cosd_r8, and _gfortran_tand_r8 */ 128 129 #define KIND 8 130 #define TINY 0x1.p-1000 /* ~= 9.33e-302 (min exp -1074) */ 131 #define COSD_SMALL 0x1.p-21 /* ~= 4.768e-7 */ 132 #define SIND_SMALL 0x1.p-19 /* ~= 9.537e-7 */ 133 #define COSD30 8.6602540378443860e-01 134 #define PIO180H 1.7453283071517944e-02 /* high 21 bits. */ 135 #define PIO180L 9.4484253514332993e-09 /* Next 53 bits. */ 136 137 #if defined(HAVE_FABS) && defined(HAVE_FMOD) && defined(HAVE_COPYSIGN) 138 139 #ifdef HAVE_SIN 140 #define ENABLE_SIND 141 #endif 142 143 #ifdef HAVE_COS 144 #define ENABLE_COSD 145 #endif 146 147 #ifdef HAVE_TAN 148 #define ENABLE_TAND 149 #endif 150 151 #endif /* HAVE_FABS && HAVE_FMOD && HAVE_COPYSIGN */ 152 153 #ifdef GFC_REAL_8_INFINITY 154 #define HAVE_INFINITY_KIND 155 #endif 156 157 #include "trigd_lib.inc" 158 159 #undef KIND 160 #undef TINY 161 #undef COSD_SMALL 162 #undef SIND_SMALL 163 #undef COSD30 164 #undef PIO180H 165 #undef PIO180L 166 #undef ENABLE_SIND 167 #undef ENABLE_COSD 168 #undef ENABLE_TAND 169 #undef HAVE_INFINITY_KIND 170 171 #endif /* HAVE_GFC_REAL_8... */ 172 173 174 #ifdef HAVE_GFC_REAL_10 175 176 /* Build _gfortran_sind_r10, _gfortran_cosd_r10, and _gfortran_tand_r10 */ 177 178 #define KIND 10 179 #define TINY 0x1.p-16400 /* ~= 1.28e-4937 (min exp -16494) */ 180 #define COSD_SMALL 0x1.p-26 /* ~= 1.490e-8 */ 181 #undef SIND_SMALL /* not precise */ 182 #define COSD30 8.66025403784438646787e-01 183 #define PIO180H 1.74532925229868851602e-02 /* high 32 bits */ 184 #define PIO180L -3.04358939097084072823e-12 /* Next 64 bits */ 185 186 #if defined(HAVE_FABSL) && defined(HAVE_FMODL) && defined(HAVE_COPYSIGNL) 187 188 #ifdef HAVE_SINL 189 #define ENABLE_SIND 190 #endif 191 192 #ifdef HAVE_COSL 193 #define ENABLE_COSD 194 #endif 195 196 #ifdef HAVE_TANL 197 #define ENABLE_TAND 198 #endif 199 200 #endif /* HAVE_FABSL && HAVE_FMODL && HAVE_COPYSIGNL */ 201 202 #ifdef GFC_REAL_10_INFINITY 203 #define HAVE_INFINITY_KIND 204 #endif 205 206 #include "trigd_lib.inc" 207 208 #undef KIND 209 #undef TINY 210 #undef COSD_SMALL 211 #undef SIND_SMALL 212 #undef COSD30 213 #undef PIO180H 214 #undef PIO180L 215 #undef ENABLE_SIND 216 #undef ENABLE_COSD 217 #undef ENABLE_TAND 218 #undef HAVE_INFINITY_KIND 219 220 #endif /* HAVE_GFC_REAL_10 */ 221 222 223 #ifdef HAVE_GFC_REAL_16 224 225 /* Build _gfortran_sind_r16, _gfortran_cosd_r16, and _gfortran_tand_r16 */ 226 227 #define KIND 16 228 #define TINY 0x1.p-16400 /* ~= 1.28e-4937 */ 229 #undef SIND_SMALL /* not precise */ 230 231 #if GFC_REAL_16_DIGITS == 64 232 /* 80 bit precision, use constants from REAL(10). */ 233 #define COSD_SMALL 0x1.p-26 /* ~= 1.490e-8 */ 234 #define COSD30 8.66025403784438646787e-01 235 #define PIO180H 1.74532925229868851602e-02 /* high 32 bits */ 236 #define PIO180L -3.04358939097084072823e-12 /* Next 64 bits */ 237 238 #else 239 /* Proper float128 precision. */ 240 #define COSD_SMALL 0x1.p-51 /* ~= 4.441e-16 */ 241 #define COSD30 8.66025403784438646763723170752936183e-01 242 #define PIO180H 1.74532925199433197605003442731685936e-02 243 #define PIO180L -2.39912634365882824665106671063098954e-17 244 #endif 245 246 #ifdef GFC_REAL_16_IS_LONG_DOUBLE 247 248 #if defined(HAVE_FABSL) && defined(HAVE_FMODL) && defined(HAVE_COPYSIGNL) 249 250 #ifdef HAVE_SINL 251 #define ENABLE_SIND 252 #endif 253 254 #ifdef HAVE_COSL 255 #define ENABLE_COSD 256 #endif 257 258 #ifdef HAVE_TANL 259 #define ENABLE_TAND 260 #endif 261 262 #endif /* HAVE_FABSL && HAVE_FMODL && HAVE_COPYSIGNL */ 263 264 #else 265 266 /* libquadmath: HAVE_*Q are never defined. They must be available. */ 267 #define ENABLE_SIND 268 #define ENABLE_COSD 269 #define ENABLE_TAND 270 271 #endif /* GFC_REAL_16_IS_LONG_DOUBLE */ 272 273 #ifdef GFC_REAL_16_INFINITY 274 #define HAVE_INFINITY_KIND 275 #endif 276 277 #include "trigd_lib.inc" 278 279 #undef KIND 280 #undef TINY 281 #undef COSD_SMALL 282 #undef SIND_SMALL 283 #undef COSD30 284 #undef PIO180H 285 #undef PIO180L 286 #undef ENABLE_SIND 287 #undef ENABLE_COSD 288 #undef ENABLE_TAND 289 #undef HAVE_INFINITY_KIND 290 291 #endif /* HAVE_GFC_REAL_16 */ 292