| Россия, Петерубрг, СПБ-ГПУ, 1998 |
Инспектор
Вы можете этот курс.
Опубликован: 19.09.2008 | Уровень: специалист | Доступ: платный
Лекция 15:
Числовые функции
- Базируется на "Быстрой и точной печати чисел с плавающей точкой" - Р.Г. Бургера и Р.К. Дайбвига - ("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]]