ghc/tests/programs/north_lias/LIAS.lhs

   1 % lias.lhs - Language Independent Arithmetic Standard in Haskell
   2
   3 % @(#)LIAS.lhs  1.11 dated 92/12/07 at 15:01:23
   4
   5 \documentstyle[a4wide,11pt,times]{article}
   6
   7 \title{Haskell and the Language Independent Arithmetic Standard}
   8 \author{N D North\\
   9         National Physical Laboratory\\
  10         Teddington, TW11 0LW, UK.\\
  11         {\tt ndn@seg.npl.co.uk}}
  12
  13 % Some macros lifted from elsewhere to make this more standalone.
  14 \makeatletter
  15 % INLINE PROGRAM CODE
  16 %
  17 % \prog{foo} sets its argument in typewriter font.
  18 \def\prog#1{\ifmmode\mbox{\tt #1}\else{\tt #1}\fi}
  19
  20 % NEWVERBATIM (from iso.sty)
  21 %
  22 % \newverbatim{foo} creates a new environment, foo, which behaves exactly
  23 % like the verbatim environment except that it is delimited by
  24 % \begin{foo} ... \end{foo}.
  25 % See the VERBATIM section of latex.tex for the inspiration behind this.
  26 %
  27 \def\newverbatim#1{\expandafter\def\csname #1\endcsname{%
  28 \@verbatim \frenchspacing\@vobeyspaces \csname @x#1verbatim\endcsname}
  29 \expandafter\let\csname end#1\endcsname=\endtrivlist
  30 \new@xverbatim{#1}}
  31
  32 \begingroup \catcode `|=0 \catcode `[= 1
  33 \catcode`]=2 \catcode `\{=12 \catcode `\}=12
  34 \catcode`\\=12
  35 |gdef|new@xverbatim#1[
  36 |expandafter|def|csname @x#1verbatim|endcsname##1\end{#1}[##1|end[#1]]]
  37 |endgroup
  38 \makeatother
  39
  40 \newverbatim{haskell}
  41
  42 % \lias{id} sets an identifier in LIAS font (italic)
  43 \def\lias#1{\mbox{\it #1}}
  44
  45 % \liass{id}{sub} sets the identifier in LIAS font, with the given
  46 %   subscript.
  47 \def\liass#1#2{\mbox{$\lias{#1}_{#2}$}}
  48
  49 % \liasss{id}{sub}{sup} sets the identifier in LIAS font, with the
  50 %   subscript and superscript.
  51 \def\liasss#1#2#3{\mbox{$\lias{#1}_{#2}^{#3}$}}
  52
  53 \begin{document}
  54 \maketitle
  55
  56 \section*{Introduction}
  57
  58 Haskell~\cite{hudak} is intended as an ``industrial strength''
  59 functional programming language and, in partial fulfillment of that
  60 aim, includes a rich set of numeric types and operators.
  61 However, the semantics of numeric operations are rather imprecise, so
  62 that determining the accuracy of numerical analysis programs is impossible
  63 in Haskell, limiting its applicability.
  64 The Language Independent Arithmetic Standard (LIAS)~\cite{lias}
  65 defines the behaviour of numerical operations precisely, yet flexibly
  66 enough that it is compatible with virtually all major arithmetic
  67 implementations including, for example IEEE~754~\cite{ieee754}.
  68
  69 This report examines the extent to which Haskell and LIAS are compatible,
  70 provides a model implementation of LIAS in Haskell, and recommends a
  71 small addition to Haskell to improve compatibility.
  72 The intention is to improve the portability of programs, both between
  73 Haskell implementations and between Haskell and other LIAS-compliant
  74 languages.
  75
  76
  77 \section{Compatibility between Haskell and LIAS}
  78
  79 \subsection{Denormalised numbers}
  80 \label{denorm}
  81 Parameters for LIAS are all available in Haskell, with the exception of
  82 \lias{denorm}, so few problems in principle.
  83
  84 \subsection{Accuracy}
  85 Haskell implementations tend to use arithmetic of underlying system, so
  86 extent to which accuracy complies depends on that of underlying system.
  87
  88 \subsection{Notification}
  89 \label{notification}
  90 Semantics of overflow etc are ``undefined'' so demanding notification is
  91 impossible.
  92 {\em Check what systems actually do.}
  93
  94 \subsection{Integers}
  95 Haskell provides a class \prog{Integral}, whose members are integer
  96 types.
  97 In particular, \prog{Integer} is the type of arbitrary-precision
  98 integers, and \prog{Int} is a type of fixed-precision integers with
  99 range at least $[-2^{29} + 1, 2^{29} - 1]$ and closed under
 100 negation.
 101 Implementations are at liberty provide other integer types.
 102
 103 Both \prog{Integer} and \prog{Int} should comply with LIAS, with the
 104 exception of notifications, as described in Section~\ref{notification}.
 105
 106 \subsection{Floating point}
 107 Haskell provides a class \prog{RealFloat}, whose members are real (as
 108 opposed to complex) floating point numbers.
 109 In particular, \prog{Float} and \prog{Double} are supposed to
 110 be at least equal in range and precision to IEEE single and double
 111 precision respectively.
 112 Implementations are at liberty provide other floating point types.
 113
 114 Both \prog{Float} and \prog{Double} should comply with LIAS, with the
 115 exceptions of notifications, as described in Section~\ref{notification},
 116 and \lias{denorm}, as described in Section~\ref{denorm}.
 117
 118
 119 \section{LIAS in Haskell}
 120
 121 This section provides a model implementation of LIAS in Haskell.
 122 Many of the required parameters and functions already exist, in which
 123 case this section just describes how to obtain them in Haskell.
 124 Others are not part of the standard language, and code is given to
 125 implement these.
 126
 127 The section begins with the module header giving the module name and
 128 exported identifiers.
 129
 130 \begin{haskell}
 131
 132 > module LIAS (
 133 >     emax, emin, denorm,
 134 >     fmax, fminN, fmin, epsilon,
 135 >     signf, succf, predf, ulpf, truncf, roundf, fractpart
 136 >     ) where
 137
 138 \end{haskell}
 139
 140
 141 \subsection{Integers}
 142 The LIAS parameters are: \lias{minint}, \lias{maxint} and \lias{bounded}.
 143 Whether an integer type is bounded or not is part of the language
 144 definition, so this parameter is not available to the user in Haskell.
 145 The minimum and maximum parameters are available for the \prog{Int}
 146 type, and are accessed as follows:
 147 \begin{tabbing}
 148 \lias{minint} \= \prog{minInt} \kill
 149 \lias{minint} \> \prog{minInt} \\
 150 \lias{maxint} \> \prog{maxInt} (\prog{= -minInt})
 151 \end{tabbing}
 152 Note that the Haskell Report (Section 6.8.2) states that
 153 \prog{maxInt = -minInt}, which is compatible with LIAS.
 154
 155 All the integer operations required by LIAS are available in Haskell,
 156 and are accessed as described in the table below:
 157 \begin{tabbing}
 158 mmmmmmmmmmmmmmm \= mmmmmmmmmm \= \kill
 159 \liass{add}{I}     \>  \prog{x + y} \\
 160 \liass{sub}{I}     \>  \prog{x - y} \\
 161 \liass{mul}{I}     \>  \prog{x * y} \\
 162 \liass{div}{I}     \>  \prog{x `div` y} (round to $-\infty$) \\
 163                    \>  \prog{x `quot` y} (round to 0) \\
 164 \liass{rem}{I}     \>  \prog{x `mod` y} (round to $-\infty$) \\
 165                    \>  \prog{x `rem` y} (round to 0) \\
 166 \liass{mod}{I}     \>  \prog{x `mod` y} (this is \liasss{mod}{I}{1}) \\
 167 \liass{neg}{I}     \>  \prog{negate x} \\
 168 \liass{abs}{I}     \>  \prog{abs x} \\
 169 \liass{eq}{I}      \>  \prog{x == y} \\
 170 \liass{neq}{I}     \>  \prog{x /= y} \\
 171 \liass{lss}{I}     \>  \prog{x < y} \\
 172 \liass{leq}{I}     \>  \prog{x <= y} \\
 173 \liass{gtr}{I}     \>  \prog{x > y} \\
 174 \liass{geq}{I}     \>  \prog{x >= y}
 175 \end{tabbing}
 176 The table shows that Haskell provides integer division with rounding
 177 towards $-\infty$ and with rounding towards 0.
 178
 179 \subsection{Floating point}
 180 Haskell provides all the parameters for floating point numbers, except
 181 for \lias{denorm}.
 182 The available parameters are determined as follows:
 183 \begin{tabbing}
 184 \lias{emax}  \=  fst (floatRange x) \= \prog{Int} \kill
 185 \lias{r}     \>  floatRadix x       \> \prog{Integer} \\
 186 \lias{p}     \>  floatDigits x      \> \prog{Int} \\
 187 \lias{emax}  \>  fst (floatRange x) \> \prog{Int} \\
 188 \lias{emin}  \>  snd (floatRange x) \> \prog{Int}
 189 \end{tabbing}
 190 In the table, \prog{x} is an expression of the type for which the
 191 parameter is required.
 192 For example, \prog{floatRadix (1.0 :: Float)} would give the radix of
 193 the \prog{Float} type.
 194 The alternative to this mechanism is to provide a separate set of
 195 identifiers for each floating point type.
 196
 197 For convenience, we provide Haskell identifiers for \lias{emax} and
 198 \lias{emin}.
 199 \begin{haskell}
 200
 201 > emax, emin :: (RealFloat a) => a -> Int
 202 > emax x  =  snd (floatRange x)
 203 > emin x  =  fst (floatRange x)
 204
 205 \end{haskell}
 206
 207 The derived constants require some coding as follows:
 208 \begin{haskell}
 209
 210 > fmax, fminN, fminD, fmin, epsilon :: (RealFloat a) => a -> a
 211
 212 > fmax x  =  encodeFloat (floatRadix x ^ floatDigits x - 1)
 213 >                        (emax x - floatDigits x)
 214
 215 > fminN x  =  encodeFloat 1 (emin x - 1)
 216
 217 > fminD x  =  encodeFloat 1 (emin x - floatDigits x)
 218
 219 > fmin x  =  if denorm x then fminD x else fminN x
 220
 221 > epsilon x  =  encodeFloat 1 (1 - floatDigits x)
 222
 223 \end{haskell}
 224
 225 The definition of \lias{denorm} assumes that the implementation gives
 226 zero on underflow.
 227 The Haskell Report leaves behaviour on underflow undefined, which
 228 makes this definition less than satisfactory and suggests that
 229 \prog{denorm} should be part of the language.
 230 \begin{haskell}
 231
 232 > denorm :: (RealFloat a) => a -> Bool
 233 > denorm x  =  fminN x / fromInteger (floatRadix x) /= 0
 234
 235 \end{haskell}
 236
 237 The floating point operations are listed below, with the syntax for
 238 invoking them.
 239 The operations marked ``$\dagger$'' are not part of Haskell and are
 240 defined later in the LIAS module.
 241 \begin{tabbing}
 242 mmmmmmmmmmmmmmm \= mmmmmmmmmm \= \kill
 243 \liass{add}{F}       \>  \prog{x + y} \\
 244 \liass{sub}{F}       \>  \prog{x - y} \\
 245 \liass{mul}{F}       \>  \prog{x * y} \\
 246 \liass{div}{F}       \>  \prog{x / y} \\
 247 \liass{neg}{F}       \>  \prog{negate x} \\
 248 \liass{abs}{F}       \>  \prog{abs x} \\
 249 \liass{sqrt}{F}      \>  \prog{sqrt x} \\
 250 \liass{sign}{F}      \>  \prog{signf x} \> $\dagger$ \\
 251 \liass{exponent}{F}  \>  \prog{exponent x} \\
 252 \liass{fraction}{F}  \>  \prog{significand x} \\
 253 \liass{scale}{F}     \>  \prog{scaleFloat n x} \\
 254 \liass{succ}{F}      \>  \prog{succf x} \> $\dagger$ \\
 255 \liass{pred}{F}      \>  \prog{predf x} \> $\dagger$ \\
 256 \liass{ulp}{F}       \>  \prog{ulpf x} \> $\dagger$ \\
 257 \liass{trunc}{F}     \>  \prog{truncf x n} \> $\dagger$ \\
 258 \liass{round}{F}     \>  \prog{roundf x n} \> $\dagger$ \\
 259 \liass{intpart}{F}   \>  \prog{truncate x} \\
 260 \liass{fractpart}{F} \>  \prog{snd (properFraction x)} \\
 261 \liass{eq}{F}        \>  \prog{x == y} \\
 262 \liass{neq}{F}       \>  \prog{x /= y} \\
 263 \liass{lss}{F}       \>  \prog{x < y} \\
 264 \liass{leq}{F}       \>  \prog{x <= y} \\
 265 \liass{gtr}{F}       \>  \prog{x > y} \\
 266 \liass{geq}{F}       \>  \prog{x >= y}
 267 \end{tabbing}
 268
 269 The code below provides definitions of the operations marked ``$\dagger$''
 270 and, for convenience, a definition of \liass{fractpart}{F}.
 271
 272 \begin{haskell}
 273
 274 > signf :: (RealFloat a) => a -> a
 275 > signf x | x >= 0  =  1
 276 >         | x <  0  =  -1
 277
 278 \end{haskell}
 279
 280 \prog{floatRadixf} is a useful utility function which gives the floating
 281 point radix as a member of the class \prog{RealFloat}.
 282 \begin{haskell}
 283
 284 > floatRadixf :: (RealFloat a) => a -> a
 285 > floatRadixf x  =  fromInteger (floatRadix x)
 286
 287 \end{haskell}
 288
 289 \begin{haskell}
 290
 291 > succf, predf :: (RealFloat a) => a -> a
 292 > succf x | x == 0          =  fmin x
 293 >         | x == -(fmin x)  =  0
 294 >         | True            =  encodeFloat (m + 1) n
 295 >                              where
 296 >                              (m, n)  =  decodeFloat x
 297
 298 > predf x  =  - succf (- x)
 299
 300 \end{haskell}
 301
 302 \begin{haskell}
 303
 304 > ulpf :: (RealFloat a) => a -> a
 305 > ulpf x | x == 0  =  error "ulpf of 0"
 306 >        | True    =  res (encodeFloat 1 (expf x - floatDigits x))
 307 >                     where
 308 >                     res 0  =  error "ulpf underflow"
 309 >                     res x  =  x
 310
 311 \end{haskell}
 312
 313 \begin{haskell}
 314
 315 > floorf :: (RealFloat a) => a -> a
 316 > floorf x  =  fromInteger (floor x)
 317
 318 > expf :: (RealFloat a) => a -> Int
 319 > expf x  =  if abs x >= fminN x then exponent x else emin x
 320
 321 > truncf :: (RealFloat a) => a -> Int -> a
 322 > truncf x n | n <= 0          =  error "truncf with n <= 0"
 323 >            | j >= eemin - p  =  encodeFloat (i `quot` (r ^ (p - n)))
 324 >                                             (j + p - n)
 325 >            | True            =  encodeFloat (i `quot` (r ^ (eemin - j - n)))
 326 >                                             (eemin -n)
 327 >                                 where
 328 >                                 (i, j)  =  decodeFloat x
 329 >                                 eemin   =  emin x
 330 >                                 r       =  floatRadix x
 331 >                                 p       =  floatDigits x
 332
 333 \end{haskell}
 334
 335 \begin{haskell}
 336
 337 > roundf ::  (RealFloat a) => a -> Int -> a
 338 > roundf x n | n <= 0              =  error "roundf with n <= 0"
 339 >            | n >= floatDigits x  =  x
 340 >            | True                =  signf x * floorf (abs x / y + 0.5) * y
 341 >                                     where
 342 >                                     y  =  encodeFloat 1 (expf x - n)
 343
 344 \end{haskell}
 345
 346 \begin{haskell}
 347
 348 > fractpart :: (RealFloat a) => a -> a
 349 > fractpart x  =  snd (properFraction x)
 350
 351 \end{haskell}
 352
 353 \section{Recommendations}
 354
 355
 356 \begin{thebibliography}{9}
 357 \bibitem{hudak} P Hudak, S Peyton Jones, P Wadler et al.
 358 {\it Report on the Functional Programming Language Haskell, Version 1.1.}
 359 Department of Computing Science, University of Glasgow, August 1991.
 360 \bibitem{ieee754} IEEE Standard for Binary Floating-Point Arithmetic.
 361     ANSI/IEEE Std 754-1985, 1985.
 362 \bibitem{lias} M~Payne, C~Schaffert, and B~A~Wichmann.
 363 {\em The Language Compatible Arithmetic Standard}.
 364  January 1990. ACM SIGPLAN Notices, Vol 25,
 365   pp59-86, and ACM SIGNUM Newsletter, Vol 25, No 1, pp2-43.
 366 \end{thebibliography}
 367 \end{document}