OSDN Git Service

* configure.ac: Check for ieeefp.h. Check for fabsf in libm.
[pf3gnuchains/gcc-fork.git] / libgfortran / intrinsics / c99_functions.c
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
7 modify it under the terms of the GNU Lesser General Public
8 License as published by the Free Software Foundation; either
9 version 2.1 of the License, or (at your option) any later version.
10
11 Libgfortran is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14 GNU Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with libgfortran; see the file COPYING.LIB.  If not,
18 write to the Free Software Foundation, Inc., 59 Temple Place - Suite 330,
19 Boston, MA 02111-1307, USA.  */
20
21 #include "config.h"
22 #include <sys/types.h>
23 #include <float.h>
24 #include <math.h>
25 #include "libgfortran.h"
26
27
28 #ifndef HAVE_ACOSF
29 float
30 acosf(float x)
31 {
32   return (float) acos(x);
33 }
34 #endif
35
36 #ifndef HAVE_ASINF
37 float
38 asinf(float x)
39 {
40   return (float) asin(x);
41 }
42 #endif
43
44 #ifndef HAVE_ATAN2F
45 float
46 atan2f(float y, float x)
47 {
48   return (float) atan2(y, x);
49 }
50 #endif
51
52 #ifndef HAVE_ATANF
53 float
54 atanf(float x)
55 {
56   return (float) atan(x);
57 }
58 #endif
59
60 #ifndef HAVE_CEILF
61 float
62 ceilf(float x)
63 {
64   return (float) ceil(x);
65 }
66 #endif
67
68 #ifndef HAVE_COPYSIGNF
69 float
70 copysignf(float x, float y)
71 {
72   return (float) copysign(x, y);
73 }
74 #endif
75
76 #ifndef HAVE_COSF
77 float
78 cosf(float x)
79 {
80   return (float) cos(x);
81 }
82 #endif
83
84 #ifndef HAVE_COSHF
85 float
86 coshf(float x)
87 {
88   return (float) cosh(x);
89 }
90 #endif
91
92 #ifndef HAVE_EXPF
93 float
94 expf(float x)
95 {
96   return (float) exp(x);
97 }
98 #endif
99
100 #ifndef HAVE_FABSF
101 float
102 fabsf(float x)
103 {
104   return (float) fabs(x);
105 }
106 #endif
107
108 #ifndef HAVE_FLOORF
109 float
110 floorf(float x)
111 {
112   return (float) floor(x);
113 }
114 #endif
115
116 #ifndef HAVE_FREXPF
117 float
118 frexpf(float x, int *exp)
119 {
120   return (float) frexp(x, exp);
121 }
122 #endif
123
124 #ifndef HAVE_HYPOTF
125 float
126 hypotf(float x, float y)
127 {
128   return (float) hypot(x, y);
129 }
130 #endif
131
132 #ifndef HAVE_LOGF
133 float
134 logf(float x)
135 {
136   return (float) log(x);
137 }
138 #endif
139
140 #ifndef HAVE_LOG10F
141 float
142 log10f(float x)
143 {
144   return (float) log10(x);
145 }
146 #endif
147
148 #ifndef HAVE_SCALBNF
149 float
150 scalbnf(float x, int y)
151 {
152   return (float) scalbn(x, y);
153 }
154 #endif
155
156 #ifndef HAVE_SINF
157 float
158 sinf(float x)
159 {
160   return (float) sin(x);
161 }
162 #endif
163
164 #ifndef HAVE_SINHF
165 float
166 sinhf(float x)
167 {
168   return (float) sinh(x);
169 }
170 #endif
171
172 #ifndef HAVE_SQRTF
173 float
174 sqrtf(float x)
175 {
176   return (float) sqrt(x);
177 }
178 #endif
179
180 #ifndef HAVE_TANF
181 float
182 tanf(float x)
183 {
184   return (float) tan(x);
185 }
186 #endif
187
188 #ifndef HAVE_TANHF
189 float
190 tanhf(float x)
191 {
192   return (float) tanh(x);
193 }
194 #endif
195
196 #ifndef HAVE_NEXTAFTERF
197 /* This is a portable implementation of nextafterf that is intended to be
198    independent of the floating point format or its in memory representation.
199    This implementation works correctly with denormalized values.  */
200 float
201 nextafterf(float x, float y)
202 {
203   /* This variable is marked volatile to avoid excess precision problems
204      on some platforms, including IA-32.  */
205   volatile float delta;
206   float absx, denorm_min;
207
208   if (isnan(x) || isnan(y))
209     return x + y;
210   if (x == y)
211     return x;
212
213   /* absx = fabsf (x);  */
214   absx = (x < 0.0) ? -x : x;
215
216   /* __FLT_DENORM_MIN__ is non-zero iff the target supports denormals.  */
217   if (__FLT_DENORM_MIN__ == 0.0f)
218     denorm_min = __FLT_MIN__;
219   else
220     denorm_min = __FLT_DENORM_MIN__;
221
222   if (absx < __FLT_MIN__)
223     delta = denorm_min;
224   else
225     {
226       float frac;
227       int exp;
228
229       /* Discard the fraction from x.  */
230       frac = frexpf (absx, &exp);
231       delta = scalbnf (0.5f, exp);
232
233       /* Scale x by the epsilon of the representation.  By rights we should
234          have been able to combine this with scalbnf, but some targets don't
235          get that correct with denormals.  */
236       delta *= __FLT_EPSILON__;
237
238       /* If we're going to be reducing the absolute value of X, and doing so
239          would reduce the exponent of X, then the delta to be applied is
240          one exponent smaller.  */
241       if (frac == 0.5f && (y < x) == (x > 0))
242         delta *= 0.5f;
243
244       /* If that underflows to zero, then we're back to the minimum.  */
245       if (delta == 0.0f)
246         delta = denorm_min;
247     }
248
249   if (y < x)
250     delta = -delta;
251
252   return x + delta;
253 }
254 #endif
255
256
257 #ifndef HAVE_POWF
258 float
259 powf(float x, float y)
260 {
261   return (float) pow(x, y);
262 }
263 #endif
264
265 /* Note that if fpclassify is not defined, then NaN is not handled */
266
267 /* Algorithm by Steven G. Kargl.  */
268
269 #ifndef HAVE_ROUND
270 /* Round to nearest integral value.  If the argument is halfway between two
271    integral values then round away from zero.  */
272
273 double
274 round(double x)
275 {
276    double t;
277 #if defined(fpclassify)
278    int i;
279    i = fpclassify(x);
280    if (i == FP_INFINITE || i == FP_NAN)
281      return (x);
282 #endif
283
284    if (x >= 0.0) 
285     {
286       t = ceil(x);
287       if (t - x > 0.5)
288         t -= 1.0;
289       return (t);
290     } 
291    else 
292     {
293       t = ceil(-x);
294       if (t + x > 0.5)
295         t -= 1.0;
296       return (-t);
297     }
298 }
299 #endif
300
301 #ifndef HAVE_ROUNDF
302 /* Round to nearest integral value.  If the argument is halfway between two
303    integral values then round away from zero.  */
304
305 float
306 roundf(float x)
307 {
308    float t;
309 #if defined(fpclassify)
310    int i;
311
312    i = fpclassify(x);
313    if (i == FP_INFINITE || i == FP_NAN)
314      return (x);
315 #endif
316
317    if (x >= 0.0) 
318     {
319       t = ceilf(x);
320       if (t - x > 0.5)
321         t -= 1.0;
322       return (t);
323     } 
324    else 
325     {
326       t = ceilf(-x);
327       if (t + x > 0.5)
328         t -= 1.0;
329       return (-t);
330     }
331 }
332 #endif