1 /* ieee-utils/fp-aix.c
2  *
3  * Copyright (C) 1996, 1997, 1998, 1999, 2000 Tim Mooney
4  *
5  * This program is free software; you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License as published by
7  * the Free Software Foundation; either version 3 of the License, or (at
8  * your option) any later version.
9  *
10  * This program is distributed in the hope that it will be useful, but
11  * WITHOUT ANY WARRANTY; without even the implied warranty of
12  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
13  * General Public License for more details.
14  *
15  * You should have received a copy of the GNU General Public License
16  * along with this program; if not, write to the Free Software
17  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
18  */
19 
20 #include <math.h>
21 #include <fptrap.h>
22 #include <float.h>
23 #include <gsl/gsl_ieee_utils.h>
24 #include <gsl/gsl_errno.h>
25 
26 /* GCC uses a common float.h for all platforms, and ignores all vendor
27 specific code in float.h. That causes problems with fprnd_t,
28 FP_RND_RZ, FP_RND_RN and FP_RND_RP on AIX.
29 
30 Personally I consider that a bug, and was advised by a GCC developer
31 to report this as a bug, which I did.
32 
33 http://gcc.gnu.org/bugzilla/show_bug.cgi?id=46155
34 
35 However, there is not universal agreement whether this is a gcc bug, so
36 the issue needs to be worked around on AIX.
37 
38 David Kirkby, 26th October 2010. */
39 
40 #if !HAVE_DECL_FPRND_T
41 typedef unsigned short fprnd_t;
42 #endif
43 
44 #ifndef FP_RND_RZ
45 #define FP_RND_RZ       0
46 #endif
47 
48 #ifndef FP_RND_RN
49 #define FP_RND_RN       1
50 #endif
51 
52 
53 #ifndef FP_RND_RP
54 #define FP_RND_RP       2
55 #endif
56 
57 #ifndef FP_RND_RM
58 #define FP_RND_RM       3
59 #endif
60 
61 int
gsl_ieee_set_mode(int precision,int rounding,int exception_mask)62 gsl_ieee_set_mode (int precision, int rounding, int exception_mask)
63 {
64   fptrap_t   mode = 0 ;
65   fprnd_t    rnd  = 0 ;
66 
67   switch (precision)
68     {
69 
70     /* I'm not positive about AIX only supporting default precision rounding,
71      * but this is the best assumption until it's proven otherwise. */
72 
73     case GSL_IEEE_SINGLE_PRECISION:
74       GSL_ERROR ("AIX only supports default precision rounding",
75                  GSL_EUNSUP) ;
76       break ;
77     case GSL_IEEE_DOUBLE_PRECISION:
78       GSL_ERROR ("AIX only supports default precision rounding",
79                  GSL_EUNSUP) ;
80       break ;
81     case GSL_IEEE_EXTENDED_PRECISION:
82       GSL_ERROR ("AIX only supports default precision rounding",
83                  GSL_EUNSUP) ;
84       break ;
85     }
86 
87   switch (rounding)
88     {
89     case GSL_IEEE_ROUND_TO_NEAREST:
90       rnd = FP_RND_RN ;
91       fp_swap_rnd (rnd) ;
92       break ;
93     case GSL_IEEE_ROUND_DOWN:
94       rnd = FP_RND_RM ;
95       fp_swap_rnd (rnd) ;
96       break ;
97     case GSL_IEEE_ROUND_UP:
98       rnd = FP_RND_RP ;
99       fp_swap_rnd (rnd) ;
100       break ;
101     case GSL_IEEE_ROUND_TO_ZERO:
102       rnd = FP_RND_RZ ;
103       fp_swap_rnd (rnd) ;
104       break ;
105     default:
106       rnd = FP_RND_RN ;
107       fp_swap_rnd (rnd) ;
108     }
109 
110   /* Turn on all the exceptions apart from 'inexact' */
111 
112   mode = TRP_INVALID | TRP_DIV_BY_ZERO | TRP_OVERFLOW | TRP_UNDERFLOW ;
113 
114   if (exception_mask & GSL_IEEE_MASK_INVALID)
115     mode &= ~ TRP_INVALID ;
116 
117   if (exception_mask & GSL_IEEE_MASK_DENORMALIZED)
118     {
119       /* do nothing */
120     }
121   else
122     {
123       GSL_ERROR ("AIX does not support the denormalized operand exception. "
124                  "Use 'mask-denormalized' to work around this.",
125                  GSL_EUNSUP) ;
126     }
127 
128   if (exception_mask & GSL_IEEE_MASK_DIVISION_BY_ZERO)
129     mode &= ~ TRP_DIV_BY_ZERO ;
130 
131   if (exception_mask & GSL_IEEE_MASK_OVERFLOW)
132     mode &= ~ TRP_OVERFLOW ;
133 
134   if (exception_mask & GSL_IEEE_MASK_UNDERFLOW)
135     mode &=  ~ TRP_UNDERFLOW ;
136 
137   if (exception_mask & GSL_IEEE_TRAP_INEXACT)
138     {
139       mode |= TRP_INEXACT ;
140     }
141   else
142     {
143       mode &= ~ TRP_INEXACT ;
144     }
145 
146   /* AIX appears to require two steps -- first enable floating point traps
147    * in general... */
148   fp_trap(FP_TRAP_SYNC);
149 
150   /* next, enable the traps we're interested in */
151   fp_enable(mode);
152 
153   return GSL_SUCCESS ;
154 
155 }
156