1(* 2 Title: Standard Basis Library: Real Signature and structure. 3 Author: David Matthews 4 Copyright David Matthews 2000, 2005, 2008, 2016-17 5 6 This library is free software; you can redistribute it and/or 7 modify it under the terms of the GNU Lesser General Public 8 License version 2.1 as published by the Free Software Foundation. 9 10 This library is distributed in the hope that it will be useful, 11 but WITHOUT ANY WARRANTY; without even the implied warranty of 12 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 13 Lesser General Public License for more details. 14 15 You should have received a copy of the GNU Lesser General Public 16 License along with this library; if not, write to the Free Software 17 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA 18*) 19 20structure Real: REAL = 21struct 22 open IEEEReal 23 val fromLargeInt: LargeInt.int -> real = Real.rtsCallFastI_F "PolyFloatArbitraryPrecision" 24 25 val fromInt: int -> real = 26 (* We have to select the appropriate conversion. This will be 27 reduced down to the appropriate function but has to be 28 type-correct whether int is arbitrary precision or fixed 29 precision. Hence the "o Large/FixedInt.fromInt". *) 30 if Bootstrap.intIsArbitraryPrecision 31 then fromLargeInt o LargeInt.fromInt 32 else Real.fromFixedInt o FixedInt.fromInt 33 34 (* These are needed because we don't yet have conversion from string 35 to real. They are filtered out by the signature. *) 36 val zero = fromInt 0 and one = fromInt 1 and four = fromInt 4 37 38 local 39 val doReal : int*real->real = RunCall.rtsCallFull2 "PolyRealGeneral" 40 in 41 fun callReal n x = doReal(n, x) 42 end 43 44 local 45 val doReal : int*(real*real)->real = RunCall.rtsCallFull2 "PolyRealGeneral" 46 in 47 fun callRealReal n p = doReal(n, p) 48 end 49 50 local 51 val doReal : int*real->bool = RunCall.rtsCallFull2 "PolyRealGeneral" 52 in 53 fun callRealToBool n x = doReal(n, x) 54 end 55 56 local 57 val doReal : int*real->int = RunCall.rtsCallFull2 "PolyRealGeneral" 58 in 59 fun callRealToInt n x = doReal(n, x) 60 end 61 62 type real = real (* Pick up from globals. *) 63 64 structure Math: MATH = 65 struct 66 type real = real (* Pick up from globals. *) 67 val sqrt = Real.rtsCallFastF_F "PolyRealSqrt" 68 and sin = Real.rtsCallFastF_F "PolyRealSin" 69 and cos = Real.rtsCallFastF_F "PolyRealCos" 70 and atan = Real.rtsCallFastF_F "PolyRealArctan" 71 and exp = Real.rtsCallFastF_F "PolyRealExp" 72 and ln = Real.rtsCallFastF_F "PolyRealLog" 73 and tan = Real.rtsCallFastF_F "PolyRealTan" 74 and asin = Real.rtsCallFastF_F "PolyRealArcSin" 75 and acos = Real.rtsCallFastF_F "PolyRealArcCos" 76 and log10 = Real.rtsCallFastF_F "PolyRealLog10" 77 and sinh = Real.rtsCallFastF_F "PolyRealSinh" 78 and cosh = Real.rtsCallFastF_F "PolyRealCosh" 79 and tanh = Real.rtsCallFastF_F "PolyRealTanh" 80 81 (* These have not yet been done. *) 82 val atan2 = callRealReal 3 83 val pow = callRealReal 4 84 85 (* Derived values. *) 86 val e = exp one 87 val pi = four * atan one 88 89 end; 90 91 92 infix 4 == != ?=; 93 94 val op == = Real.== 95 val op != : real * real -> bool = not o op == 96 97 val radix : int = callRealToInt 11 zero 98 val precision : int = callRealToInt 12 zero 99 val maxFinite : real = callReal 13 zero 100 val minNormalPos : real = callReal 14 zero 101 102 val posInf : real = one/zero; 103 val negInf : real = ~one/zero; 104 105 (* We only implement this sort of real. *) 106 fun toLarge (x: real) : (*LargeReal.*)real =x 107 fun fromLarge (_ : IEEEReal.rounding_mode) (x: (*LargeReal.*)real): real = x 108 109 local 110 open Real 111 in 112 (* NAN values fail any test including equality with themselves. *) 113 fun isNan x = x != x 114 (* NAN values do not match and infinities when multiplied by 0 produce NAN. *) 115 fun isFinite x = x * zero == zero 116 117 val signBit : real -> bool = callRealToBool 17 118 val copySign : (real * real) -> real = callRealReal 18 119 120 (* If we assume that all functions produce normalised results where 121 possible, the only subnormal values will be those smaller than 122 minNormalPos. *) 123 fun isNormal x = isFinite x andalso abs x >= minNormalPos 124 125 fun class x = 126 if isFinite x then if x == zero then ZERO 127 else if abs x >= minNormalPos then NORMAL 128 else SUBNORMAL 129 else if isNan x then NAN 130 else (* not finite and not Nan *) INF 131 132 fun sign x = 133 if isNan x then raise General.Domain 134 else if x == zero then 0 else if x < zero then ~1 else 1 135 end 136 137 fun sameSign (x, y) = signBit x = signBit y 138 139 fun unordered (x, y) = isNan x orelse isNan y 140 141 (* Returns the minimum. In the case where one is a NaN it returns the 142 other. In that case the comparison will be false. *) 143 fun min (a: real, b: real): real = if a < b orelse isNan b then a else b 144 (* Similarly for max. *) 145 fun max (a: real, b: real): real = if a > b orelse isNan b then a else b 146 147 fun checkFloat x = 148 if isFinite x then x 149 else if isNan x then raise General.Div else raise General.Overflow 150 151 val radixAsReal (* Not exported *) = fromInt radix 152 val epsilon (* Not exported *) = Math.pow(radixAsReal, fromInt (Int.-(1, precision))) 153 154 val minPos : real = minNormalPos*epsilon; 155 156 local 157 val toMantissa : real->real = callReal 24 158 and toExponent : real->int = callRealToInt 25 159 160 val doReal : int*(real*int)->real = RunCall.rtsCallFull2 "PolyRealGeneral" 161 162 fun fromManAndExp (ri: real*int): real = doReal(23, ri) 163 164 open Real 165 in 166 fun toManExp r = 167 if not (isFinite r) orelse r == zero 168 (* Nan, infinities and +/-0 all return r in the mantissa. 169 We include 0 to preserve its sign. *) 170 then {man=r, exp=0} 171 else {man=toMantissa r, exp=toExponent r} 172 173 fun fromManExp {man, exp} = 174 if not (isFinite man) orelse man == zero 175 (* Nan, infinities and +/-0 in the mantissa all return 176 their argument. *) 177 then man 178 else fromManAndExp(man, exp) 179 end 180 181 (* Convert to integer. Ideally this would do the rounding/truncation as part of the 182 conversion but it doesn't seem to be possible to detect overflow properly. 183 Instead we use the real rounding/truncation, convert to arbitrary 184 precision and then check for overflow if necessary. *) 185 local 186 (* The RTS function converts to at most a 64-bit value (even on 187 32-bits). That will convert all the bits of the mantissa 188 but if the exponent is large we may have to multiply by 189 some power of two. *) 190 val realToInt: real -> LargeInt.int = RunCall.rtsCallFull1 "PolyRealBoxedToLongInt" 191 in 192 val realFloor = Real.rtsCallFastF_F "PolyRealFloor" 193 and realCeil = Real.rtsCallFastF_F "PolyRealCeil" 194 and realTrunc = Real.rtsCallFastF_F "PolyRealTrunc" 195 and realRound = Real.rtsCallFastF_F "PolyRealRound" 196 197 fun toArbitrary x = 198 if isNan x then raise General.Domain 199 else if not (isFinite x) then raise General.Overflow 200 else 201 let 202 val { man, exp } = toManExp x 203 in 204 if exp <= precision 205 then realToInt x 206 else IntInf.<< (realToInt(fromManExp{man=man, exp=precision}), Word.fromInt(exp - precision)) 207 end 208 209 fun floor x = toArbitrary(realFloor x) 210 (* Returns the largest integer <= x. *) 211 212 fun ceil x = toArbitrary(realCeil x) 213 (* Returns the smallest integer >= x. *) 214 215 fun trunc x = toArbitrary(realTrunc x) 216 (* Truncate towards zero. *) 217 218 fun round x = toArbitrary(realRound x) 219 (* Return the nearest integer, returning an even value if equidistant. *) 220 221 fun toLargeInt IEEEReal.TO_NEGINF r = floor r 222 | toLargeInt IEEEReal.TO_POSINF r = ceil r 223 | toLargeInt IEEEReal.TO_ZERO r = trunc r 224 | toLargeInt IEEEReal.TO_NEAREST r = round r 225 226 fun toInt mode x = LargeInt.toInt(toLargeInt mode x) 227 228 val floor = LargeInt.toInt o floor 229 and ceil = LargeInt.toInt o ceil 230 and trunc = LargeInt.toInt o trunc 231 and round = LargeInt.toInt o round 232 end; 233 234 local 235 val realConv: string->real = RunCall.rtsCallFull1 "PolyRealBoxedFromString" 236 237 val posNan = abs(zero / zero) 238 val negNan = ~posNan 239 in 240 fun fromDecimal { class = INF, sign=true, ...} = SOME negInf 241 | fromDecimal { class = INF, sign=false, ...} = SOME posInf 242 | fromDecimal { class = ZERO, sign=true, ...} = SOME (~ zero) 243 | fromDecimal { class = ZERO, sign=false, ...} = SOME zero 244 (* Generate signed Nans ignoring the digits and mantissa. There 245 was code here to set the mantissa but there's no reference to 246 that in the current version of the Basis library. *) 247 | fromDecimal { class = NAN, sign=true, ... } = SOME negNan 248 | fromDecimal { class = NAN, sign=false, ... } = SOME posNan 249 250 | fromDecimal { class = _ (* NORMAL or SUBNORMAL *), sign, digits, exp} = 251 (let 252 fun toChar x = 253 if x < 0 orelse x > 9 then raise General.Domain 254 else Char.chr (x + Char.ord #"0") 255 (* Turn the number into a string. *) 256 val str = "0." ^ String.implode(List.map toChar digits) ^"E" ^ 257 Int.toString exp 258 (* Convert it to a real using the RTS conversion function. 259 Change any Conversion exceptions into Domain. *) 260 val result = realConv str handle RunCall.Conversion _ => raise General.Domain 261 in 262 if sign then SOME (~result) else SOME result 263 end 264 handle General.Domain => NONE 265 ) 266 end 267 268 local 269 val dtoa: real*int*int -> string*int*int = RunCall.rtsCallFull3 "PolyRealBoxedToString" 270 open StringCvt 271 272 fun addZeros n = 273 if n <= 0 then "" else "0" ^ addZeros (n-1) 274 275 fun fixFmt ndigs r = 276 if isNan r then "nan" 277 else if not (isFinite r) 278 then if r < zero then "~inf" else "inf" 279 else 280 let 281 (* Try to get ndigs past the decimal point. *) 282 val (str, exp, sign) = dtoa(r, 3, ndigs) 283 val strLen = String.size str 284 (* If the exponents is negative or zero we need to put a zero 285 before the decimal point. If the exponent is positive and 286 less than the number of digits we can take that 287 many characters off, otherwise we have to pad with zeros. *) 288 val numb = 289 if exp <= 0 290 then (* Exponent is zero or negative - all significant digits are 291 after the decimal point. Put in any zeros before 292 the significant digits, then the significant digits 293 and then any trailing zeros. *) 294 if ndigs = 0 then "0" 295 else "0." ^ addZeros(~exp) ^ str ^ addZeros(ndigs-strLen+exp) 296 else if strLen <= exp 297 then (* Exponent is not less than the length of the string - 298 all significant digits are before the decimal point. Add 299 any extra zeros before the decimal point then zeros after it. *) 300 str ^ addZeros(exp-strLen) ^ 301 (if ndigs = 0 then "" else "." ^ addZeros ndigs) 302 else (* Significant digits straddle the decimal point - insert the 303 decimal point and add any trailing zeros. *) 304 String.substring(str, 0, exp) ^ "." ^ 305 String.substring(str, exp, strLen-exp) ^ 306 addZeros(ndigs-strLen+exp) 307 in 308 if sign <> 0 then "~" ^ numb else numb 309 end 310 311 fun sciFmt ndigs r = 312 if isNan r then "nan" 313 else if not (isFinite r) 314 then if r < zero then "~inf" else "inf" 315 else 316 let 317 (* Try to get ndigs+1 digits. 1 before the decimal point and ndigs after. *) 318 val (str, exp, sign) = dtoa(r, 2, ndigs+1) 319 val strLen = String.size str 320 fun addZeros n = 321 if n <= 0 then "" else "0" ^ addZeros (n-1) 322 val numb = 323 if strLen = 0 324 then "0" ^ (if ndigs = 0 then "" else "." ^ addZeros ndigs) ^ "E0" 325 else 326 (if strLen = 1 327 then str ^ (if ndigs = 0 then "" else "." ^ addZeros ndigs) 328 else String.substring(str, 0, 1) ^ "." ^ 329 String.substring(str, 1, strLen-1) ^ addZeros (ndigs-strLen+1) 330 ) ^ "E" ^ Int.toString (exp-1) 331 in 332 if sign <> 0 then "~" ^ numb else numb 333 end 334 335 fun genFmt ndigs r = 336 if isNan r then "nan" 337 else if not (isFinite r) 338 then if r < zero then "~inf" else "inf" 339 else 340 let 341 (* Try to get ndigs digits. *) 342 val (str, exp, sign) = dtoa(r, 2, ndigs) 343 val strLen = String.size str 344 val numb = 345 (* Have to use scientific notation if exp > ndigs. Also use it 346 if the exponent is small (TODO: adjust this) *) 347 if exp > ndigs orelse exp < ~5 348 then (* Scientific format *) 349 (if strLen = 1 then str 350 else String.substring(str, 0, 1) ^ "." ^ 351 String.substring(str, 1, strLen-1) 352 ) ^ "E" ^ Int.toString (exp-1) 353 354 else (* Fixed format (N.B. no trailing zeros are added after the 355 decimal point apart from one if necessary) *) 356 if exp <= 0 357 then (* Exponent is zero or negative - all significant digits are 358 after the decimal point. Put in any zeros before 359 the significant digits, then the significant digits 360 and then any trailing zeros. *) 361 "0." ^ addZeros(~exp) ^ str 362 else if strLen <= exp 363 then (* Exponent is not less than the length of the string - 364 all significant digits are before the decimal point. Add 365 any extra zeros before the decimal point. Insert .0 at the 366 end to make it a valid real number. *) 367 str ^ addZeros(exp-strLen) ^ ".0" 368 else (* Significant digits straddle the decimal point - insert the 369 decimal point. *) 370 String.substring(str, 0, exp) ^ "." ^ 371 String.substring(str, exp, strLen-exp) 372 in 373 if sign <> 0 then "~" ^ numb else numb 374 end 375 376 fun strToDigitList str = 377 let 378 fun getDigs i l = 379 if i < 0 then l 380 else getDigs (i-1) 381 ((Char.ord(String.sub(str, i)) - Char.ord #"0") :: l) 382 in 383 getDigs (String.size str - 1) [] 384 end 385 in 386 fun toDecimal r = 387 let 388 val sign = signBit r 389 val kind = class r 390 in 391 case kind of 392 ZERO => { class = ZERO, sign = sign, digits=[], exp = 0 } 393 | INF => { class = INF, sign = sign, digits=[], exp = 0 } 394 | NAN => { class = NAN, sign = sign, digits=[], exp = 0 } 395 | _ => (* NORMAL/SUBNORMAL *) 396 let 397 val (str, exp, sign) = dtoa(r, 0, 0) 398 val digits = strToDigitList str 399 in 400 { class = kind, sign = sign <> 0, digits = digits, exp = exp } 401 end 402 end 403 404 (* Note: The definition says, reasonably, that negative values 405 for the number of digits raises Size. The tests also check 406 for a very large value for the number of digits and seem to 407 expect Size to be raised in that case. Note that the exception 408 is raised when fmt spec is evaluated and before it is applied 409 to an actual real argument. 410 In all cases, even EXACT format, this should produce "nan" for a NaN 411 and ignore the sign bit. *) 412 fun fmt (SCI NONE) = sciFmt 6 413 | fmt (SCI (SOME d) ) = 414 if d < 0 orelse d > 200 then raise General.Size 415 else sciFmt d 416 | fmt (FIX NONE) = fixFmt 6 417 | fmt (FIX (SOME d) ) = 418 if d < 0 orelse d > 200 then raise General.Size 419 else fixFmt d 420 | fmt (GEN NONE) = genFmt 12 421 | fmt (GEN (SOME d) ) = 422 if d < 1 orelse d > 200 then raise General.Size 423 else genFmt d 424 | fmt EXACT = (fn r => if isNan r then "nan" else IEEEReal.toString(toDecimal r)) 425 426 val toString = fmt (GEN NONE) 427 end 428 429 430 fun scan getc src = 431 let 432 (* Return a list of digits. *) 433 fun getdigits inp src = 434 case getc src of 435 NONE => (List.rev inp, src) 436 | SOME(ch, src') => 437 if ch >= #"0" andalso ch <= #"9" 438 then getdigits ((Char.ord ch - Char.ord #"0") :: inp) src' 439 else (List.rev inp, src) 440 441 (* Read an unsigned integer. Returns NONE if no digits have been read. *) 442 fun getNumber sign digits acc src = 443 case getc src of 444 NONE => if digits = 0 then NONE else SOME(if sign then ~acc else acc, src) 445 | SOME(ch, src') => 446 if ch >= #"0" andalso ch <= #"9" 447 then getNumber sign (digits+1) (acc*10 + Char.ord ch - Char.ord #"0") src' 448 else if digits = 0 then NONE else SOME(if sign then ~acc else acc, src') 449 450 (* Return the signed exponent. *) 451 fun getExponent src = 452 case getc src of 453 NONE => NONE 454 | SOME(ch, src') => 455 if ch = #"+" 456 then getNumber false 0 0 src' 457 else if ch = #"-" orelse ch = #"~" 458 then getNumber true 0 0 src' 459 else getNumber false 0 0 src 460 461 fun read_number sign src = 462 case getc src of 463 NONE => NONE 464 | SOME(ch, _) => 465 if not (ch >= #"0" andalso ch <= #"9" orelse ch = #".") 466 then NONE (* Bad *) 467 else (* Digits or decimal. *) 468 let 469 (* Get the digits before the decimal point (if any) *) 470 val (intPart, srcAfterDigs) = getdigits [] src 471 (* Get the digits after the decimal point (if any). 472 If there is a decimal point we only accept it if 473 there is at least one digit after it. *) 474 val (decimals, srcAfterMant) = 475 case getc srcAfterDigs of 476 NONE => ([], srcAfterDigs) 477 | SOME (#".", srcAfterDP) => 478 ( (* Check that the next character is a digit. *) 479 case getc srcAfterDP of 480 NONE => ([], srcAfterDigs) 481 | SOME(ch, _) => 482 if ch >= #"0" andalso ch <= #"9" 483 then getdigits [] srcAfterDP 484 else ([], srcAfterDigs) 485 ) 486 | SOME (_, _) => ([], srcAfterDigs) 487 (* The exponent is optional. If it doesn't form a valid 488 exponent we return zero as the value and the 489 continuation is the beginning of the "exponent". *) 490 val (exponent, srcAfterExp) = 491 case getc srcAfterMant of 492 NONE => (0, srcAfterMant) 493 | SOME (ch, src'''') => 494 if ch = #"e" orelse ch = #"E" 495 then 496 ( 497 case getExponent src'''' of 498 NONE => (0, srcAfterMant) 499 | SOME x => x 500 ) 501 else (0, srcAfterMant) 502 (* Generate a decimal representation ready for conversion. 503 We don't bother to strip off leading or trailing zeros. *) 504 val decimalRep = {class=NORMAL, sign=sign, 505 digits=List.@(intPart, decimals), 506 exp=exponent + List.length intPart} 507 in 508 case fromDecimal decimalRep of 509 SOME r => SOME(r, srcAfterExp) 510 | NONE => NONE 511 end 512 in 513 case getc src of 514 NONE => NONE 515 | SOME(ch, src') => 516 if Char.isSpace ch (* Skip white space. *) 517 then scan getc src' (* Recurse *) 518 else if ch = #"+" (* Remove the + sign *) 519 then read_number false src' 520 else if ch = #"-" orelse ch = #"~" 521 then read_number true src' 522 else (* See if it's a valid digit. *) 523 read_number false src 524 end 525 526 val fromString = StringCvt.scanString scan 527 528 (* Converter to real values. This replaces the basic conversion 529 function for reals installed in the bootstrap process. 530 For more information see convInt in Int. *) 531 local 532 fun convReal (s: string) : real = 533 let 534 (* Set the rounding mode to TO_NEAREST whatever the current 535 rounding mode. Otherwise the result of compiling a piece of 536 code with a literal constant could depend on what the rounding 537 mode was set to. *) 538 val oldRounding = IEEEReal.getRoundingMode() handle Fail _ => IEEEReal.TO_NEAREST 539 val () = IEEEReal.setRoundingMode IEEEReal.TO_NEAREST handle Fail _ => () 540 val scanResult = StringCvt.scanString scan s 541 val () = IEEEReal.setRoundingMode oldRounding handle Fail _ => () 542 in 543 case scanResult of 544 NONE => raise RunCall.Conversion "Invalid real constant" 545 | SOME res => res 546 end 547 in 548 (* Install this as a conversion function for real literals. *) 549 val (): unit = RunCall.addOverload convReal "convReal" 550 end 551 552 open Real (* Get the other definitions. *) 553 554 fun compare (r1, r2) = 555 if r1 == r2 then General.EQUAL 556 else if r1 < r2 then General.LESS 557 else if r1 > r2 then General.GREATER 558 else raise Unordered 559 560 fun compareReal (r1, r2) = 561 if r1 == r2 then EQUAL 562 else if r1 < r2 then LESS 563 else if r1 > r2 then GREATER 564 else UNORDERED 565 566 (* Question: The definition says "bitwise equal, ignoring signs on zeros". 567 If we assume that all numbers are normalised, is that the same as "equal"?*) 568 fun op ?= (x, y) = 569 isNan x orelse isNan y orelse x == y 570 571 (* Although these may be built in in some architectures it's 572 probably not worth treating them specially at the moment. *) 573 fun *+ (x: real, y: real, z: real): real = x*y+z 574 and *- (x: real, y: real, z: real): real = x*y-z 575 576 fun rem (x, y) = 577 if not (isFinite y) andalso not (isNan y) then x 578 else x - realTrunc(x / y)*y 579 580 (* Split a real into whole and fractional parts. The fractional part must have 581 the same sign as the number even if it is zero. *) 582 fun split r = 583 let 584 val whole = realTrunc r 585 val frac = r - whole 586 in 587 { whole = whole, 588 frac = 589 if not (isFinite r) 590 then if isNan r then r else (* Infinity *) if r < zero then ~zero else zero 591 else if frac == zero then if signBit r then ~zero else zero 592 else frac } 593 end 594 595 (* Get the fractional part of a real. *) 596 fun realMod r = #frac(split r) 597 598 (* nextAfter: This was previously implemented in ML but, at the very least, 599 needed to work with rounding to something other than TO_NEAREST. This should 600 be implemented as a fast call but we don't currently support fast calls for 601 real * real -> real. *) 602 val nextAfter = callRealReal 26 603 604end; 605 606structure Math = Real.Math; 607 608structure LargeReal: REAL = Real; 609 610(* Values available unqualified at the top-level. *) 611val real : int -> real = Real.fromInt 612val trunc : real -> int = Real.trunc 613val floor : real -> int = Real.floor 614val ceil : real -> int = Real.ceil 615val round : real -> int =Real.round; 616 617(* Install print function. *) 618local 619 fun print_real _ _ (r: real) = 620 PolyML.PrettyString(Real.fmt (StringCvt.GEN(SOME 10)) r) 621in 622 val () = PolyML.addPrettyPrinter print_real; 623end; 624