Описание библиотеки Haskell 98: Числовые функции Описание Haskell 98
наверх | назад | вперед | содержание | предметный указатель функций

14  Числовые функции


module Numeric(fromRat,
               showSigned, showIntAtBase,
               showInt, showOct, showHex,
               readSigned, readInt,
               readDec, readOct, readHex, 
               floatToDigits,
               showEFloat, showFFloat, showGFloat, showFloat, 
               readFloat, lexDigits) where

fromRat        :: (RealFloat a) => Rational -> a

showSigned     :: (Real a) => (a -> ShowS) -> Int -> a -> ShowS
showIntAtBase  :: Integral a => a -> (Int -> Char) -> a -> ShowS
showInt        :: Integral a => a -> ShowS
showOct        :: Integral a => a -> ShowS
showHex        :: Integral a => a -> ShowS

readSigned     :: (Real a) => ReadS a -> ReadS a
readInt        :: (Integral a) =>
                    a -> (Char -> Bool) -> (Char -> Int) -> ReadS a
readDec        :: (Integral a) => ReadS a
readOct        :: (Integral a) => ReadS a
readHex        :: (Integral a) => ReadS a

showEFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showFFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showGFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showFloat      :: (RealFloat a) => a -> ShowS

floatToDigits  :: (RealFloat a) => Integer -> a -> ([Int], Int)

readFloat      :: (RealFrac a) => ReadS a
lexDigits      :: ReadS String 

Эта библиотека содержит числовые функции разных сортов, многие из которых используются в стандартном Prelude.

Далее, напомним следующие определения типов из Prelude:
 
  type ShowS = String -> String
  type ReadS = String -> [(a,String)]

14.1  Функции преобразования величин в строки

14.2  Функции преобразования строк в другие величины

(NB: readInt является "двойственной" для showIntAtBase, а readDec является "двойственной" для showInt. Противоречивые имена этих функций сложились исторически.)

14.3  Прочие функции

14.4  Библиотека Numeric


module Numeric(fromRat,
               showSigned, showIntAtBase,
               showInt, showOct, showHex,
               readSigned, readInt,
               readDec, readOct, readHex, 
               floatToDigits,
               showEFloat, showFFloat, showGFloat, showFloat, 
               readFloat, lexDigits) where

import Char   ( isDigit, isOctDigit, isHexDigit
              , digitToInt, intToDigit )
import Ratio  ( (%), numerator, denominator )
import Array  ( (!), Array, array )

-- Эта функция выполняет преобразование рационального числа в число с плавающей точкой. 
-- Ее следует использовать в экземплярах Fractional классов Float и Double.

fromRat :: (RealFloat a) => Rational -> a
fromRat x = 
    if x == 0 then encodeFloat 0 0              -- Сначала обрабатывает исключительные ситуации
    else if x < 0 then - fromRat' (-x)         
    else fromRat' x

-- Процесс преобразования:
-- Перевести (масштабировать) рациональное число в систему счисления с основанием RealFloat, пока 
-- оно не будет лежать в диапазоне мантиссы (используется функциями decodeFloat/encodeFloat).
-- Затем округлить рациональное число до Integer и закодировать его с помощью экспоненты,
-- полученной при переводе числа в систему счисления.
-- Для того чтобы ускорить процесс масштабирования, мы вычисляем log2 числа, чтобы получить
-- первое приближение экспоненты.
fromRat' :: (RealFloat a) => Rational -> a
fromRat' x = r
  where b = floatRadix r
        p = floatDigits r
        (minExp0, _) = floatRange r
        minExp = minExp0 - p            -- действительная минимальная экспонента
        xMin = toRational (expt b (p-1))
        xMax = toRational (expt b p)
        p0 = (integerLogBase b (numerator x) -
              integerLogBase b (denominator x) - p) `max` minExp
        f = if p0 < 0 then 1 % expt b (-p0) else expt b p0 % 1
        (x', p') = scaleRat (toRational b) minExp xMin xMax p0 (x / f)
        r = encodeFloat (round x') p'

-- Масштабировать x, пока не выполнится условие xMin <= x < xMax или p (экспонета) <= minExp.
scaleRat :: Rational -> Int -> Rational -> Rational -> 
             Int -> Rational -> (Rational, Int)
scaleRat b minExp xMin xMax p x =
    if p <= minExp then
        (x, p)
    else if x >= xMax then
        scaleRat b minExp xMin xMax (p+1) (x/b)
    else if x < xMin  then
        scaleRat b minExp xMin xMax (p-1) (x*b)
    else
        (x, p)

-- Возведение в степень с помощью кэша наиболее частых чисел.
minExpt = 0::Int
maxExpt = 1100::Int
expt :: Integer -> Int -> Integer
expt base n =
    if base == 2 && n >= minExpt && n <= maxExpt then
        expts!n
    else
        base^n

expts :: Array Int Integer
expts = array (minExpt,maxExpt) [(n,2^n) | n <- [minExpt .. maxExpt]]

-- Вычисляет (нижнюю границу) log i по основанию b.
-- Наиболее простой способ --- просто делить i на b, пока оно не станет меньше b,
-- но это было бы очень медленно!  Мы просто немного более сообразительны.
integerLogBase :: Integer -> Integer -> Int
integerLogBase b i =
     if i < b then
        0
     else
        -- Пытается сначала возвести в квадрат основание, чтобы  сократить число делений.
        let l = 2 * integerLogBase (b*b) i
            doDiv :: Integer -> Int -> Int
            doDiv i l = if i < b then l else doDiv (i `div` b) (l+1)
        in  doDiv (i `div` (b^l)) l


-- Разные утилиты для отображения целых чисел и чисел с плавающей точкой 

showSigned :: Real a => (a -> ShowS) -> Int -> a -> ShowS
showSigned showPos p x 
  | x < 0     = showParen (p > 6) (showChar '-' . showPos (-x))
  | otherwise = showPos x

-- showInt, showOct, showHex используются только для положительных чисел
showInt, showOct, showHex :: Integral a => a -> ShowS
showOct = showIntAtBase  8 intToDigit
showInt = showIntAtBase 10 intToDigit
showHex = showIntAtBase 16 intToDigit

showIntAtBase :: Integral a 
      => a              -- основание
      -> (Int -> Char)  -- цифра для символа
      -> a              -- число для отображения
      -> ShowS
showIntAtBase base intToDig n rest
  | n < 0     = error "Numeric.showIntAtBase: не могу отображать отрицательные числа"
  | n' == 0   = rest'
  | otherwise = showIntAtBase base intToDig n' rest'
  where
    (n',d) = quotRem n base
    rest'  = intToDig (fromIntegral d) : rest


readSigned :: (Real a) => ReadS a -> ReadS a
readSigned readPos = readParen False read'
                     where read' r  = read'' r ++
                                      [(-x,t) | ("-",s) <- lex r,
                                                (x,t)   <- read'' s]
                           read'' r = [(n,s)  | (str,s) <- lex r,
                                                (n,"")  <- readPos str]


-- readInt считывает строку цифр, используя произвольное основание.  
-- Знаки минус в начале строки должны обрабатываться где-то в другом месте.

readInt :: (Integral a) => a -> (Char -> Bool) -> (Char -> Int) -> ReadS a
readInt radix isDig digToInt s =
   [(foldl1 (\n d -> n * radix + d) (map (fromIntegral . digToInt) ds), r)
          | (ds,r) <- nonnull isDig s ]

-- Функции для считывания беззнаковых чисел с различными основаниями
readDec, readOct, readHex :: (Integral a) => ReadS a
readDec = readInt 10 isDigit    digitToInt
readOct = readInt  8 isOctDigit digitToInt
readHex = readInt 16 isHexDigit digitToInt


showEFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showFFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showGFloat     :: (RealFloat a) => Maybe Int -> a -> ShowS
showFloat      :: (RealFloat a) => a -> ShowS

showEFloat d x =  showString (formatRealFloat FFExponent d x)
showFFloat d x =  showString (formatRealFloat FFFixed d x)
showGFloat d x =  showString (formatRealFloat FFGeneric d x)
showFloat      =  showGFloat Nothing 

-- Это типы форматов.  Этот тип не экспортируется.

data FFFormat = FFExponent | FFFixed | FFGeneric

formatRealFloat :: (RealFloat a) => FFFormat -> Maybe Int -> a -> String
formatRealFloat fmt decs x 
  = s
  where 
    base = 10
    s = if isNaN x then 
            "NaN"
        else if isInfinite x then 
            if x < 0 then "-Infinity" else "Infinity"
        else if x < 0 || isNegativeZero x then 
            '-' : doFmt fmt (floatToDigits (toInteger base) (-x))
        else 
            doFmt fmt (floatToDigits (toInteger base) x)
    
    doFmt fmt (is, e)
      = let 
           ds = map intToDigit is
        in  
        case fmt of
          FFGeneric -> 
              doFmt (if e < 0 || e > 7 then FFExponent else FFFixed)
                    (is, e)
          FFExponent ->
            case decs of
              Nothing ->
                case ds of
                   []    -> "0.0e0"
                   [d]   -> d : ".0e" ++ show (e-1)
                   d:ds  -> d : '.' : ds ++ 'e':show (e-1)
          
              Just dec ->
                let dec' = max dec 1 in
                case is of
                  [] -> '0':'.':take dec' (repeat '0') ++ "e0"
                  _ ->
                    let (ei, is') = roundTo base (dec'+1) is
                        d:ds = map intToDigit
                                   (if ei > 0 then init is' else is')
                    in d:'.':ds  ++ "e" ++ show (e-1+ei)
          
          FFFixed ->
            case decs of
               Nothing  -- Всегда печатает десятичную точку
                 | e > 0     -> take e (ds ++ repeat '0')
                                ++ '.' : mk0 (drop e ds)
                 | otherwise -> "0." ++ mk0 (replicate (-e) '0' ++ ds)
              
               Just dec ->  -- Печатает десятичную точку, если dec > 0
                 let dec' = max dec 0 in
                 if e >= 0 then
                   let (ei, is') = roundTo base (dec' + e) is
                       (ls, rs)  = splitAt (e+ei) 
                                              (map intToDigit is')
                   in  mk0 ls ++ mkdot0 rs
                 else
                   let (ei, is') = roundTo base dec' 
                                           (replicate (-e) 0 ++ is)
                       d : ds = map intToDigit 
                                    (if ei > 0 then is' else 0:is')
                   in  d : mkdot0 ds
            where   
              mk0 "" = "0"        -- Печатает 0.34, а не .34
              mk0 s  = s  
    
              mkdot0 "" = ""       -- Печатает 34, а не 34.
              mkdot0 s  = '.' : s  -- когда формат задает отсутствие
           -- цифр после десятичной точки
    

roundTo :: Int -> Int -> [Int] -> (Int, [Int])
roundTo base d is = case f d is of
                (0, is) -> (0, is)
                (1, is) -> (1, 1 : is)
  where b2 = base `div` 2
        f n [] = (0, replicate n 0)
        f 0 (i:_) = (if i >= b2 then 1 else 0, [])
        f d (i:is) = 
            let (c, ds) = f (d-1) is
                i' = c + i
            in  if i' == base then (1, 0:ds) else (0, i':ds)

--
-- Базируется на "Быстрой и точной печати чисел с плавающей точкой"
-- Р.Г. Бургера и Р.К. Дайбвига
-- ("Printing Floating-Point Numbers Quickly and Accurately"
-- R.G. Burger и R. K. Dybvig)
-- в PLDI 96.
-- Версия, приведенная здесь, использует намного более медленную оценку алгоритма.  
-- Ее следует усовершенствовать.

-- Эта функция возвращает непустой список цифр (целые числа в диапазоне [0..base-1])
-- и экспоненту.  В общем случае, если
--      floatToDigits r = ([a, b, ... z], e)
-- то
--      r = 0.ab..z * base^e
-- 

floatToDigits :: (RealFloat a) => Integer -> a -> ([Int], Int)

floatToDigits _ 0 = ([], 0)
floatToDigits base x =
    let (f0, e0) = decodeFloat x
        (minExp0, _) = floatRange x
        p = floatDigits x
        b = floatRadix x
        minExp = minExp0 - p            -- действительная минимальная экспонента

        -- В Haskell требуется, чтобы f было скорректировано так, чтобы  денормализационные числа
        -- имели невозможно низкую экспоненту.  Для этого используется коррекция.
        f :: Integer
        e :: Int
        (f, e) = let n = minExp - e0
                 in  if n > 0 then (f0 `div` (b^n), e0+n) else (f0, e0)

        (r, s, mUp, mDn) =
           if e >= 0 then
               let be = b^e in
               if f == b^(p-1) then
                   (f*be*b*2, 2*b, be*b, b)
               else
                   (f*be*2, 2, be, be)
           else
               if e > minExp && f == b^(p-1) then
                   (f*b*2, b^(-e+1)*2, b, 1)
               else
                   (f*2, b^(-e)*2, 1, 1)
        k = 
            let k0 =
                    if b==2 && base==10 then
                        -- logBase 10 2 немного больше, чем 3/10, поэтому
                        -- следующее вызовет ошибку на нижней стороне.  Игнорирование
                        -- дроби создаст эту ошибку даже больше.
                        -- Haskell обещает, что p-1 <= logBase b f < p.
                        (p - 1 + e0) * 3 `div` 10
                    else
                        ceiling ((log (fromInteger (f+1)) + 
                                 fromIntegral e * log (fromInteger b)) / 
                                  log (fromInteger base))
                fixup n =
                    if n >= 0 then
                        if r + mUp <= expt base n * s then n else fixup (n+1)
                    else
                        if expt base (-n) * (r + mUp) <= s then n
                                                           else fixup (n+1)
            in  fixup k0

        gen ds rn sN mUpN mDnN =
            let (dn, rn') = (rn * base) `divMod` sN
                mUpN' = mUpN * base
                mDnN' = mDnN * base
            in  case (rn' < mDnN', rn' + mUpN' > sN) of
                (True,  False) -> dn : ds
                (False, True)  -> dn+1 : ds
                (True,  True)  -> if rn' * 2 < sN then dn : ds else dn+1 : ds
                (False, False) -> gen (dn:ds) rn' sN mUpN' mDnN'
        rds =
            if k >= 0 then
                gen [] r (s * expt base k) mUp mDn
            else
                let bk = expt base (-k)
                in  gen [] (r * bk) s (mUp * bk) (mDn * bk)
    in  (map fromIntegral (reverse rds), k)



-- Эта функция для считывания чисел с плавающей точкой использует менее  ограничивающий синтаксис 
-- для чисел с плавающей точкой, чем лексический анализатор Haskell. `.' является необязательной.

readFloat     :: (RealFrac a) => ReadS a
readFloat r    = [(fromRational ((n%1)*10^^(k-d)),t) | (n,d,s) <- readFix r,
                                                       (k,t)   <- readExp s] ++
                 [ (0/0, t) | ("NaN",t)      <- lex r] ++
                 [ (1/0, t) | ("Infinity",t) <- lex r]
               where 
                 readFix r = [(read (ds++ds'), length ds', t)
                             | (ds,d) <- lexDigits r,
                               (ds',t) <- lexFrac d ]
               
                 lexFrac ('.':ds) = lexDigits ds
                 lexFrac s        = [("",s)]        
                 
                 readExp (e:s) | e `elem` "eE" = readExp' s
                 readExp s                     = [(0,s)]
                 
                 readExp' ('-':s) = [(-k,t) | (k,t) <- readDec s]
                 readExp' ('+':s) = readDec s
                 readExp' s       = readDec s

lexDigits        :: ReadS String 
lexDigits        =  nonnull isDigit

nonnull          :: (Char -> Bool) -> ReadS String
nonnull p s      =  [(cs,t) | (cs@(_:_),t) <- [span p s]]



Описание Haskell 98
наверх | назад | вперед | содержание | предметный указатель функций
Декабрь 2002