Import Upstream version 20180207
[hcoop/debian/mlton.git] / basis-library / real / real.sml
CommitLineData
7f918cf1
CE
1(* Copyright (C) 2011-2014,2017 Matthew Fluet.
2 * Copyright (C) 2003-2007 Henry Cejtin, Matthew Fluet, Suresh
3 * Jagannathan, and Stephen Weeks.
4 *
5 * MLton is released under a BSD-style license.
6 * See the file MLton-LICENSE for details.
7 *)
8
9functor Real (structure W: WORD_EXTRA
10 structure R:
11 sig
12 include PRE_REAL
13 val castToWord: real -> W.word
14 val castFromWord: W.word -> real
15 end): REAL_EXTRA =
16 struct
17 structure MLton = Primitive.MLton
18 structure Prim = R
19 local
20 open IEEEReal
21 in
22 datatype float_class = datatype float_class
23 datatype rounding_mode = datatype rounding_mode
24 end
25 infix 4 == != ?=
26 type real = R.real
27
28 local
29 open Prim
30 in
31 val realSize = Int32.toInt realSize
32 val exponentBias = Int32.toInt exponentBias
33 val precision = Int32.toInt precision
34 val radix = Int32.toInt radix
35 end
36
37 val signBits = Word.one
38 val exponentSignificandBits = Word.- (Word.fromInt realSize, signBits)
39 val significandBits = Word.- (Word.fromInt precision, Word.one)
40 val exponentBits = Word.- (exponentSignificandBits, significandBits)
41
42 local
43 val mkMask : Word.word -> W.word =
44 fn b => W.notb (W.<< (W.notb W.zero, b))
45 in
46 val signMask =
47 W.<< (mkMask signBits, exponentSignificandBits)
48 val exponentMask =
49 W.<< (mkMask exponentBits, significandBits)
50 val significandMask =
51 mkMask significandBits
52 end
53
54 val class : real -> float_class =
55 fn r =>
56 let
57 val w = R.castToWord r
58 in
59 if W.andb (w, exponentMask) = exponentMask
60 then if W.andb (w, significandMask) = W.zero
61 then IEEEReal.INF
62 else IEEEReal.NAN
63 else if W.andb (w, exponentMask) = W.zero
64 then if W.andb (w, significandMask) = W.zero
65 then IEEEReal.ZERO
66 else IEEEReal.SUBNORMAL
67 else IEEEReal.NORMAL
68 end
69
70 val toBits : real -> {sign: bool, exponent: W.word, significand: W.word} =
71 fn r =>
72 let
73 val w = R.castToWord r
74 val significand =
75 W.andb (w, significandMask)
76 val exponent =
77 W.>> (W.andb (w, exponentMask), significandBits)
78 val sign =
79 W.andb (w, signMask) = signMask
80 in
81 {sign = sign,
82 exponent = exponent,
83 significand = significand}
84 end
85
86 val fromBits : {sign: bool, exponent: W.word, significand: W.word} -> real =
87 fn {sign, exponent, significand} =>
88 let
89 val w =
90 W.orb (if sign then W.<< (W.one, exponentSignificandBits) else W.zero,
91 W.orb (W.andb (W.<< (exponent, significandBits), exponentMask),
92 W.andb (significand, significandMask)))
93 val r = R.castFromWord w
94 in
95 r
96 end
97
98 local
99 open Prim
100 in
101 val op *+ = op *+
102 val op *- = op *-
103 val op * = op *
104 val op + = op +
105 val op - = op -
106 val op / = op /
107 val op < = op <
108 val op <= = op <=
109 val op > = op >
110 val op >= = op >=
111 val ~ = ~
112 val abs = abs
113 end
114
115 local
116 fun 'a make {fromRealUnsafe: 'a -> real,
117 toRealUnsafe: real -> 'a,
118 other : {precision: Primitive.Int32.int}} =
119 if R.precision = #precision other
120 then (fn (_: rounding_mode) => fromRealUnsafe,
121 toRealUnsafe)
122 else (fn (m: rounding_mode) => fn r =>
123 IEEEReal.withRoundingMode (m, fn () => fromRealUnsafe r),
124 toRealUnsafe)
125 in
126 val (fromReal32,toReal32) =
127 make {fromRealUnsafe = R.fromReal32Unsafe,
128 toRealUnsafe = R.toReal32Unsafe,
129 other = {precision = Primitive.Real32.precision}}
130 val (fromReal64,toReal64) =
131 make {fromRealUnsafe = R.fromReal64Unsafe,
132 toRealUnsafe = R.toReal64Unsafe,
133 other = {precision = Primitive.Real64.precision}}
134 end
135 local
136 structure S =
137 LargeReal_ChooseRealN
138 (type 'a t = real -> 'a
139 val fReal32 = toReal32
140 val fReal64 = toReal64)
141 in
142 val toLarge = S.f
143 end
144 local
145 structure S =
146 LargeReal_ChooseRealN
147 (type 'a t = rounding_mode -> 'a -> real
148 val fReal32 = fromReal32
149 val fReal64 = fromReal64)
150 in
151 val fromLarge = S.f
152 end
153
154 val negInf = R.castFromWord (W.orb (signMask, exponentMask))
155 val negOne = R.castFromWord (W.orb (signMask, W.<< (W.fromInt exponentBias, significandBits)))
156 val negZero = R.castFromWord signMask
157 val zero = R.castFromWord W.zero
158 val minPos = R.castFromWord W.one
159 val minNormalPos = R.castFromWord (W.<< (W.one, significandBits))
160 val half = R.castFromWord (W.<< (W.- (W.fromInt exponentBias, W.one), significandBits))
161 val one = R.castFromWord (W.<< (W.fromInt exponentBias, significandBits))
162 val two = R.castFromWord (W.<< (W.+ (W.fromInt exponentBias, W.one), significandBits))
163 val maxFinite = R.castFromWord (W.- (exponentMask, W.one))
164 val posInf = R.castFromWord exponentMask
165
166 val nan = posInf + negInf
167 val posNan = R.castFromWord (W.andb (R.castToWord nan, W.notb signMask))
168 val negNan = R.castFromWord (W.orb (R.castToWord nan, signMask))
169
170
171 fun isFinite r =
172 abs r <= maxFinite
173
174 val op == = Prim.==
175
176 val op != = not o op ==
177
178 fun isNan r = r != r
179
180 fun isNormal r = class r = NORMAL
181
182 val op ?= =
183 if MLton.Codegen.isAMD64 orelse MLton.Codegen.isLLVM orelse MLton.Codegen.isX86
184 then R.?=
185 else
186 fn (x, y) =>
187 case (class x, class y) of
188 (NAN, _) => true
189 | (_, NAN) => true
190 | (ZERO, ZERO) => true
191 | _ => R.== (x, y)
192
193 fun min (x, y) =
194 if x <= y then x
195 else if x > y then y
196 else if isNan y then x
197 else y
198
199 fun max (x, y) =
200 if x >= y then x
201 else if x < y then y
202 else if isNan y then x
203 else y
204
205 fun sign (x: real): int =
206 if x > zero then 1
207 else if x < zero then ~1
208 else if x == zero then 0
209 else raise Domain
210
211 val signBit = #sign o toBits
212
213 fun sameSign (x, y) = signBit x = signBit y
214
215 fun copySign (x, y) =
216 if sameSign (x, y)
217 then x
218 else ~ x
219
220 local
221 structure I = IEEEReal
222 in
223 fun compareReal (x, y) =
224 if x < y then I.LESS
225 else if x > y then I.GREATER
226 else if x == y then I.EQUAL
227 else I.UNORDERED
228 end
229
230 local
231 structure I = IEEEReal
232 structure G = General
233 in
234 fun compare (x, y) =
235 case compareReal (x, y) of
236 I.EQUAL => G.EQUAL
237 | I.GREATER => G.GREATER
238 | I.LESS => G.LESS
239 | I.UNORDERED => raise IEEEReal.Unordered
240 end
241
242 fun unordered (x, y) = isNan x orelse isNan y
243
244 (* nextAfter for subnormal and normal values works by converting
245 * the real to a word of equivalent size and doing an increment
246 * or decrement on the word. Because of the way IEEE floating
247 * point numbers are represented, word {de,in}crement
248 * automatically does the right thing at the boundary between
249 * normals and denormals. Also, convienently,
250 * maxFinite+1 = posInf and minFinite-1 = negInf.
251 *)
252 val nextAfter: real * real -> real =
253 fn (r, t) =>
254 case (class r, class t) of
255 (NAN, _) => nan
256 | (_, NAN) => nan
257 | (INF, _) => r
258 | (ZERO, ZERO) => t (* want "t", not "r", to get the sign right *)
259 | (ZERO, _) => if t > zero then minPos else ~minPos
260 | _ =>
261 if r == t then
262 r
263 else if (r > t) = (r > zero) then
264 R.castFromWord (W.- (R.castToWord r, W.one))
265 else
266 R.castFromWord (W.+ (R.castToWord r, W.one))
267
268 local
269 val one = One.make (fn () => ref (0 : C_Int.t))
270 in
271 fun toManExp x =
272 case class x of
273 INF => {exp = 0, man = x}
274 | NAN => {exp = 0, man = nan}
275 | ZERO => {exp = 0, man = x}
276 | _ => One.use (one, fn r =>
277 let
278 val man = R.frexp (x, r)
279 in
280 {exp = C_Int.toInt (!r), man = man}
281 end)
282 end
283
284 fun fromManExp {exp, man} =
285 (R.ldexp (man, C_Int.fromInt exp))
286 handle Overflow =>
287 man * (if Int.< (exp, 0) then zero else posInf)
288
289 val fromManExp =
290 if MLton.Codegen.isX86
291 then fromManExp
292 else
293 fn {exp, man} =>
294 case class man of
295 INF => man
296 | NAN => man
297 | ZERO => man
298 | _ => fromManExp {exp = exp, man = man}
299
300 local
301 val one = One.make (fn () => ref zero)
302 in
303 fun split x =
304 case class x of
305 INF => {frac = if x > zero then zero else ~zero,
306 whole = x}
307 | NAN => {frac = nan, whole = nan}
308 | _ =>
309 let
310 val (frac, whole) =
311 One.use (one, fn int =>
312 (R.modf (x, int), !int))
313 (* Some platforms' C libraries don't get sign of
314 * zero right.
315 *)
316 fun fix y =
317 if class y = ZERO andalso not (sameSign (x, y))
318 then ~ y
319 else y
320 in
321 {frac = fix frac,
322 whole = fix whole}
323 end
324 end
325
326 val realMod = #frac o split
327
328 fun checkFloat x =
329 if isFinite x then x
330 else if isNan x then raise Div
331 else raise Overflow
332
333 val realCeil = R.realCeil
334 val realFloor = R.realFloor
335 val realTrunc = R.realTrunc
336
337 (* Unfortunately, libc round ties to zero instead of even values. *)
338 (* Fortunately, if any rounding mode is supported, it's TO_NEAREST. *)
339 val realRound = fn r => IEEEReal.withRoundingMode (TO_NEAREST, fn () => R.round r)
340
341 fun rem (x, y) =
342 (case class x of
343 INF => nan
344 | NAN => nan
345 | ZERO => zero
346 | _ => (case class y of
347 INF => x
348 | NAN => nan
349 | ZERO => nan
350 | _ => x - realTrunc (x/y) * y))
351
352 (* fromDecimal, scan, fromString: decimal -> binary conversions *)
353 fun strtor (str: NullString.t,
354 rounding_mode: IEEEReal.rounding_mode) =
355 let
356 val rounding : C_Int.int =
357 case rounding_mode of
358 TO_NEAREST => 1
359 | TO_NEGINF => 3
360 | TO_POSINF => 2
361 | TO_ZERO => 0
362 in
363 Prim.strtor (str, rounding)
364 end
365 exception Bad
366 fun fromDecimalWithRoundingMode
367 ({class, digits, exp, sign}: IEEEReal.decimal_approx,
368 rounding_mode: IEEEReal.rounding_mode) =
369 let
370 fun doit () =
371 let
372 val exp =
373 if Int.< (exp, 0)
374 then concat ["-", Int.toString (Int.~ exp)]
375 else Int.toString exp
376(*
377 val str = concat [if sign then "-" else "",
378 "0.", digits,
379 "E", exp, "\000"]
380*)
381 val n = Int.+ (if sign then 1 else 0,
382 Int.+ (4 (* "0." + "E" + "\000" *),
383 Int.+ (List.length digits,
384 String.size exp)))
385 val a = Array.alloc n
386 fun upd (i, c) = (Array.update (a, i, c); Int.+ (i, 1))
387 val i = 0
388 val i = if sign then upd (i, #"-") else i
389 val i = upd (i, #"0")
390 val i = upd (i, #".")
391 val i =
392 List.foldl
393 (fn (d, i) =>
394 if Int.< (d, 0) orelse Int.> (d, 9)
395 then raise Bad
396 else upd (i, Char.chr (Int.+ (d, Char.ord #"0"))))
397 i digits
398 val i = upd (i, #"E")
399 val i = CharVector.foldl (fn (c, i) => upd (i, c)) i exp
400 val _ = upd (i, #"\000")
401 val str = Vector.unsafeFromArray a
402 val x = strtor (NullString.fromString str, rounding_mode)
403 in
404 x
405 end
406 in
407 SOME (case class of
408 INF => if sign then negInf else posInf
409 | NAN => if sign then negNan else posNan
410 | NORMAL => doit ()
411 | SUBNORMAL => doit ()
412 | ZERO => if sign then negZero else zero)
413 handle Bad => NONE
414 end
415
416 fun fromDecimal da = fromDecimalWithRoundingMode (da, TO_NEAREST)
417
418 fun scan reader state =
419 case IEEEReal.scan reader state of
420 NONE => NONE
421 | SOME (da, state) =>
422 SOME (valOf (fromDecimalWithRoundingMode
423 (da, IEEEReal.getRoundingMode ())),
424 state)
425
426 val fromString = StringCvt.scanString scan
427
428 (* toDecimal, fmt, toString: binary -> decimal conversions. *)
429 datatype mode = Fix | Gen | Sci
430 local
431 val one = One.make (fn () => ref (0: C_Int.int))
432 in
433 fun gdtoa (x: real, mode: mode, ndig: int,
434 rounding_mode: IEEEReal.rounding_mode) =
435 let
436 val mode : C_Int.int =
437 case mode of
438 Fix => 3
439 | Gen => 0
440 | Sci => 2
441 val ndig : C_Int.int = C_Int.fromInt ndig
442 val rounding : C_Int.int =
443 case rounding_mode of
444 TO_NEAREST => 1
445 | TO_NEGINF => 3
446 | TO_POSINF => 2
447 | TO_ZERO => 0
448 in
449 One.use (one, fn decpt =>
450 (Prim.gdtoa (x, mode, ndig, rounding, decpt),
451 C_Int.toInt (!decpt)))
452 end
453 end
454
455 fun toDecimal (x: real): IEEEReal.decimal_approx =
456 case class x of
457 INF => {class = INF,
458 digits = [],
459 exp = 0,
460 sign = x < zero}
461 | NAN => {class = NAN,
462 digits = [],
463 exp = 0,
464 sign = signBit x}
465 | ZERO => {class = ZERO,
466 digits = [],
467 exp = 0,
468 sign = signBit x}
469 | c =>
470 let
471 val (cs, exp) = gdtoa (x, Gen, 0, TO_NEAREST)
472 fun loop (i, ac) =
473 if Int.< (i, 0)
474 then ac
475 else loop (Int.- (i, 1),
476 (Int.- (Char.ord (CUtil.C_String.sub (cs, i)),
477 Char.ord #"0"))
478 :: ac)
479 val digits = loop (Int.- (CUtil.C_String.length cs, 1), [])
480 in
481 {class = c,
482 digits = digits,
483 exp = exp,
484 sign = x < zero}
485 end
486
487 datatype realfmt = datatype StringCvt.realfmt
488
489 local
490 fun fix (sign: string, cs: CUtil.C_String.t, decpt: int, ndig: int): string =
491 let
492 val length = CUtil.C_String.length cs
493 in
494 if Int.< (decpt, 0)
495 then
496 concat [sign,
497 "0.",
498 String.new (Int.~ decpt, #"0"),
499 CUtil.C_String.toString cs,
500 String.new (Int.+ (Int.- (ndig, length),
501 decpt),
502 #"0")]
503 else
504 let
505 val whole =
506 if decpt = 0
507 then "0"
508 else
509 String.tabulate (decpt, fn i =>
510 if Int.< (i, length)
511 then CUtil.C_String.sub (cs, i)
512 else #"0")
513 in
514 if 0 = ndig
515 then concat [sign, whole]
516 else
517 let
518 val frac =
519 String.tabulate
520 (ndig, fn i =>
521 let
522 val j = Int.+ (i, decpt)
523 in
524 if Int.< (j, length)
525 then CUtil.C_String.sub (cs, j)
526 else #"0"
527 end)
528 in
529 concat [sign, whole, ".", frac]
530 end
531 end
532 end
533 fun sci (x: real, ndig: int): string =
534 let
535 val sign = if x < zero then "~" else ""
536 val (cs, decpt) =
537 gdtoa (x, Sci, Int.+ (1, ndig), IEEEReal.getRoundingMode ())
538 val length = CUtil.C_String.length cs
539 val whole = String.tabulate (1, fn _ => CUtil.C_String.sub (cs, 0))
540 val frac =
541 if 0 = ndig
542 then ""
543 else concat [".",
544 String.tabulate
545 (ndig, fn i =>
546 let
547 val j = Int.+ (i, 1)
548 in
549 if Int.< (j, length)
550 then CUtil.C_String.sub (cs, j)
551 else #"0"
552 end)]
553 val exp = Int.- (decpt, 1)
554 val exp =
555 let
556 val (exp, sign) =
557 if Int.< (exp, 0)
558 then (Int.~ exp, "~")
559 else (exp, "")
560 in
561 concat [sign, Int.toString exp]
562 end
563 in
564 concat [sign, whole, frac, "E", exp]
565 end
566 fun gen (x: real, n: int): string =
567 let
568 val (prefix, x) =
569 if x < zero
570 then ("~", ~ x)
571 else ("", x)
572 val ss = Substring.full (sci (x, Int.- (n, 1)))
573 fun isE c = c = #"E"
574 fun isZero c = c = #"0"
575 val expS =
576 Substring.string (Substring.taker (not o isE) ss)
577 val exp = valOf (Int.fromString expS)
578 val man =
579 String.translate
580 (fn #"." => "" | c => str c)
581 (Substring.string (Substring.dropr isZero
582 (Substring.takel (not o isE) ss)))
583 val manSize = String.size man
584 fun zeros i = CharVector.tabulate (i, fn _ => #"0")
585 fun dotAt i =
586 concat [String.substring (man, 0, i),
587 ".", String.extract (man, i, NONE)]
588 fun sci () = concat [prefix,
589 if manSize = 1 then man else dotAt 1,
590 "E", expS]
591 val op - = Int.-
592 val op + = Int.+
593 val ~ = Int.~
594 val op >= = Int.>=
595 in
596 if exp >= (if manSize = 1 then 3 else manSize + 3)
597 then sci ()
598 else if exp >= manSize - 1
599 then concat [prefix, man, zeros (exp - (manSize - 1))]
600 else if exp >= 0
601 then concat [prefix, dotAt (exp + 1)]
602 else if exp >= (if manSize = 1 then ~2 else ~3)
603 then concat [prefix, "0.", zeros (~exp - 1), man]
604 else sci ()
605 end
606 in
607 fun fmt spec =
608 let
609 val doit =
610 case spec of
611 EXACT => IEEEReal.toString o toDecimal
612 | FIX opt =>
613 let
614 val n =
615 case opt of
616 NONE => 6
617 | SOME n =>
618 if Primitive.Controls.safe andalso Int.< (n, 0)
619 then raise Size
620 else n
621 in
622 fn x =>
623 let
624 val sign = if x < zero then "~" else ""
625 val (cs, decpt) =
626 gdtoa (x, Fix, n, IEEEReal.getRoundingMode ())
627 in
628 fix (sign, cs, decpt, n)
629 end
630 end
631 | GEN opt =>
632 let
633 val n =
634 case opt of
635 NONE => 12
636 | SOME n =>
637 if Primitive.Controls.safe andalso Int.< (n, 1)
638 then raise Size
639 else n
640 in
641 fn x => gen (x, n)
642 end
643 | SCI opt =>
644 let
645 val n =
646 case opt of
647 NONE => 6
648 | SOME n =>
649 if Primitive.Controls.safe andalso Int.< (n, 0)
650 then raise Size
651 else n
652 in
653 fn x => sci (x, n)
654 end
655 in
656 fn x =>
657 case class x of
658 NAN => (* if signBit x then "~nan" else *) "nan"
659 | INF => if x > zero then "inf" else "~inf"
660 | _ => doit x
661 end
662 end
663
664 val toString = fmt (StringCvt.GEN NONE)
665
666 (* Not all devices support all rounding modes.
667 * However, every device has ceil/floor/round/trunc.
668 *)
669 fun safeConvert (m, cvt, x) =
670 case m of
671 TO_POSINF => cvt (realCeil x)
672 | TO_NEGINF => cvt (realFloor x)
673 | TO_NEAREST => cvt (realRound x)
674 | TO_ZERO => cvt (realTrunc x)
675
676 local
677 fun 'a make {fromIntUnsafe: 'a -> real,
678 toIntUnsafe: real -> 'a,
679 other : {maxInt': Word.word -> 'a,
680 minInt': 'a,
681 precision': int}} =
682 (fromIntUnsafe,
683 if Int.< (precision, #precision' other) then
684 let
685 val trim = Int.- (Int.- (#precision' other, precision), 1)
686 val maxInt' = (#maxInt' other) (Word.fromInt trim)
687 val minInt' = #minInt' other
688 val maxInt = fromIntUnsafe maxInt'
689 val minInt = fromIntUnsafe minInt'
690 in
691 fn (m: rounding_mode) => fn x =>
692 if minInt <= x then
693 if x <= maxInt then
694 safeConvert (m, toIntUnsafe, x)
695 else
696 raise Overflow
697 else
698 if x < minInt then
699 raise Overflow
700 else
701 raise Domain (* NaN *)
702 end
703 else
704 let
705 val maxInt' = (#maxInt' other) 0w0
706 val minInt' = #minInt' other
707 val maxInt = fromIntUnsafe maxInt'
708 val minInt = fromIntUnsafe minInt'
709 in
710 fn (m: rounding_mode) => fn x =>
711 if minInt <= x then
712 if x <= maxInt then
713 safeConvert (m, toIntUnsafe, x)
714 else
715 if x < maxInt + one then
716 (case m of
717 TO_NEGINF => maxInt'
718 | TO_POSINF => raise Overflow
719 | TO_ZERO => maxInt'
720 | TO_NEAREST =>
721 (* Depends on maxInt being odd. *)
722 if x - maxInt >= half then
723 raise Overflow
724 else
725 maxInt')
726 else
727 raise Overflow
728 else
729 if x < minInt then
730 if minInt - one < x then
731 (case m of
732 TO_NEGINF => raise Overflow
733 | TO_POSINF => minInt'
734 | TO_ZERO => minInt'
735 | TO_NEAREST =>
736 (* Depends on minInt being even. *)
737 if x - minInt < ~half then
738 raise Overflow
739 else
740 minInt')
741 else
742 raise Overflow
743 else
744 raise Domain (* NaN *)
745 end)
746 in
747 val (fromInt8,toInt8) =
748 make {fromIntUnsafe = R.fromInt8Unsafe,
749 toIntUnsafe = R.toInt8Unsafe,
750 other = {maxInt' = fn w => Int8.<< (Int8.>> (Int8.maxInt', w), w),
751 minInt' = Int8.minInt',
752 precision' = Int8.precision'}}
753 val (fromInt16,toInt16) =
754 make {fromIntUnsafe = R.fromInt16Unsafe,
755 toIntUnsafe = R.toInt16Unsafe,
756 other = {maxInt' = fn w => Int16.<< (Int16.>> (Int16.maxInt', w), w),
757 minInt' = Int16.minInt',
758 precision' = Int16.precision'}}
759 val (fromInt32,toInt32) =
760 make {fromIntUnsafe = R.fromInt32Unsafe,
761 toIntUnsafe = R.toInt32Unsafe,
762 other = {maxInt' = fn w => Int32.<< (Int32.>> (Int32.maxInt', w), w),
763 minInt' = Int32.minInt',
764 precision' = Int32.precision'}}
765 val (fromInt64,toInt64) =
766 make {fromIntUnsafe = R.fromInt64Unsafe,
767 toIntUnsafe = R.toInt64Unsafe,
768 other = {maxInt' = fn w => Int64.<< (Int64.>> (Int64.maxInt', w), w),
769 minInt' = Int64.minInt',
770 precision' = Int64.precision'}}
771 end
772
773 val fromIntInf: IntInf.int -> real =
774 fn i =>
775 let
776 val str =
777 if IntInf.< (i, 0)
778 then "-" ^ (IntInf.toString (IntInf.~ i))
779 else IntInf.toString i
780 val x = strtor (NullString.nullTerm str,
781 IEEEReal.getRoundingMode ())
782 in
783 x
784 end
785
786 val toIntInf: rounding_mode -> real -> LargeInt.int =
787 fn mode => fn x =>
788 case class x of
789 INF => raise Overflow
790 | NAN => raise Domain
791 | ZERO => (0 : LargeInt.int)
792 | _ =>
793 let
794 (* This round may turn x into an INF, so we need to check the
795 * class again.
796 *)
797 val x =
798 case mode of
799 TO_POSINF => realCeil x
800 | TO_NEGINF => realFloor x
801 | TO_NEAREST => realRound x
802 | TO_ZERO => realTrunc x
803 in
804 case class x of
805 INF => raise Overflow
806 | _ => valOf (IntInf.fromString (fmt (StringCvt.FIX (SOME 0)) x))
807 end
808
809 local
810 structure S =
811 Int_ChooseInt
812 (type 'a t = 'a -> real
813 val fInt8 = fromInt8
814 val fInt16 = fromInt16
815 val fInt32 = fromInt32
816 val fInt64 = fromInt64
817 val fIntInf = fromIntInf)
818 in
819 val fromInt = S.f
820 end
821 local
822 structure S =
823 LargeInt_ChooseInt
824 (type 'a t = 'a -> real
825 val fInt8 = fromInt8
826 val fInt16 = fromInt16
827 val fInt32 = fromInt32
828 val fInt64 = fromInt64
829 val fIntInf = fromIntInf)
830 in
831 val fromLargeInt = S.f
832 end
833 local
834 structure S =
835 Int_ChooseInt
836 (type 'a t = rounding_mode -> real -> 'a
837 val fInt8 = toInt8
838 val fInt16 = toInt16
839 val fInt32 = toInt32
840 val fInt64 = toInt64
841 val fIntInf = toIntInf)
842 in
843 val toInt = S.f
844 end
845 local
846 structure S =
847 LargeInt_ChooseInt
848 (type 'a t = rounding_mode -> real -> 'a
849 val fInt8 = toInt8
850 val fInt16 = toInt16
851 val fInt32 = toInt32
852 val fInt64 = toInt64
853 val fIntInf = toIntInf)
854 in
855 val toLargeInt = S.f
856 end
857
858 val floor = toInt TO_NEGINF
859 val ceil = toInt TO_POSINF
860 val trunc = toInt TO_ZERO
861 val round = toInt TO_NEAREST
862
863 local
864 fun 'a make {fromWordUnsafe: 'a -> real,
865 toWordUnsafe: real -> 'a,
866 other : {maxWord': Word.word -> 'a,
867 wordSize: int,
868 zeroWord: 'a}} =
869 (fromWordUnsafe,
870 if Int.<= (precision, #wordSize other)
871 then let
872 val trim = Int.- (#wordSize other, precision)
873 val maxWord' = (#maxWord' other) (Word.fromInt trim)
874 val maxWord = fromWordUnsafe maxWord'
875 val zeroWord = #zeroWord other
876 in
877 fn (m: rounding_mode) => fn x =>
878 case class x of
879 INF => raise Overflow
880 | NAN => raise Domain
881 | _ => if zero <= x
882 then if x <= maxWord
883 then safeConvert (m, toWordUnsafe, x)
884 else raise Overflow
885 else if x > ~one
886 then (case m of
887 TO_NEGINF => raise Overflow
888 | TO_POSINF => zeroWord
889 | TO_ZERO => zeroWord
890 | TO_NEAREST =>
891 if x < ~half
892 then raise Overflow
893 else zeroWord)
894 else raise Overflow
895 end
896 else let
897 val maxWord' = (#maxWord' other) 0w0
898 val maxWord = fromWordUnsafe maxWord'
899 val zeroWord = #zeroWord other
900 in
901 fn (m: rounding_mode) => fn x =>
902 case class x of
903 INF => raise Overflow
904 | NAN => raise Domain
905 | _ => if zero <= x
906 then if x <= maxWord
907 then safeConvert (m, toWordUnsafe, x)
908 else if x < maxWord + one
909 then (case m of
910 TO_NEGINF => maxWord'
911 | TO_POSINF => raise Overflow
912 | TO_ZERO => maxWord'
913 | TO_NEAREST =>
914 (* Depends on maxWord being odd. *)
915 if x - maxWord >= half
916 then raise Overflow
917 else maxWord')
918 else raise Overflow
919 else if x > ~one
920 then (case m of
921 TO_NEGINF => raise Overflow
922 | TO_POSINF => zeroWord
923 | TO_ZERO => zeroWord
924 | TO_NEAREST =>
925 if x < ~half
926 then raise Overflow
927 else zeroWord)
928 else raise Overflow
929 end)
930 in
931 val (fromWord8,toWord8) =
932 make {fromWordUnsafe = R.fromWord8Unsafe,
933 toWordUnsafe = R.toWord8Unsafe,
934 other = {maxWord' = fn w => Word8.<< (Word8.>> (Word8.maxWord', w), w),
935 wordSize = Word8.wordSize,
936 zeroWord = Word8.zero}}
937 val (fromWord16,toWord16) =
938 make {fromWordUnsafe = R.fromWord16Unsafe,
939 toWordUnsafe = R.toWord16Unsafe,
940 other = {maxWord' = fn w => Word16.<< (Word16.>> (Word16.maxWord', w), w),
941 wordSize = Word16.wordSize,
942 zeroWord = Word16.zero}}
943 val (fromWord32,toWord32) =
944 make {fromWordUnsafe = R.fromWord32Unsafe,
945 toWordUnsafe = R.toWord32Unsafe,
946 other = {maxWord' = fn w => Word32.<< (Word32.>> (Word32.maxWord', w), w),
947 wordSize = Word32.wordSize,
948 zeroWord = Word32.zero}}
949 val (fromWord64,toWord64) =
950 make {fromWordUnsafe = R.fromWord64Unsafe,
951 toWordUnsafe = R.toWord64Unsafe,
952 other = {maxWord' = fn w => Word64.<< (Word64.>> (Word64.maxWord', w), w),
953 wordSize = Word64.wordSize,
954 zeroWord = Word64.zero}}
955 end
956
957 local
958 structure S =
959 Word_ChooseWordN
960 (type 'a t = 'a -> real
961 val fWord8 = fromWord8
962 val fWord16 = fromWord16
963 val fWord32 = fromWord32
964 val fWord64 = fromWord64)
965 in
966 val fromWord = S.f
967 end
968 local
969 structure S =
970 LargeWord_ChooseWordN
971 (type 'a t = 'a -> real
972 val fWord8 = fromWord8
973 val fWord16 = fromWord16
974 val fWord32 = fromWord32
975 val fWord64 = fromWord64)
976 in
977 val fromLargeWord = S.f
978 end
979 local
980 structure S =
981 Word_ChooseWordN
982 (type 'a t = rounding_mode -> real -> 'a
983 val fWord8 = toWord8
984 val fWord16 = toWord16
985 val fWord32 = toWord32
986 val fWord64 = toWord64)
987 in
988 val toWord = S.f
989 end
990 local
991 structure S =
992 LargeWord_ChooseWordN
993 (type 'a t = rounding_mode -> real -> 'a
994 val fWord8 = toWord8
995 val fWord16 = toWord16
996 val fWord32 = toWord32
997 val fWord64 = toWord64)
998 in
999 val toLargeWord = S.f
1000 end
1001
1002 structure Math =
1003 struct
1004 open Prim.Math
1005
1006 (* Patch functions to handle out-of-range args. Many C math
1007 * libraries do not do what the SML Basis Spec requires.
1008 *)
1009
1010 local
1011 fun patch f x =
1012 if x < ~one orelse x > one
1013 then nan
1014 else f x
1015 in
1016 val acos = patch acos
1017 val asin = patch asin
1018 end
1019
1020 local
1021 fun patch f x = if x < zero then nan else f x
1022 in
1023 val ln = patch ln
1024 val log10 = patch log10
1025 end
1026
1027 (* The x86 doesn't get exp right on infs. *)
1028 val exp =
1029 if MLton.Codegen.isX86
1030 andalso let open MLton.Platform.Arch in host = X86 end
1031 then (fn x =>
1032 case class x of
1033 INF => if x > zero then posInf else zero
1034 | _ => exp x)
1035 else exp
1036
1037 (* The Cygwin math library doesn't get pow right on some exceptional
1038 * cases.
1039 *
1040 * The Linux math library doesn't get pow (x, y) right when x < 0
1041 * and y is large (but finite).
1042 *
1043 * So, we define a pow function that gives the correct result on
1044 * exceptional cases, and only calls the C pow with x > 0.
1045 *)
1046 fun isInt (x: real): bool = x == realFloor x
1047
1048 (* isEven x assumes isInt x. *)
1049 fun isEven (x: real): bool = isInt (x / two)
1050
1051 fun isOddInt x = isInt x andalso not (isEven x)
1052
1053 fun isNeg x = x < zero
1054
1055 fun pow (x, y) =
1056 case class y of
1057 INF =>
1058 if class x = NAN
1059 then nan
1060 else if x < negOne orelse x > one
1061 then if isNeg y then zero else posInf
1062 else if negOne < x andalso x < one
1063 then if isNeg y then posInf else zero
1064 else (* x = 1 orelse x = ~1 *)
1065 nan
1066 | NAN => nan
1067 | ZERO => one
1068 | _ =>
1069 (case class x of
1070 INF =>
1071 if isNeg x
1072 then if isNeg y
1073 then if isOddInt y
1074 then negZero
1075 else zero
1076 else if isOddInt y
1077 then negInf
1078 else posInf
1079 else (* x = posInf *)
1080 if isNeg y then zero else posInf
1081 | NAN => nan
1082 | ZERO =>
1083 if isNeg y
1084 then if isOddInt y
1085 then copySign (posInf, x)
1086 else posInf
1087 else if isOddInt y
1088 then x
1089 else zero
1090 | _ =>
1091 if isNeg x
1092 then if isInt y
1093 then if isEven y
1094 then Prim.Math.pow (~ x, y)
1095 else negOne * Prim.Math.pow (~ x, y)
1096 else nan
1097 else Prim.Math.pow (x, y))
1098
1099 fun cosh x =
1100 case class x of
1101 INF => x
1102 | ZERO => one
1103 | _ => R.Math.cosh x
1104
1105 fun sinh x =
1106 case class x of
1107 INF => x
1108 | ZERO => x
1109 | _ => R.Math.sinh x
1110
1111 fun tanh x =
1112 case class x of
1113 INF => if x > zero then one else negOne
1114 | ZERO => x
1115 | _ => R.Math.tanh x
1116 end
1117 end
1118
1119structure Real32 = Real (structure W = Word32
1120 structure R =
1121 struct
1122 open Primitive.Real32
1123 local open Primitive.PackReal32 in
1124 val castToWord = castToWord
1125 val castFromWord = castFromWord
1126 end
1127 end)
1128structure Real64 = Real (structure W = Word64
1129 structure R =
1130 struct
1131 open Primitive.Real64
1132 local open Primitive.PackReal64 in
1133 val castToWord = castToWord
1134 val castFromWord = castFromWord
1135 end
1136 end)