]>
Commit | Line | Data |
---|---|---|
a2a2059f BD |
1 | /* Implementation of various C99 functions |
2 | Copyright (C) 2004 Free Software Foundation, Inc. | |
3 | ||
4 | This file is part of the GNU Fortran 95 runtime library (libgfortran). | |
5 | ||
6 | Libgfortran is free software; you can redistribute it and/or | |
57dea9f6 | 7 | modify it under the terms of the GNU General Public |
a2a2059f | 8 | License as published by the Free Software Foundation; either |
57dea9f6 TM |
9 | version 2 of the License, or (at your option) any later version. |
10 | ||
11 | In addition to the permissions in the GNU General Public License, the | |
12 | Free Software Foundation gives you unlimited permission to link the | |
13 | compiled version of this file into combinations with other programs, | |
14 | and to distribute those combinations without any restriction coming | |
15 | from the use of this file. (The General Public License restrictions | |
16 | do apply in other respects; for example, they cover modification of | |
17 | the file, and distribution when not linked into a combine | |
18 | executable.) | |
a2a2059f BD |
19 | |
20 | Libgfortran is distributed in the hope that it will be useful, | |
21 | but WITHOUT ANY WARRANTY; without even the implied warranty of | |
22 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the | |
57dea9f6 | 23 | GNU General Public License for more details. |
a2a2059f | 24 | |
57dea9f6 TM |
25 | You should have received a copy of the GNU General Public |
26 | License along with libgfortran; see the file COPYING. If not, | |
fe2ae685 KC |
27 | write to the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor, |
28 | Boston, MA 02110-1301, USA. */ | |
a2a2059f BD |
29 | |
30 | #include "config.h" | |
31 | #include <sys/types.h> | |
453310d8 | 32 | #include <float.h> |
a2a2059f | 33 | #include <math.h> |
1409cd0b FXC |
34 | |
35 | #define C99_PROTOS_H WE_DONT_WANT_PROTOS_NOW | |
a2a2059f BD |
36 | #include "libgfortran.h" |
37 | ||
38 | ||
453310d8 | 39 | #ifndef HAVE_ACOSF |
1409cd0b | 40 | #define HAVE_ACOSF |
453310d8 RS |
41 | float |
42 | acosf(float x) | |
43 | { | |
44 | return (float) acos(x); | |
45 | } | |
46 | #endif | |
47 | ||
48 | #ifndef HAVE_ASINF | |
1409cd0b | 49 | #define HAVE_ASINF |
453310d8 RS |
50 | float |
51 | asinf(float x) | |
52 | { | |
53 | return (float) asin(x); | |
54 | } | |
55 | #endif | |
56 | ||
57 | #ifndef HAVE_ATAN2F | |
1409cd0b | 58 | #define HAVE_ATAN2F |
453310d8 RS |
59 | float |
60 | atan2f(float y, float x) | |
61 | { | |
62 | return (float) atan2(y, x); | |
63 | } | |
64 | #endif | |
65 | ||
66 | #ifndef HAVE_ATANF | |
1409cd0b | 67 | #define HAVE_ATANF |
453310d8 RS |
68 | float |
69 | atanf(float x) | |
70 | { | |
71 | return (float) atan(x); | |
72 | } | |
73 | #endif | |
74 | ||
75 | #ifndef HAVE_CEILF | |
1409cd0b | 76 | #define HAVE_CEILF |
453310d8 RS |
77 | float |
78 | ceilf(float x) | |
79 | { | |
80 | return (float) ceil(x); | |
81 | } | |
82 | #endif | |
83 | ||
84 | #ifndef HAVE_COPYSIGNF | |
1409cd0b | 85 | #define HAVE_COPYSIGNF |
453310d8 RS |
86 | float |
87 | copysignf(float x, float y) | |
88 | { | |
89 | return (float) copysign(x, y); | |
90 | } | |
91 | #endif | |
92 | ||
93 | #ifndef HAVE_COSF | |
1409cd0b | 94 | #define HAVE_COSF |
453310d8 RS |
95 | float |
96 | cosf(float x) | |
97 | { | |
98 | return (float) cos(x); | |
99 | } | |
100 | #endif | |
101 | ||
102 | #ifndef HAVE_COSHF | |
1409cd0b | 103 | #define HAVE_COSHF |
453310d8 RS |
104 | float |
105 | coshf(float x) | |
106 | { | |
107 | return (float) cosh(x); | |
108 | } | |
109 | #endif | |
110 | ||
111 | #ifndef HAVE_EXPF | |
1409cd0b | 112 | #define HAVE_EXPF |
453310d8 RS |
113 | float |
114 | expf(float x) | |
115 | { | |
116 | return (float) exp(x); | |
117 | } | |
118 | #endif | |
119 | ||
6e4d9244 | 120 | #ifndef HAVE_FABSF |
1409cd0b | 121 | #define HAVE_FABSF |
6e4d9244 EB |
122 | float |
123 | fabsf(float x) | |
124 | { | |
125 | return (float) fabs(x); | |
126 | } | |
127 | #endif | |
128 | ||
453310d8 | 129 | #ifndef HAVE_FLOORF |
1409cd0b | 130 | #define HAVE_FLOORF |
453310d8 RS |
131 | float |
132 | floorf(float x) | |
133 | { | |
134 | return (float) floor(x); | |
135 | } | |
136 | #endif | |
137 | ||
138 | #ifndef HAVE_FREXPF | |
1409cd0b | 139 | #define HAVE_FREXPF |
453310d8 RS |
140 | float |
141 | frexpf(float x, int *exp) | |
142 | { | |
143 | return (float) frexp(x, exp); | |
144 | } | |
145 | #endif | |
146 | ||
147 | #ifndef HAVE_HYPOTF | |
1409cd0b | 148 | #define HAVE_HYPOTF |
453310d8 RS |
149 | float |
150 | hypotf(float x, float y) | |
151 | { | |
152 | return (float) hypot(x, y); | |
153 | } | |
154 | #endif | |
155 | ||
156 | #ifndef HAVE_LOGF | |
1409cd0b | 157 | #define HAVE_LOGF |
453310d8 RS |
158 | float |
159 | logf(float x) | |
160 | { | |
161 | return (float) log(x); | |
162 | } | |
163 | #endif | |
164 | ||
165 | #ifndef HAVE_LOG10F | |
1409cd0b | 166 | #define HAVE_LOG10F |
453310d8 RS |
167 | float |
168 | log10f(float x) | |
169 | { | |
170 | return (float) log10(x); | |
171 | } | |
172 | #endif | |
173 | ||
ae973d6a | 174 | #ifndef HAVE_SCALBN |
1409cd0b | 175 | #define HAVE_SCALBN |
ae973d6a FXC |
176 | double |
177 | scalbn(double x, int y) | |
178 | { | |
179 | return x * pow(FLT_RADIX, y); | |
180 | } | |
181 | #endif | |
182 | ||
453310d8 | 183 | #ifndef HAVE_SCALBNF |
1409cd0b | 184 | #define HAVE_SCALBNF |
453310d8 RS |
185 | float |
186 | scalbnf(float x, int y) | |
187 | { | |
188 | return (float) scalbn(x, y); | |
189 | } | |
190 | #endif | |
191 | ||
192 | #ifndef HAVE_SINF | |
1409cd0b | 193 | #define HAVE_SINF |
453310d8 RS |
194 | float |
195 | sinf(float x) | |
196 | { | |
197 | return (float) sin(x); | |
198 | } | |
199 | #endif | |
200 | ||
201 | #ifndef HAVE_SINHF | |
1409cd0b | 202 | #define HAVE_SINHF |
453310d8 RS |
203 | float |
204 | sinhf(float x) | |
205 | { | |
206 | return (float) sinh(x); | |
207 | } | |
208 | #endif | |
209 | ||
210 | #ifndef HAVE_SQRTF | |
1409cd0b | 211 | #define HAVE_SQRTF |
453310d8 RS |
212 | float |
213 | sqrtf(float x) | |
214 | { | |
215 | return (float) sqrt(x); | |
216 | } | |
217 | #endif | |
218 | ||
219 | #ifndef HAVE_TANF | |
1409cd0b | 220 | #define HAVE_TANF |
453310d8 RS |
221 | float |
222 | tanf(float x) | |
223 | { | |
224 | return (float) tan(x); | |
225 | } | |
226 | #endif | |
227 | ||
228 | #ifndef HAVE_TANHF | |
1409cd0b | 229 | #define HAVE_TANHF |
453310d8 RS |
230 | float |
231 | tanhf(float x) | |
232 | { | |
233 | return (float) tanh(x); | |
234 | } | |
235 | #endif | |
236 | ||
69a2d125 | 237 | #ifndef HAVE_TRUNC |
1409cd0b | 238 | #define HAVE_TRUNC |
69a2d125 EB |
239 | double |
240 | trunc(double x) | |
241 | { | |
242 | if (!isfinite (x)) | |
243 | return x; | |
244 | ||
245 | if (x < 0.0) | |
246 | return - floor (-x); | |
247 | else | |
248 | return floor (x); | |
249 | } | |
250 | #endif | |
251 | ||
252 | #ifndef HAVE_TRUNCF | |
1409cd0b | 253 | #define HAVE_TRUNCF |
69a2d125 EB |
254 | float |
255 | truncf(float x) | |
256 | { | |
257 | return (float) trunc (x); | |
258 | } | |
259 | #endif | |
260 | ||
453310d8 | 261 | #ifndef HAVE_NEXTAFTERF |
1409cd0b | 262 | #define HAVE_NEXTAFTERF |
453310d8 RS |
263 | /* This is a portable implementation of nextafterf that is intended to be |
264 | independent of the floating point format or its in memory representation. | |
067a5735 | 265 | This implementation works correctly with denormalized values. */ |
453310d8 RS |
266 | float |
267 | nextafterf(float x, float y) | |
268 | { | |
067a5735 RS |
269 | /* This variable is marked volatile to avoid excess precision problems |
270 | on some platforms, including IA-32. */ | |
271 | volatile float delta; | |
272 | float absx, denorm_min; | |
453310d8 RS |
273 | |
274 | if (isnan(x) || isnan(y)) | |
067a5735 | 275 | return x + y; |
453310d8 RS |
276 | if (x == y) |
277 | return x; | |
74421469 EB |
278 | if (!isfinite (x)) |
279 | return x > 0 ? __FLT_MAX__ : - __FLT_MAX__; | |
453310d8 | 280 | |
067a5735 RS |
281 | /* absx = fabsf (x); */ |
282 | absx = (x < 0.0) ? -x : x; | |
453310d8 | 283 | |
067a5735 RS |
284 | /* __FLT_DENORM_MIN__ is non-zero iff the target supports denormals. */ |
285 | if (__FLT_DENORM_MIN__ == 0.0f) | |
286 | denorm_min = __FLT_MIN__; | |
287 | else | |
288 | denorm_min = __FLT_DENORM_MIN__; | |
289 | ||
290 | if (absx < __FLT_MIN__) | |
291 | delta = denorm_min; | |
453310d8 RS |
292 | else |
293 | { | |
067a5735 RS |
294 | float frac; |
295 | int exp; | |
296 | ||
297 | /* Discard the fraction from x. */ | |
298 | frac = frexpf (absx, &exp); | |
299 | delta = scalbnf (0.5f, exp); | |
300 | ||
301 | /* Scale x by the epsilon of the representation. By rights we should | |
302 | have been able to combine this with scalbnf, but some targets don't | |
303 | get that correct with denormals. */ | |
304 | delta *= __FLT_EPSILON__; | |
305 | ||
306 | /* If we're going to be reducing the absolute value of X, and doing so | |
307 | would reduce the exponent of X, then the delta to be applied is | |
308 | one exponent smaller. */ | |
309 | if (frac == 0.5f && (y < x) == (x > 0)) | |
310 | delta *= 0.5f; | |
311 | ||
312 | /* If that underflows to zero, then we're back to the minimum. */ | |
313 | if (delta == 0.0f) | |
314 | delta = denorm_min; | |
453310d8 | 315 | } |
067a5735 RS |
316 | |
317 | if (y < x) | |
318 | delta = -delta; | |
319 | ||
320 | return x + delta; | |
453310d8 RS |
321 | } |
322 | #endif | |
323 | ||
bf4d99cf TS |
324 | |
325 | #ifndef HAVE_POWF | |
1409cd0b | 326 | #define HAVE_POWF |
bf4d99cf TS |
327 | float |
328 | powf(float x, float y) | |
329 | { | |
330 | return (float) pow(x, y); | |
331 | } | |
332 | #endif | |
333 | ||
69d3c9a4 | 334 | /* Note that if fpclassify is not defined, then NaN is not handled */ |
bc20e36d | 335 | |
a2a2059f BD |
336 | /* Algorithm by Steven G. Kargl. */ |
337 | ||
338 | #ifndef HAVE_ROUND | |
1409cd0b | 339 | #define HAVE_ROUND |
a2a2059f BD |
340 | /* Round to nearest integral value. If the argument is halfway between two |
341 | integral values then round away from zero. */ | |
342 | ||
343 | double | |
344 | round(double x) | |
345 | { | |
346 | double t; | |
118ea208 | 347 | if (!isfinite (x)) |
a2a2059f BD |
348 | return (x); |
349 | ||
350 | if (x >= 0.0) | |
351 | { | |
352 | t = ceil(x); | |
353 | if (t - x > 0.5) | |
354 | t -= 1.0; | |
355 | return (t); | |
356 | } | |
357 | else | |
358 | { | |
359 | t = ceil(-x); | |
360 | if (t + x > 0.5) | |
361 | t -= 1.0; | |
362 | return (-t); | |
363 | } | |
364 | } | |
365 | #endif | |
366 | ||
367 | #ifndef HAVE_ROUNDF | |
1409cd0b | 368 | #define HAVE_ROUNDF |
a2a2059f BD |
369 | /* Round to nearest integral value. If the argument is halfway between two |
370 | integral values then round away from zero. */ | |
371 | ||
372 | float | |
373 | roundf(float x) | |
374 | { | |
375 | float t; | |
118ea208 | 376 | if (!isfinite (x)) |
a2a2059f BD |
377 | return (x); |
378 | ||
379 | if (x >= 0.0) | |
380 | { | |
381 | t = ceilf(x); | |
382 | if (t - x > 0.5) | |
383 | t -= 1.0; | |
384 | return (t); | |
385 | } | |
386 | else | |
387 | { | |
388 | t = ceilf(-x); | |
389 | if (t + x > 0.5) | |
390 | t -= 1.0; | |
391 | return (-t); | |
392 | } | |
393 | } | |
394 | #endif | |
32aa3bff FXC |
395 | |
396 | #ifndef HAVE_LOG10L | |
1409cd0b | 397 | #define HAVE_LOG10L |
32aa3bff FXC |
398 | /* log10 function for long double variables. The version provided here |
399 | reduces the argument until it fits into a double, then use log10. */ | |
400 | long double | |
401 | log10l(long double x) | |
402 | { | |
403 | #if LDBL_MAX_EXP > DBL_MAX_EXP | |
404 | if (x > DBL_MAX) | |
405 | { | |
406 | double val; | |
407 | int p2_result = 0; | |
408 | if (x > 0x1p16383L) { p2_result += 16383; x /= 0x1p16383L; } | |
409 | if (x > 0x1p8191L) { p2_result += 8191; x /= 0x1p8191L; } | |
410 | if (x > 0x1p4095L) { p2_result += 4095; x /= 0x1p4095L; } | |
411 | if (x > 0x1p2047L) { p2_result += 2047; x /= 0x1p2047L; } | |
412 | if (x > 0x1p1023L) { p2_result += 1023; x /= 0x1p1023L; } | |
413 | val = log10 ((double) x); | |
414 | return (val + p2_result * .30102999566398119521373889472449302L); | |
415 | } | |
416 | #endif | |
417 | #if LDBL_MIN_EXP < DBL_MIN_EXP | |
418 | if (x < DBL_MIN) | |
419 | { | |
420 | double val; | |
421 | int p2_result = 0; | |
422 | if (x < 0x1p-16380L) { p2_result += 16380; x /= 0x1p-16380L; } | |
423 | if (x < 0x1p-8189L) { p2_result += 8189; x /= 0x1p-8189L; } | |
424 | if (x < 0x1p-4093L) { p2_result += 4093; x /= 0x1p-4093L; } | |
425 | if (x < 0x1p-2045L) { p2_result += 2045; x /= 0x1p-2045L; } | |
426 | if (x < 0x1p-1021L) { p2_result += 1021; x /= 0x1p-1021L; } | |
427 | val = fabs(log10 ((double) x)); | |
428 | return (- val - p2_result * .30102999566398119521373889472449302L); | |
429 | } | |
430 | #endif | |
431 | return log10 (x); | |
432 | } | |
433 | #endif | |
1409cd0b FXC |
434 | |
435 | ||
436 | #if !defined(HAVE_CABSF) | |
437 | #define HAVE_CABSF | |
438 | float | |
439 | cabsf (float complex z) | |
440 | { | |
441 | return hypotf (REALPART (z), IMAGPART (z)); | |
442 | } | |
443 | #endif | |
444 | ||
445 | #if !defined(HAVE_CABS) | |
446 | #define HAVE_CABS | |
447 | double | |
448 | cabs (double complex z) | |
449 | { | |
450 | return hypot (REALPART (z), IMAGPART (z)); | |
451 | } | |
452 | #endif | |
453 | ||
454 | #if !defined(HAVE_CABSL) && defined(HAVE_HYPOTL) | |
455 | #define HAVE_CABSL | |
456 | long double | |
457 | cabsl (long double complex z) | |
458 | { | |
459 | return hypotl (REALPART (z), IMAGPART (z)); | |
460 | } | |
461 | #endif | |
462 | ||
463 | ||
464 | #if !defined(HAVE_CARGF) | |
465 | #define HAVE_CARGF | |
466 | float | |
467 | cargf (float complex z) | |
468 | { | |
469 | return atan2f (IMAGPART (z), REALPART (z)); | |
470 | } | |
471 | #endif | |
472 | ||
473 | #if !defined(HAVE_CARG) | |
474 | #define HAVE_CARG | |
475 | double | |
476 | carg (double complex z) | |
477 | { | |
478 | return atan2 (IMAGPART (z), REALPART (z)); | |
479 | } | |
480 | #endif | |
481 | ||
482 | #if !defined(HAVE_CARGL) && defined(HAVE_ATAN2L) | |
483 | #define HAVE_CARGL | |
484 | long double | |
485 | cargl (long double complex z) | |
486 | { | |
487 | return atan2l (IMAGPART (z), REALPART (z)); | |
488 | } | |
489 | #endif | |
490 | ||
491 | ||
492 | /* exp(z) = exp(a)*(cos(b) + i sin(b)) */ | |
493 | #if !defined(HAVE_CEXPF) | |
494 | #define HAVE_CEXPF | |
495 | float complex | |
496 | cexpf (float complex z) | |
497 | { | |
498 | float a, b; | |
499 | float complex v; | |
500 | ||
501 | a = REALPART (z); | |
502 | b = IMAGPART (z); | |
503 | COMPLEX_ASSIGN (v, cosf (b), sinf (b)); | |
504 | return expf (a) * v; | |
505 | } | |
506 | #endif | |
507 | ||
508 | #if !defined(HAVE_CEXP) | |
509 | #define HAVE_CEXP | |
510 | double complex | |
511 | cexp (double complex z) | |
512 | { | |
513 | double a, b; | |
514 | double complex v; | |
515 | ||
516 | a = REALPART (z); | |
517 | b = IMAGPART (z); | |
518 | COMPLEX_ASSIGN (v, cos (b), sin (b)); | |
519 | return exp (a) * v; | |
520 | } | |
521 | #endif | |
522 | ||
523 | #if !defined(HAVE_CEXPL) && defined(HAVE_COSL) && defined(HAVE_SINL) && defined(EXPL) | |
524 | #define HAVE_CEXPL | |
525 | long double complex | |
526 | cexpl (long double complex z) | |
527 | { | |
528 | long double a, b; | |
529 | long double complex v; | |
530 | ||
531 | a = REALPART (z); | |
532 | b = IMAGPART (z); | |
533 | COMPLEX_ASSIGN (v, cosl (b), sinl (b)); | |
534 | return expl (a) * v; | |
535 | } | |
536 | #endif | |
537 | ||
538 | ||
539 | /* log(z) = log (cabs(z)) + i*carg(z) */ | |
540 | #if !defined(HAVE_CLOGF) | |
541 | #define HAVE_CLOGF | |
542 | float complex | |
543 | clogf (float complex z) | |
544 | { | |
545 | float complex v; | |
546 | ||
547 | COMPLEX_ASSIGN (v, logf (cabsf (z)), cargf (z)); | |
548 | return v; | |
549 | } | |
550 | #endif | |
551 | ||
552 | #if !defined(HAVE_CLOG) | |
553 | #define HAVE_CLOG | |
554 | double complex | |
555 | clog (double complex z) | |
556 | { | |
557 | double complex v; | |
558 | ||
559 | COMPLEX_ASSIGN (v, log (cabs (z)), carg (z)); | |
560 | return v; | |
561 | } | |
562 | #endif | |
563 | ||
564 | #if !defined(HAVE_CLOGL) && defined(HAVE_LOGL) && defined(HAVE_CABSL) && defined(HAVE_CARGL) | |
565 | #define HAVE_CLOGL | |
566 | long double complex | |
567 | clogl (long double complex z) | |
568 | { | |
569 | long double complex v; | |
570 | ||
571 | COMPLEX_ASSIGN (v, logl (cabsl (z)), cargl (z)); | |
572 | return v; | |
573 | } | |
574 | #endif | |
575 | ||
576 | ||
577 | /* log10(z) = log10 (cabs(z)) + i*carg(z) */ | |
578 | #if !defined(HAVE_CLOG10F) | |
579 | #define HAVE_CLOG10F | |
580 | float complex | |
581 | clog10f (float complex z) | |
582 | { | |
583 | float complex v; | |
584 | ||
585 | COMPLEX_ASSIGN (v, log10f (cabsf (z)), cargf (z)); | |
586 | return v; | |
587 | } | |
588 | #endif | |
589 | ||
590 | #if !defined(HAVE_CLOG10) | |
591 | #define HAVE_CLOG10 | |
592 | double complex | |
593 | clog10 (double complex z) | |
594 | { | |
595 | double complex v; | |
596 | ||
597 | COMPLEX_ASSIGN (v, log10 (cabs (z)), carg (z)); | |
598 | return v; | |
599 | } | |
600 | #endif | |
601 | ||
602 | #if !defined(HAVE_CLOG10L) && defined(HAVE_LOG10L) && defined(HAVE_CABSL) && defined(HAVE_CARGL) | |
603 | #define HAVE_CLOG10L | |
604 | long double complex | |
605 | clog10l (long double complex z) | |
606 | { | |
607 | long double complex v; | |
608 | ||
609 | COMPLEX_ASSIGN (v, log10l (cabsl (z)), cargl (z)); | |
610 | return v; | |
611 | } | |
612 | #endif | |
613 | ||
614 | ||
615 | /* pow(base, power) = cexp (power * clog (base)) */ | |
616 | #if !defined(HAVE_CPOWF) | |
617 | #define HAVE_CPOWF | |
618 | float complex | |
619 | cpowf (float complex base, float complex power) | |
620 | { | |
621 | return cexpf (power * clogf (base)); | |
622 | } | |
623 | #endif | |
624 | ||
625 | #if !defined(HAVE_CPOW) | |
626 | #define HAVE_CPOW | |
627 | double complex | |
628 | cpow (double complex base, double complex power) | |
629 | { | |
630 | return cexp (power * clog (base)); | |
631 | } | |
632 | #endif | |
633 | ||
634 | #if !defined(HAVE_CPOWL) && defined(HAVE_CEXPL) && defined(HAVE_CLOGL) | |
635 | #define HAVE_CPOWL | |
636 | long double complex | |
637 | cpowl (long double complex base, long double complex power) | |
638 | { | |
639 | return cexpl (power * clogl (base)); | |
640 | } | |
641 | #endif | |
642 | ||
643 | ||
644 | /* sqrt(z). Algorithm pulled from glibc. */ | |
645 | #if !defined(HAVE_CSQRTF) | |
646 | #define HAVE_CSQRTF | |
647 | float complex | |
648 | csqrtf (float complex z) | |
649 | { | |
650 | float re, im; | |
651 | float complex v; | |
652 | ||
653 | re = REALPART (z); | |
654 | im = IMAGPART (z); | |
655 | if (im == 0) | |
656 | { | |
657 | if (re < 0) | |
658 | { | |
659 | COMPLEX_ASSIGN (v, 0, copysignf (sqrtf (-re), im)); | |
660 | } | |
661 | else | |
662 | { | |
663 | COMPLEX_ASSIGN (v, fabsf (sqrtf (re)), copysignf (0, im)); | |
664 | } | |
665 | } | |
666 | else if (re == 0) | |
667 | { | |
668 | float r; | |
669 | ||
670 | r = sqrtf (0.5 * fabsf (im)); | |
671 | ||
a2694f68 | 672 | COMPLEX_ASSIGN (v, r, copysignf (r, im)); |
1409cd0b FXC |
673 | } |
674 | else | |
675 | { | |
676 | float d, r, s; | |
677 | ||
678 | d = hypotf (re, im); | |
679 | /* Use the identity 2 Re res Im res = Im x | |
680 | to avoid cancellation error in d +/- Re x. */ | |
681 | if (re > 0) | |
682 | { | |
683 | r = sqrtf (0.5 * d + 0.5 * re); | |
684 | s = (0.5 * im) / r; | |
685 | } | |
686 | else | |
687 | { | |
688 | s = sqrtf (0.5 * d - 0.5 * re); | |
689 | r = fabsf ((0.5 * im) / s); | |
690 | } | |
691 | ||
692 | COMPLEX_ASSIGN (v, r, copysignf (s, im)); | |
693 | } | |
694 | return v; | |
695 | } | |
696 | #endif | |
697 | ||
698 | #if !defined(HAVE_CSQRT) | |
699 | #define HAVE_CSQRT | |
700 | double complex | |
701 | csqrt (double complex z) | |
702 | { | |
703 | double re, im; | |
704 | double complex v; | |
705 | ||
706 | re = REALPART (z); | |
707 | im = IMAGPART (z); | |
708 | if (im == 0) | |
709 | { | |
710 | if (re < 0) | |
711 | { | |
712 | COMPLEX_ASSIGN (v, 0, copysign (sqrt (-re), im)); | |
713 | } | |
714 | else | |
715 | { | |
716 | COMPLEX_ASSIGN (v, fabs (sqrt (re)), copysign (0, im)); | |
717 | } | |
718 | } | |
719 | else if (re == 0) | |
720 | { | |
721 | double r; | |
722 | ||
723 | r = sqrt (0.5 * fabs (im)); | |
724 | ||
a2694f68 | 725 | COMPLEX_ASSIGN (v, r, copysign (r, im)); |
1409cd0b FXC |
726 | } |
727 | else | |
728 | { | |
729 | double d, r, s; | |
730 | ||
731 | d = hypot (re, im); | |
732 | /* Use the identity 2 Re res Im res = Im x | |
733 | to avoid cancellation error in d +/- Re x. */ | |
734 | if (re > 0) | |
735 | { | |
736 | r = sqrt (0.5 * d + 0.5 * re); | |
737 | s = (0.5 * im) / r; | |
738 | } | |
739 | else | |
740 | { | |
741 | s = sqrt (0.5 * d - 0.5 * re); | |
742 | r = fabs ((0.5 * im) / s); | |
743 | } | |
744 | ||
745 | COMPLEX_ASSIGN (v, r, copysign (s, im)); | |
746 | } | |
747 | return v; | |
748 | } | |
749 | #endif | |
750 | ||
751 | #if !defined(HAVE_CSQRTL) && defined(HAVE_COPYSIGNL) && defined(HAVE_SQRTL) && defined(HAVE_FABSL) && defined(HAVE_HYPOTL) | |
752 | #define HAVE_CSQRTL | |
753 | long double complex | |
754 | csqrtl (long double complex z) | |
755 | { | |
756 | long double re, im; | |
757 | long double complex v; | |
758 | ||
759 | re = REALPART (z); | |
760 | im = IMAGPART (z); | |
761 | if (im == 0) | |
762 | { | |
763 | if (re < 0) | |
764 | { | |
765 | COMPLEX_ASSIGN (v, 0, copysignl (sqrtl (-re), im)); | |
766 | } | |
767 | else | |
768 | { | |
769 | COMPLEX_ASSIGN (v, fabsl (sqrtl (re)), copysignl (0, im)); | |
770 | } | |
771 | } | |
772 | else if (re == 0) | |
773 | { | |
774 | long double r; | |
775 | ||
776 | r = sqrtl (0.5 * fabsl (im)); | |
777 | ||
778 | COMPLEX_ASSIGN (v, copysignl (r, im), r); | |
779 | } | |
780 | else | |
781 | { | |
782 | long double d, r, s; | |
783 | ||
784 | d = hypotl (re, im); | |
785 | /* Use the identity 2 Re res Im res = Im x | |
786 | to avoid cancellation error in d +/- Re x. */ | |
787 | if (re > 0) | |
788 | { | |
789 | r = sqrtl (0.5 * d + 0.5 * re); | |
790 | s = (0.5 * im) / r; | |
791 | } | |
792 | else | |
793 | { | |
794 | s = sqrtl (0.5 * d - 0.5 * re); | |
795 | r = fabsl ((0.5 * im) / s); | |
796 | } | |
797 | ||
798 | COMPLEX_ASSIGN (v, r, copysignl (s, im)); | |
799 | } | |
800 | return v; | |
801 | } | |
802 | #endif | |
803 | ||
804 | ||
805 | /* sinh(a + i b) = sinh(a) cos(b) + i cosh(a) sin(b) */ | |
806 | #if !defined(HAVE_CSINHF) | |
807 | #define HAVE_CSINHF | |
808 | float complex | |
809 | csinhf (float complex a) | |
810 | { | |
811 | float r, i; | |
812 | float complex v; | |
813 | ||
814 | r = REALPART (a); | |
815 | i = IMAGPART (a); | |
816 | COMPLEX_ASSIGN (v, sinhf (r) * cosf (i), coshf (r) * sinf (i)); | |
817 | return v; | |
818 | } | |
819 | #endif | |
820 | ||
821 | #if !defined(HAVE_CSINH) | |
822 | #define HAVE_CSINH | |
823 | double complex | |
824 | csinh (double complex a) | |
825 | { | |
826 | double r, i; | |
827 | double complex v; | |
828 | ||
829 | r = REALPART (a); | |
830 | i = IMAGPART (a); | |
831 | COMPLEX_ASSIGN (v, sinh (r) * cos (i), cosh (r) * sin (i)); | |
832 | return v; | |
833 | } | |
834 | #endif | |
835 | ||
836 | #if !defined(HAVE_CSINHL) && defined(HAVE_COSL) && defined(HAVE_COSHL) && defined(HAVE_SINL) && defined(HAVE_SINHL) | |
837 | #define HAVE_CSINHL | |
838 | long double complex | |
839 | csinhl (long double complex a) | |
840 | { | |
841 | long double r, i; | |
842 | long double complex v; | |
843 | ||
844 | r = REALPART (a); | |
845 | i = IMAGPART (a); | |
846 | COMPLEX_ASSIGN (v, sinhl (r) * cosl (i), coshl (r) * sinl (i)); | |
847 | return v; | |
848 | } | |
849 | #endif | |
850 | ||
851 | ||
852 | /* cosh(a + i b) = cosh(a) cos(b) - i sinh(a) sin(b) */ | |
853 | #if !defined(HAVE_CCOSHF) | |
854 | #define HAVE_CCOSHF | |
855 | float complex | |
856 | ccoshf (float complex a) | |
857 | { | |
858 | float r, i; | |
859 | float complex v; | |
860 | ||
861 | r = REALPART (a); | |
862 | i = IMAGPART (a); | |
863 | COMPLEX_ASSIGN (v, coshf (r) * cosf (i), - (sinhf (r) * sinf (i))); | |
864 | return v; | |
865 | } | |
866 | #endif | |
867 | ||
868 | #if !defined(HAVE_CCOSH) | |
869 | #define HAVE_CCOSH | |
870 | double complex | |
871 | ccosh (double complex a) | |
872 | { | |
873 | double r, i; | |
874 | double complex v; | |
875 | ||
876 | r = REALPART (a); | |
877 | i = IMAGPART (a); | |
878 | COMPLEX_ASSIGN (v, cosh (r) * cos (i), - (sinh (r) * sin (i))); | |
879 | return v; | |
880 | } | |
881 | #endif | |
882 | ||
883 | #if !defined(HAVE_CCOSHL) && defined(HAVE_COSL) && defined(HAVE_COSHL) && defined(HAVE_SINL) && defined(HAVE_SINHL) | |
884 | #define HAVE_CCOSHL | |
885 | long double complex | |
886 | ccoshl (long double complex a) | |
887 | { | |
888 | long double r, i; | |
889 | long double complex v; | |
890 | ||
891 | r = REALPART (a); | |
892 | i = IMAGPART (a); | |
893 | COMPLEX_ASSIGN (v, coshl (r) * cosl (i), - (sinhl (r) * sinl (i))); | |
894 | return v; | |
895 | } | |
896 | #endif | |
897 | ||
898 | ||
899 | /* tanh(a + i b) = (tanh(a) + i tan(b)) / (1 - i tanh(a) tan(b)) */ | |
900 | #if !defined(HAVE_CTANHF) | |
901 | #define HAVE_CTANHF | |
902 | float complex | |
903 | ctanhf (float complex a) | |
904 | { | |
905 | float rt, it; | |
906 | float complex n, d; | |
907 | ||
908 | rt = tanhf (REALPART (a)); | |
909 | it = tanf (IMAGPART (a)); | |
910 | COMPLEX_ASSIGN (n, rt, it); | |
911 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
912 | ||
913 | return n / d; | |
914 | } | |
915 | #endif | |
916 | ||
917 | #if !defined(HAVE_CTANH) | |
918 | #define HAVE_CTANH | |
919 | double complex | |
920 | ctanh (double complex a) | |
921 | { | |
922 | double rt, it; | |
923 | double complex n, d; | |
924 | ||
925 | rt = tanh (REALPART (a)); | |
926 | it = tan (IMAGPART (a)); | |
927 | COMPLEX_ASSIGN (n, rt, it); | |
928 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
929 | ||
930 | return n / d; | |
931 | } | |
932 | #endif | |
933 | ||
934 | #if !defined(HAVE_CTANHL) && defined(HAVE_TANL) && defined(HAVE_TANHL) | |
935 | #define HAVE_CTANHL | |
936 | long double complex | |
937 | ctanhl (long double complex a) | |
938 | { | |
939 | long double rt, it; | |
940 | long double complex n, d; | |
941 | ||
942 | rt = tanhl (REALPART (a)); | |
943 | it = tanl (IMAGPART (a)); | |
944 | COMPLEX_ASSIGN (n, rt, it); | |
945 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
946 | ||
947 | return n / d; | |
948 | } | |
949 | #endif | |
950 | ||
951 | ||
952 | /* sin(a + i b) = sin(a) cosh(b) + i cos(a) sinh(b) */ | |
953 | #if !defined(HAVE_CSINF) | |
954 | #define HAVE_CSINF | |
955 | float complex | |
956 | csinf (float complex a) | |
957 | { | |
958 | float r, i; | |
959 | float complex v; | |
960 | ||
961 | r = REALPART (a); | |
962 | i = IMAGPART (a); | |
963 | COMPLEX_ASSIGN (v, sinf (r) * coshf (i), cosf (r) * sinhf (i)); | |
964 | return v; | |
965 | } | |
966 | #endif | |
967 | ||
968 | #if !defined(HAVE_CSIN) | |
969 | #define HAVE_CSIN | |
970 | double complex | |
971 | csin (double complex a) | |
972 | { | |
973 | double r, i; | |
974 | double complex v; | |
975 | ||
976 | r = REALPART (a); | |
977 | i = IMAGPART (a); | |
978 | COMPLEX_ASSIGN (v, sin (r) * cosh (i), cos (r) * sinh (i)); | |
979 | return v; | |
980 | } | |
981 | #endif | |
982 | ||
983 | #if !defined(HAVE_CSINL) && defined(HAVE_COSL) && defined(HAVE_COSHL) && defined(HAVE_SINL) && defined(HAVE_SINHL) | |
984 | #define HAVE_CSINL | |
985 | long double complex | |
986 | csinl (long double complex a) | |
987 | { | |
988 | long double r, i; | |
989 | long double complex v; | |
990 | ||
991 | r = REALPART (a); | |
992 | i = IMAGPART (a); | |
993 | COMPLEX_ASSIGN (v, sinl (r) * coshl (i), cosl (r) * sinhl (i)); | |
994 | return v; | |
995 | } | |
996 | #endif | |
997 | ||
998 | ||
999 | /* cos(a + i b) = cos(a) cosh(b) - i sin(a) sinh(b) */ | |
1000 | #if !defined(HAVE_CCOSF) | |
1001 | #define HAVE_CCOSF | |
1002 | float complex | |
1003 | ccosf (float complex a) | |
1004 | { | |
1005 | float r, i; | |
1006 | float complex v; | |
1007 | ||
1008 | r = REALPART (a); | |
1009 | i = IMAGPART (a); | |
1010 | COMPLEX_ASSIGN (v, cosf (r) * coshf (i), - (sinf (r) * sinhf (i))); | |
1011 | return v; | |
1012 | } | |
1013 | #endif | |
1014 | ||
1015 | #if !defined(HAVE_CCOS) | |
1016 | #define HAVE_CCOS | |
1017 | double complex | |
1018 | ccos (double complex a) | |
1019 | { | |
1020 | double r, i; | |
1021 | double complex v; | |
1022 | ||
1023 | r = REALPART (a); | |
1024 | i = IMAGPART (a); | |
1025 | COMPLEX_ASSIGN (v, cos (r) * cosh (i), - (sin (r) * sinh (i))); | |
1026 | return v; | |
1027 | } | |
1028 | #endif | |
1029 | ||
1030 | #if !defined(HAVE_CCOSL) && defined(HAVE_COSL) && defined(HAVE_COSHL) && defined(HAVE_SINL) && defined(HAVE_SINHL) | |
1031 | #define HAVE_CCOSL | |
1032 | long double complex | |
1033 | ccosl (long double complex a) | |
1034 | { | |
1035 | long double r, i; | |
1036 | long double complex v; | |
1037 | ||
1038 | r = REALPART (a); | |
1039 | i = IMAGPART (a); | |
1040 | COMPLEX_ASSIGN (v, cosl (r) * coshl (i), - (sinl (r) * sinhl (i))); | |
1041 | return v; | |
1042 | } | |
1043 | #endif | |
1044 | ||
1045 | ||
1046 | /* tan(a + i b) = (tan(a) + i tanh(b)) / (1 - i tan(a) tanh(b)) */ | |
1047 | #if !defined(HAVE_CTANF) | |
1048 | #define HAVE_CTANF | |
1049 | float complex | |
1050 | ctanf (float complex a) | |
1051 | { | |
1052 | float rt, it; | |
1053 | float complex n, d; | |
1054 | ||
1055 | rt = tanf (REALPART (a)); | |
1056 | it = tanhf (IMAGPART (a)); | |
1057 | COMPLEX_ASSIGN (n, rt, it); | |
1058 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
1059 | ||
1060 | return n / d; | |
1061 | } | |
1062 | #endif | |
1063 | ||
1064 | #if !defined(HAVE_CTAN) | |
1065 | #define HAVE_CTAN | |
1066 | double complex | |
1067 | ctan (double complex a) | |
1068 | { | |
1069 | double rt, it; | |
1070 | double complex n, d; | |
1071 | ||
1072 | rt = tan (REALPART (a)); | |
1073 | it = tanh (IMAGPART (a)); | |
1074 | COMPLEX_ASSIGN (n, rt, it); | |
1075 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
1076 | ||
1077 | return n / d; | |
1078 | } | |
1079 | #endif | |
1080 | ||
1081 | #if !defined(HAVE_CTANL) && defined(HAVE_TANL) && defined(HAVE_TANHL) | |
1082 | #define HAVE_CTANL | |
1083 | long double complex | |
1084 | ctanl (long double complex a) | |
1085 | { | |
1086 | long double rt, it; | |
1087 | long double complex n, d; | |
1088 | ||
1089 | rt = tanl (REALPART (a)); | |
1090 | it = tanhl (IMAGPART (a)); | |
1091 | COMPLEX_ASSIGN (n, rt, it); | |
1092 | COMPLEX_ASSIGN (d, 1, - (rt * it)); | |
1093 | ||
1094 | return n / d; | |
1095 | } | |
1096 | #endif | |
1097 |