OpenCores
URL https://opencores.org/ocsvn/openrisc/openrisc/trunk

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-dev/] [or1k-gcc/] [libgcc/] [config/] [rs6000/] [ibm-ldouble.c] - Blame information for rev 737

Go to most recent revision | Details | Compare with Previous | View Log

Line No. Rev Author Line
1 734 jeremybenn
/* 128-bit long double support routines for Darwin.
2
   Copyright (C) 1993, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2011
3
   Free Software Foundation, Inc.
4
 
5
This file is part of GCC.
6
 
7
GCC is free software; you can redistribute it and/or modify it under
8
the terms of the GNU General Public License as published by the Free
9
Software Foundation; either version 3, or (at your option) any later
10
version.
11
 
12
GCC is distributed in the hope that it will be useful, but WITHOUT ANY
13
WARRANTY; without even the implied warranty of MERCHANTABILITY or
14
FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
15
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
 
27
/* Implementations of floating-point long double basic arithmetic
28
   functions called by the IBM C compiler when generating code for
29
   PowerPC platforms.  In particular, the following functions are
30
   implemented: __gcc_qadd, __gcc_qsub, __gcc_qmul, and __gcc_qdiv.
31
   Double-double algorithms are based on the paper "Doubled-Precision
32
   IEEE Standard 754 Floating-Point Arithmetic" by W. Kahan, February 26,
33
   1987.  An alternative published reference is "Software for
34
   Doubled-Precision Floating-Point Computations", by Seppo Linnainmaa,
35
   ACM TOMS vol 7 no 3, September 1981, pages 272-283.  */
36
 
37
/* Each long double is made up of two IEEE doubles.  The value of the
38
   long double is the sum of the values of the two parts.  The most
39
   significant part is required to be the value of the long double
40
   rounded to the nearest double, as specified by IEEE.  For Inf
41
   values, the least significant part is required to be one of +0.0 or
42
   -0.0.  No other requirements are made; so, for example, 1.0 may be
43
   represented as (1.0, +0.0) or (1.0, -0.0), and the low part of a
44
   NaN is don't-care.
45
 
46
   This code currently assumes big-endian.  */
47
 
48
#if (!defined (__LITTLE_ENDIAN__) \
49
     && (defined (__MACH__) || defined (__powerpc__) || defined (_AIX)))
50
 
51
#define fabs(x) __builtin_fabs(x)
52
#define isless(x, y) __builtin_isless (x, y)
53
#define inf() __builtin_inf()
54
 
55
#define unlikely(x) __builtin_expect ((x), 0)
56
 
57
#define nonfinite(a) unlikely (! isless (fabs (a), inf ()))
58
 
59
/* Define ALIASNAME as a strong alias for NAME.  */
60
# define strong_alias(name, aliasname) _strong_alias(name, aliasname)
61
# define _strong_alias(name, aliasname) \
62
  extern __typeof (name) aliasname __attribute__ ((alias (#name)));
63
 
64
/* All these routines actually take two long doubles as parameters,
65
   but GCC currently generates poor code when a union is used to turn
66
   a long double into a pair of doubles.  */
67
 
68
long double __gcc_qadd (double, double, double, double);
69
long double __gcc_qsub (double, double, double, double);
70
long double __gcc_qmul (double, double, double, double);
71
long double __gcc_qdiv (double, double, double, double);
72
 
73
#if defined __ELF__ && defined SHARED \
74
    && (defined __powerpc64__ || !(defined __linux__ || defined __gnu_hurd__))
75
/* Provide definitions of the old symbol names to satisfy apps and
76
   shared libs built against an older libgcc.  To access the _xlq
77
   symbols an explicit version reference is needed, so these won't
78
   satisfy an unadorned reference like _xlqadd.  If dot symbols are
79
   not needed, the assembler will remove the aliases from the symbol
80
   table.  */
81
__asm__ (".symver __gcc_qadd,_xlqadd@GCC_3.4\n\t"
82
         ".symver __gcc_qsub,_xlqsub@GCC_3.4\n\t"
83
         ".symver __gcc_qmul,_xlqmul@GCC_3.4\n\t"
84
         ".symver __gcc_qdiv,_xlqdiv@GCC_3.4\n\t"
85
         ".symver .__gcc_qadd,._xlqadd@GCC_3.4\n\t"
86
         ".symver .__gcc_qsub,._xlqsub@GCC_3.4\n\t"
87
         ".symver .__gcc_qmul,._xlqmul@GCC_3.4\n\t"
88
         ".symver .__gcc_qdiv,._xlqdiv@GCC_3.4");
89
#endif
90
 
91
typedef union
92
{
93
  long double ldval;
94
  double dval[2];
95
} longDblUnion;
96
 
97
/* Add two 'long double' values and return the result.  */
98
long double
99
__gcc_qadd (double a, double aa, double c, double cc)
100
{
101
  longDblUnion x;
102
  double z, q, zz, xh;
103
 
104
  z = a + c;
105
 
106
  if (nonfinite (z))
107
    {
108
      z = cc + aa + c + a;
109
      if (nonfinite (z))
110
        return z;
111
      x.dval[0] = z;  /* Will always be DBL_MAX.  */
112
      zz = aa + cc;
113
      if (fabs(a) > fabs(c))
114
        x.dval[1] = a - z + c + zz;
115
      else
116
        x.dval[1] = c - z + a + zz;
117
    }
118
  else
119
    {
120
      q = a - z;
121
      zz = q + c + (a - (q + z)) + aa + cc;
122
 
123
      /* Keep -0 result.  */
124
      if (zz == 0.0)
125
        return z;
126
 
127
      xh = z + zz;
128
      if (nonfinite (xh))
129
        return xh;
130
 
131
      x.dval[0] = xh;
132
      x.dval[1] = z - xh + zz;
133
    }
134
  return x.ldval;
135
}
136
 
137
long double
138
__gcc_qsub (double a, double b, double c, double d)
139
{
140
  return __gcc_qadd (a, b, -c, -d);
141
}
142
 
143
#ifdef __NO_FPRS__
144
static double fmsub (double, double, double);
145
#endif
146
 
147
long double
148
__gcc_qmul (double a, double b, double c, double d)
149
{
150
  longDblUnion z;
151
  double t, tau, u, v, w;
152
 
153
  t = a * c;                    /* Highest order double term.  */
154
 
155
  if (unlikely (t == 0)          /* Preserve -0.  */
156
      || nonfinite (t))
157
    return t;
158
 
159
  /* Sum terms of two highest orders. */
160
 
161
  /* Use fused multiply-add to get low part of a * c.  */
162
#ifndef __NO_FPRS__
163
  asm ("fmsub %0,%1,%2,%3" : "=f"(tau) : "f"(a), "f"(c), "f"(t));
164
#else
165
  tau = fmsub (a, c, t);
166
#endif
167
  v = a*d;
168
  w = b*c;
169
  tau += v + w;     /* Add in other second-order terms.  */
170
  u = t + tau;
171
 
172
  /* Construct long double result.  */
173
  if (nonfinite (u))
174
    return u;
175
  z.dval[0] = u;
176
  z.dval[1] = (t - u) + tau;
177
  return z.ldval;
178
}
179
 
180
long double
181
__gcc_qdiv (double a, double b, double c, double d)
182
{
183
  longDblUnion z;
184
  double s, sigma, t, tau, u, v, w;
185
 
186
  t = a / c;                    /* highest order double term */
187
 
188
  if (unlikely (t == 0)          /* Preserve -0.  */
189
      || nonfinite (t))
190
    return t;
191
 
192
  /* Finite nonzero result requires corrections to the highest order term.  */
193
 
194
  s = c * t;                    /* (s,sigma) = c*t exactly.  */
195
  w = -(-b + d * t);    /* Written to get fnmsub for speed, but not
196
                           numerically necessary.  */
197
 
198
  /* Use fused multiply-add to get low part of c * t.    */
199
#ifndef __NO_FPRS__
200
  asm ("fmsub %0,%1,%2,%3" : "=f"(sigma) : "f"(c), "f"(t), "f"(s));
201
#else
202
  sigma = fmsub (c, t, s);
203
#endif
204
  v = a - s;
205
 
206
  tau = ((v-sigma)+w)/c;   /* Correction to t.  */
207
  u = t + tau;
208
 
209
  /* Construct long double result.  */
210
  if (nonfinite (u))
211
    return u;
212
  z.dval[0] = u;
213
  z.dval[1] = (t - u) + tau;
214
  return z.ldval;
215
}
216
 
217
#if defined (_SOFT_DOUBLE) && defined (__LONG_DOUBLE_128__)
218
 
219
long double __gcc_qneg (double, double);
220
int __gcc_qeq (double, double, double, double);
221
int __gcc_qne (double, double, double, double);
222
int __gcc_qge (double, double, double, double);
223
int __gcc_qle (double, double, double, double);
224
long double __gcc_stoq (float);
225
long double __gcc_dtoq (double);
226
float __gcc_qtos (double, double);
227
double __gcc_qtod (double, double);
228
int __gcc_qtoi (double, double);
229
unsigned int __gcc_qtou (double, double);
230
long double __gcc_itoq (int);
231
long double __gcc_utoq (unsigned int);
232
 
233
extern int __eqdf2 (double, double);
234
extern int __ledf2 (double, double);
235
extern int __gedf2 (double, double);
236
 
237
/* Negate 'long double' value and return the result.    */
238
long double
239
__gcc_qneg (double a, double aa)
240
{
241
  longDblUnion x;
242
 
243
  x.dval[0] = -a;
244
  x.dval[1] = -aa;
245
  return x.ldval;
246
}
247
 
248
/* Compare two 'long double' values for equality.  */
249
int
250
__gcc_qeq (double a, double aa, double c, double cc)
251
{
252
  if (__eqdf2 (a, c) == 0)
253
    return __eqdf2 (aa, cc);
254
  return 1;
255
}
256
 
257
strong_alias (__gcc_qeq, __gcc_qne);
258
 
259
/* Compare two 'long double' values for less than or equal.  */
260
int
261
__gcc_qle (double a, double aa, double c, double cc)
262
{
263
  if (__eqdf2 (a, c) == 0)
264
    return __ledf2 (aa, cc);
265
  return __ledf2 (a, c);
266
}
267
 
268
strong_alias (__gcc_qle, __gcc_qlt);
269
 
270
/* Compare two 'long double' values for greater than or equal.  */
271
int
272
__gcc_qge (double a, double aa, double c, double cc)
273
{
274
  if (__eqdf2 (a, c) == 0)
275
    return __gedf2 (aa, cc);
276
  return __gedf2 (a, c);
277
}
278
 
279
strong_alias (__gcc_qge, __gcc_qgt);
280
 
281
/* Convert single to long double.  */
282
long double
283
__gcc_stoq (float a)
284
{
285
  longDblUnion x;
286
 
287
  x.dval[0] = (double) a;
288
  x.dval[1] = 0.0;
289
 
290
  return x.ldval;
291
}
292
 
293
/* Convert double to long double.  */
294
long double
295
__gcc_dtoq (double a)
296
{
297
  longDblUnion x;
298
 
299
  x.dval[0] = a;
300
  x.dval[1] = 0.0;
301
 
302
  return x.ldval;
303
}
304
 
305
/* Convert long double to single.  */
306
float
307
__gcc_qtos (double a, double aa __attribute__ ((__unused__)))
308
{
309
  return (float) a;
310
}
311
 
312
/* Convert long double to double.  */
313
double
314
__gcc_qtod (double a, double aa __attribute__ ((__unused__)))
315
{
316
  return a;
317
}
318
 
319
/* Convert long double to int.  */
320
int
321
__gcc_qtoi (double a, double aa)
322
{
323
  double z = a + aa;
324
  return (int) z;
325
}
326
 
327
/* Convert long double to unsigned int.  */
328
unsigned int
329
__gcc_qtou (double a, double aa)
330
{
331
  double z = a + aa;
332
  return (unsigned int) z;
333
}
334
 
335
/* Convert int to long double.  */
336
long double
337
__gcc_itoq (int a)
338
{
339
  return __gcc_dtoq ((double) a);
340
}
341
 
342
/* Convert unsigned int to long double.  */
343
long double
344
__gcc_utoq (unsigned int a)
345
{
346
  return __gcc_dtoq ((double) a);
347
}
348
 
349
#endif
350
 
351
#ifdef __NO_FPRS__
352
 
353
int __gcc_qunord (double, double, double, double);
354
 
355
extern int __eqdf2 (double, double);
356
extern int __unorddf2 (double, double);
357
 
358
/* Compare two 'long double' values for unordered.  */
359
int
360
__gcc_qunord (double a, double aa, double c, double cc)
361
{
362
  if (__eqdf2 (a, c) == 0)
363
    return __unorddf2 (aa, cc);
364
  return __unorddf2 (a, c);
365
}
366
 
367
#include "soft-fp/soft-fp.h"
368
#include "soft-fp/double.h"
369
#include "soft-fp/quad.h"
370
 
371
/* Compute floating point multiply-subtract with higher (quad) precision.  */
372
static double
373
fmsub (double a, double b, double c)
374
{
375
    FP_DECL_EX;
376
    FP_DECL_D(A);
377
    FP_DECL_D(B);
378
    FP_DECL_D(C);
379
    FP_DECL_Q(X);
380
    FP_DECL_Q(Y);
381
    FP_DECL_Q(Z);
382
    FP_DECL_Q(U);
383
    FP_DECL_Q(V);
384
    FP_DECL_D(R);
385
    double r;
386
    long double u, x, y, z;
387
 
388
    FP_INIT_ROUNDMODE;
389
    FP_UNPACK_RAW_D (A, a);
390
    FP_UNPACK_RAW_D (B, b);
391
    FP_UNPACK_RAW_D (C, c);
392
 
393
    /* Extend double to quad.  */
394
#if (2 * _FP_W_TYPE_SIZE) < _FP_FRACBITS_Q
395
    FP_EXTEND(Q,D,4,2,X,A);
396
    FP_EXTEND(Q,D,4,2,Y,B);
397
    FP_EXTEND(Q,D,4,2,Z,C);
398
#else
399
    FP_EXTEND(Q,D,2,1,X,A);
400
    FP_EXTEND(Q,D,2,1,Y,B);
401
    FP_EXTEND(Q,D,2,1,Z,C);
402
#endif
403
    FP_PACK_RAW_Q(x,X);
404
    FP_PACK_RAW_Q(y,Y);
405
    FP_PACK_RAW_Q(z,Z);
406
    FP_HANDLE_EXCEPTIONS;
407
 
408
    /* Multiply.  */
409
    FP_INIT_ROUNDMODE;
410
    FP_UNPACK_Q(X,x);
411
    FP_UNPACK_Q(Y,y);
412
    FP_MUL_Q(U,X,Y);
413
    FP_PACK_Q(u,U);
414
    FP_HANDLE_EXCEPTIONS;
415
 
416
    /* Subtract.  */
417
    FP_INIT_ROUNDMODE;
418
    FP_UNPACK_SEMIRAW_Q(U,u);
419
    FP_UNPACK_SEMIRAW_Q(Z,z);
420
    FP_SUB_Q(V,U,Z);
421
 
422
    /* Truncate quad to double.  */
423
#if (2 * _FP_W_TYPE_SIZE) < _FP_FRACBITS_Q
424
    V_f[3] &= 0x0007ffff;
425
    FP_TRUNC(D,Q,2,4,R,V);
426
#else
427
    V_f1 &= 0x0007ffffffffffffL;
428
    FP_TRUNC(D,Q,1,2,R,V);
429
#endif
430
    FP_PACK_SEMIRAW_D(r,R);
431
    FP_HANDLE_EXCEPTIONS;
432
 
433
    return r;
434
}
435
 
436
#endif
437
 
438
#endif

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.