lsfloat.lsp

   1 ;;; -*- mode: lisp -*-
   2 ;;; Copyright (c) 2005--2007, by A.J. Rossini <blindglobe@gmail.com>
   3 ;;; See COPYRIGHT file for any additional restrictions (BSD license).
   4 ;;; Since 1991, ANSI was finally finished.  Edited for ANSI Common Lisp.
   5
   6 ;;;; lsfloat -- Floating point specs and transcendental functions
   7 ;;;;
   8 ;;;; Copyright (c) 1991, by Luke Tierney. Permission is granted for
   9 ;;;; unrestricted use.
  10 ;;;;
  11 ;;;; Common Lisp allows for four different floating point types that need
  12 ;;;; not be distinct. For statistical work, the type I prefer to use is
  13 ;;;; the one that is closest to a C double. This type is named stat-float.
  14 ;;;; By setting the variable *read-default-float-format* to this type, you
  15 ;;;; insure that data entered as floating point data is read in with this
  16 ;;;; type. The problem arises with data read as integers that is passed to
  17 ;;;; a transcendental, like sqrt. Floating point contagion rules say these
  18 ;;;; integers are to be converted to type single-float. Unless single-float
  19 ;;;; is equivalent to C double, as it is in Mac CL and KCL, this is not
  20 ;;;; what I want. Hence this file redefines the transcendentals to first
  21 ;;;; coerce their arguments to stat-float before applying the built-in
  22 ;;;; functions.
  23 ;;;;
  24 ;;;; No actual modifications to the transcendentals are needed if
  25 ;;;; single-float is the same as stat-float. The fearure
  26 ;;;; :stat-float-is-double-float is used to indicate this.
  27 ;;;;
  28 ;;;; KCL NOTE:
  29 ;;;; In (A)KCL the type short-float corresponds to C float and the types
  30 ;;;; single-float, double-float and long-float correspond to C double.
  31 ;;;; But in the implementation of the transcendentals (A)KCL coerces
  32 ;;;; rationals to short-float, not single-float. CLtL1 is a little vague
  33 ;;;; on this (it talks about "single precision") but CLtL2 clarifies that
  34 ;;;; rationals should produce single-float results. So (A)KCL is wrong, at
  35 ;;;; least relative to the clarification in CLtL2. I therefore decided
  36 ;;;; to fix (A)KCL in files c/num_sfun.c and lsp/numlib.lsp. If these
  37 ;;;; fixes are applied, the feature :stat-float-is-double-float should be
  38 ;;;; defined.
  39
  40
  41 ;;; ======== From the CLHS ==========
  42 ;;
  43 ;; Type SHORT-FLOAT, SINGLE-FLOAT, DOUBLE-FLOAT, LONG-FLOAT
  44 ;;
  45 ;; Supertypes:
  46 ;;
  47 ;; short-float: short-float, float, real, number, t
  48 ;;
  49 ;;
  50 ;; single-float: single-float, float, real, number, t
  51 ;;
  52 ;; double-float: double-float, float, real, number, t
  53 ;;
  54 ;; long-float: long-float, float, real, number, t
  55 ;;
  56 ;; Description: For the four defined subtypes of type float, it is true
  57 ;; that intermediate between the type short-float and the type long-float
  58 ;; are the type single-float and the type double-float. The precise
  59 ;; definition of these categories is implementation-defined. The
  60 ;; precision (measured in ``bits'', computed as p log 2b) and the
  61 ;; exponent size (also measured in ``bits,'' computed as log 2(n+1),
  62 ;; where n is the maximum exponent value) is recommended to be at least
  63 ;; as great as the values in the next figure. Each of the defined
  64 ;; subtypes of type float might or might not have a minus zero.
  65 ;;
  66 ;;
  67 ;; Format  Minimum Precision  Minimum Exponent Size
  68 ;; ----------
  69 ;; Short   13 bits            5 bits
  70 ;; Single  24 bits            8 bits
  71 ;; Double  50 bits            8 bits
  72 ;; Long    50 bits            8 bits
  73 ;;
  74 ;; Figure 12-12. Recommended Minimum Floating-Point Precision and Exponent Size
  75 ;;
  76 ;;  There can be fewer than four internal representations for floats. If
  77 ;; there are fewer distinct representations, the following rules apply:
  78 ;;
  79 ;;  If there is only one, it is the type single-float. In this
  80 ;; representation, an object is simultaneously of types single-float,
  81 ;; double-float, short-float, and long-float.
  82 ;;
  83 ;;  Two internal representations can be arranged in either of the
  84 ;; following ways:
  85 ;;
  86 ;;  Two types are provided: single-float and short-float. An object is
  87 ;; simultaneously of types single-float, double-float, and long-float.
  88 ;;
  89 ;;  Two types are provided: single-float and double-float. An object is
  90 ;; simultaneously of types single-float and short-float, or double-float
  91 ;; and long-float.
  92 ;;
  93 ;;  Three internal representations can be arranged in either of the
  94 ;; following ways:
  95 ;;
  96 ;;  Three types are provided: short-float, single-float, and
  97 ;; double-float. An object can simultaneously be of type double-float and
  98 ;; long-float.
  99 ;;
 100 ;;  Three types are provided: single-float, double-float, and
 101 ;; long-float. An object can simultaneously be of types single-float and
 102 ;; short-float.
 103
 104 ;;; AJR: experimental branch ansifloat:
 105 ;;; what is the problem with just using ANSI CL's floating defs?
 106
 107 ;;; Package Setup
 108
 109 (defpackage :lisp-stat-float
 110   (:use :common-lisp)
 111   (:export machine-epsilon))
 112
 113 (in-package #:lisp-stat-float)
 114
 115 (defparameter machine-epsilon
 116   (do ((epsilon (float 1.0 +stat-float-template+)
 117                 (/ epsilon 2.0)))
 118       ((= (+ 1.0 (/ epsilon 2.0)) 1.0) epsilon)))
 119
 120 (setf *read-default-float-format* 'long-float)
 121