htswanalysis/MACS/lib/gsl/gsl-1.11/ieee-utils/fp-tru64.c

   1 /* ieee-utils/fp-tru64.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
  21 /*
  22  * Under Compaq's Unix with the silly name, read the man pages for read_rnd,
  23  * write_rnd, and ieee(3) for more information on the functions used here.
  24  *
  25  * Note that enabling control of dynamic rounding mode (via write_rnd) requires
  26  * that you pass a special flag to your C compiler.  For Compaq's C compiler
  27  * the flag is `-fprm d', for gcc it's `-mfp-rounding-mode=d'.
  28  *
  29  * Enabling the trap control (via ieee_set_fp_control) also requires a
  30  * flag be passed to the C compiler.  The flag for Compaq's C compiler
  31  * is `-ieee' and for gcc it's `-mieee'.
  32
  33  * We have not implemented the `inexact' case, since it is rarely used
  34  * and requires the library being built with an additional compiler
  35  * flag that can degrade performance for everything else. If you need
  36  * to add support for `inexact' the relevant flag for Compaq's
  37  * compiler is `-ieee_with_inexact', and the flag for gcc is
  38  * `-mieee-with-inexact'.
  39  *
  40  * Problem have been reported with the "fixed" float.h installed with
  41  * gcc-2.95 lacking some of the definitions in the system float.h (the
  42  * symptoms are errors like: `FP_RND_RN' undeclared). To work around
  43  * this we can include the system float.h before the gcc version, e.g.
  44  *
  45  *  #include "/usr/include/float.h"
  46  *  #include <float.h>
  47  */
  48
  49 #include <float.h>
  50
  51 #ifndef FP_RND_RN
  52 #  undef _FLOAT_H_
  53 #  include "/usr/include/float.h"
  54 #  undef _FLOAT_H_
  55 #  include <float.h>
  56 #endif
  57
  58 #include <machine/fpu.h>
  59 #include <stdio.h>
  60 #include <gsl/gsl_ieee_utils.h>
  61 #include <gsl/gsl_errno.h>
  62
  63 int
  64 gsl_ieee_set_mode (int precision, int rounding, int exception_mask)
  65 {
  66   unsigned long int mode = 0 ;
  67   unsigned int    rnd  = 0 ;
  68
  69 /* I'm actually not completely sure that the alpha only supports default
  70  * precisions rounding, but I couldn't find any information regarding this, so
  71  * it seems safe to assume this for now until it's proven otherwise.
  72  */
  73
  74   switch (precision)
  75     {
  76     case GSL_IEEE_SINGLE_PRECISION:
  77       GSL_ERROR ("Tru64 Unix on the alpha only supports default precision rounding",
  78                  GSL_EUNSUP) ;
  79       break ;
  80     case GSL_IEEE_DOUBLE_PRECISION:
  81       GSL_ERROR ("Tru64 Unix on the alpha only supports default precision rounding",
  82                  GSL_EUNSUP) ;
  83       break ;
  84     case GSL_IEEE_EXTENDED_PRECISION:
  85       GSL_ERROR ("Tru64 Unix on the alpha only supports default precision rounding",
  86                  GSL_EUNSUP) ;
  87       break ;
  88     }
  89
  90
  91   switch (rounding)
  92     {
  93     case GSL_IEEE_ROUND_TO_NEAREST:
  94       rnd = FP_RND_RN ;
  95       write_rnd (rnd) ;
  96       break ;
  97     case GSL_IEEE_ROUND_DOWN:
  98       rnd = FP_RND_RM ;
  99       write_rnd (rnd) ;
 100       break ;
 101     case GSL_IEEE_ROUND_UP:
 102       rnd = FP_RND_RP ;
 103       write_rnd (rnd) ;
 104       break ;
 105     case GSL_IEEE_ROUND_TO_ZERO:
 106       rnd = FP_RND_RZ ;
 107       write_rnd (rnd) ;
 108       break ;
 109     default:
 110       rnd = FP_RND_RN ;
 111       write_rnd (rnd) ;
 112     }
 113
 114   /* Turn on all the exceptions apart from 'inexact' */
 115
 116   /* from the ieee(3) man page:
 117    * IEEE_TRAP_ENABLE_INV       ->      Invalid operation
 118    * IEEE_TRAP_ENABLE_DZE       ->      Divide by 0
 119    * IEEE_TRAP_ENABLE_OVF       ->      Overflow
 120    * IEEE_TRAP_ENABLE_UNF       ->      Underflow
 121    * IEEE_TRAP_ENABLE_INE       ->      Inexact (requires special option to C compiler)
 122    * IEEE_TRAP_ENABLE_DNO       ->      denormal operand
 123    * Note: IEEE_TRAP_ENABLE_DNO is not supported on OSF 3.x or Digital Unix
 124    * 4.0 - 4.0d(?).
 125    * IEEE_TRAP_ENABLE_MASK      ->      mask of all the trap enables
 126    * IEEE_MAP_DMZ                       ->      map denormal inputs to zero
 127    * IEEE_MAP_UMZ                       ->      map underflow results to zero
 128    */
 129
 130   mode = IEEE_TRAP_ENABLE_INV | IEEE_TRAP_ENABLE_DZE | IEEE_TRAP_ENABLE_OVF
 131                 | IEEE_TRAP_ENABLE_UNF ;
 132
 133   if (exception_mask & GSL_IEEE_MASK_INVALID)
 134     mode &= ~ IEEE_TRAP_ENABLE_INV ;
 135
 136   if (exception_mask & GSL_IEEE_MASK_DENORMALIZED)
 137     {
 138 #ifdef IEEE_TRAP_ENABLE_DNO
 139      mode &= ~ IEEE_TRAP_ENABLE_DNO ;
 140 #else
 141      GSL_ERROR ("Sorry, this version of Digital Unix does not support denormalized operands", GSL_EUNSUP) ;
 142 #endif
 143     }
 144
 145   if (exception_mask & GSL_IEEE_MASK_DIVISION_BY_ZERO)
 146     mode &= ~ IEEE_TRAP_ENABLE_DZE ;
 147
 148   if (exception_mask & GSL_IEEE_MASK_OVERFLOW)
 149     mode &= ~ IEEE_TRAP_ENABLE_OVF ;
 150
 151   if (exception_mask & GSL_IEEE_MASK_UNDERFLOW)
 152     mode &=  ~ IEEE_TRAP_ENABLE_UNF ;
 153
 154   if (exception_mask & GSL_IEEE_TRAP_INEXACT)
 155     {
 156       /* To implement this would require a special flag to the C
 157        compiler which can cause degraded performance */
 158
 159       GSL_ERROR ("Sorry, GSL does not implement trap-inexact for Tru64 Unix on the alpha - see fp-tru64.c for details", GSL_EUNSUP) ;
 160
 161       /* In case you need to add it, the appropriate line would be
 162        *
 163        *  mode |= IEEE_TRAP_ENABLE_INE ;
 164        *
 165        */
 166
 167     }
 168   else
 169     {
 170       mode &= ~ IEEE_TRAP_ENABLE_INE ;
 171     }
 172
 173   ieee_set_fp_control (mode) ;
 174
 175   return GSL_SUCCESS ;
 176 }