libgfortran/intrinsics/trigd.c

   1 /* Implementation of the degree trignometric functions COSD, SIND, TAND.
   2    Copyright (C) 2020-2023 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 #elif defined(GFC_REAL_16_USE_IEC_60559)
 265
 266 #if defined(HAVE_FABSF128) && defined(HAVE_FMODF128) && defined(HAVE_COPYSIGNF128)
 267
 268 #ifdef HAVE_SINF128
 269 #define ENABLE_SIND
 270 #endif
 271
 272 #ifdef HAVE_COSF128
 273 #define ENABLE_COSD
 274 #endif
 275
 276 #ifdef HAVE_TANF128
 277 #define ENABLE_TAND
 278 #endif
 279
 280 #endif /* HAVE_FABSF128 && HAVE_FMODF128 && HAVE_COPYSIGNF128 */
 281
 282 #else
 283
 284 /* libquadmath: HAVE_*Q are never defined.  They must be available.  */
 285 #define ENABLE_SIND
 286 #define ENABLE_COSD
 287 #define ENABLE_TAND
 288
 289 #endif /* GFC_REAL_16_IS_LONG_DOUBLE */
 290
 291 #ifdef GFC_REAL_16_INFINITY
 292 #define HAVE_INFINITY_KIND
 293 #endif
 294
 295 #include "trigd_lib.inc"
 296
 297 #undef KIND
 298 #undef TINY
 299 #undef COSD_SMALL
 300 #undef SIND_SMALL
 301 #undef COSD30
 302 #undef PIO180H
 303 #undef PIO180L
 304 #undef ENABLE_SIND
 305 #undef ENABLE_COSD
 306 #undef ENABLE_TAND
 307 #undef HAVE_INFINITY_KIND
 308
 309 #endif /* HAVE_GFC_REAL_16 */
 310
 311 #ifdef HAVE_GFC_REAL_17
 312
 313 /* Build _gfortran_sind_r17, _gfortran_cosd_r17, and _gfortran_tand_r17  */
 314
 315 #define KIND    17
 316 #define TINY    0x1.p-16400     /* ~= 1.28e-4937 */
 317 #undef  SIND_SMALL              /* not precise */
 318
 319 /* Proper float128 precision.  */
 320 #define COSD_SMALL  0x1.p-51    /* ~= 4.441e-16 */
 321 #define COSD30      8.66025403784438646763723170752936183e-01
 322 #define PIO180H     1.74532925199433197605003442731685936e-02
 323 #define PIO180L     -2.39912634365882824665106671063098954e-17
 324
 325 /* libquadmath or glibc 2.32+: HAVE_*Q are never defined.  They must be available.  */
 326 #define ENABLE_SIND
 327 #define ENABLE_COSD
 328 #define ENABLE_TAND
 329
 330 #ifdef GFC_REAL_17_INFINITY
 331 #define HAVE_INFINITY_KIND
 332 #endif
 333
 334 #ifdef POWER_IEEE128
 335 #define COPYSIGN __copysignieee128
 336 #define FMOD __fmodieee128
 337 #define FABS __fabsieee128
 338 #define FMA __fmaieee128
 339 #define SIN __sinieee128
 340 #define COS __cosieee128
 341 #define TAN __tanieee128
 342 #endif
 343
 344 #include "trigd_lib.inc"
 345
 346 #undef KIND
 347 #undef TINY
 348 #undef COSD_SMALL
 349 #undef SIND_SMALL
 350 #undef COSD30
 351 #undef PIO180H
 352 #undef PIO180L
 353 #undef ENABLE_SIND
 354 #undef ENABLE_COSD
 355 #undef ENABLE_TAND
 356 #undef HAVE_INFINITY_KIND
 357
 358 #endif /* HAVE_GFC_REAL_17 */