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 */