Update Trac ticket URLs to point to GitLab
[ghc.git] / libraries / base / GHC / Real.hs
1 {-# LANGUAGE Trustworthy #-}
2 {-# LANGUAGE CPP, NoImplicitPrelude, MagicHash, UnboxedTuples, BangPatterns #-}
3 {-# OPTIONS_GHC -Wno-orphans #-}
4 {-# OPTIONS_HADDOCK not-home #-}
5
6 -----------------------------------------------------------------------------
7 -- |
8 -- Module : GHC.Real
9 -- Copyright : (c) The University of Glasgow, 1994-2002
10 -- License : see libraries/base/LICENSE
11 --
12 -- Maintainer : cvs-ghc@haskell.org
13 -- Stability : internal
14 -- Portability : non-portable (GHC Extensions)
15 --
16 -- The types 'Ratio' and 'Rational', and the classes 'Real', 'Fractional',
17 -- 'Integral', and 'RealFrac'.
18 --
19 -----------------------------------------------------------------------------
20
21 module GHC.Real where
22
23 #include "MachDeps.h"
24
25 import GHC.Base
26 import GHC.Num
27 import GHC.List
28 import GHC.Enum
29 import GHC.Show
30 import {-# SOURCE #-} GHC.Exception( divZeroException, overflowException
31 , underflowException
32 , ratioZeroDenomException )
33
34 #if defined(MIN_VERSION_integer_gmp)
35 import GHC.Integer.GMP.Internals
36 #endif
37
38 infixr 8 ^, ^^
39 infixl 7 /, `quot`, `rem`, `div`, `mod`
40 infixl 7 %
41
42 default () -- Double isn't available yet,
43 -- and we shouldn't be using defaults anyway
44
45 ------------------------------------------------------------------------
46 -- Divide by zero and arithmetic overflow
47 ------------------------------------------------------------------------
48
49 -- We put them here because they are needed relatively early
50 -- in the libraries before the Exception type has been defined yet.
51
52 {-# NOINLINE divZeroError #-}
53 divZeroError :: a
54 divZeroError = raise# divZeroException
55
56 {-# NOINLINE ratioZeroDenominatorError #-}
57 ratioZeroDenominatorError :: a
58 ratioZeroDenominatorError = raise# ratioZeroDenomException
59
60 {-# NOINLINE overflowError #-}
61 overflowError :: a
62 overflowError = raise# overflowException
63
64 {-# NOINLINE underflowError #-}
65 underflowError :: a
66 underflowError = raise# underflowException
67
68
69 --------------------------------------------------------------
70 -- The Ratio and Rational types
71 --------------------------------------------------------------
72
73 -- | Rational numbers, with numerator and denominator of some 'Integral' type.
74 --
75 -- Note that `Ratio`'s instances inherit the deficiencies from the type
76 -- parameter's. For example, @Ratio Natural@'s 'Num' instance has similar
77 -- problems to `Numeric.Natural.Natural`'s.
78 data Ratio a = !a :% !a deriving Eq -- ^ @since 2.01
79
80 -- | Arbitrary-precision rational numbers, represented as a ratio of
81 -- two 'Integer' values. A rational number may be constructed using
82 -- the '%' operator.
83 type Rational = Ratio Integer
84
85 ratioPrec, ratioPrec1 :: Int
86 ratioPrec = 7 -- Precedence of ':%' constructor
87 ratioPrec1 = ratioPrec + 1
88
89 infinity, notANumber :: Rational
90 infinity = 1 :% 0
91 notANumber = 0 :% 0
92
93 -- Use :%, not % for Inf/NaN; the latter would
94 -- immediately lead to a runtime error, because it normalises.
95
96 -- | Forms the ratio of two integral numbers.
97 {-# SPECIALISE (%) :: Integer -> Integer -> Rational #-}
98 (%) :: (Integral a) => a -> a -> Ratio a
99
100 -- | Extract the numerator of the ratio in reduced form:
101 -- the numerator and denominator have no common factor and the denominator
102 -- is positive.
103 numerator :: Ratio a -> a
104
105 -- | Extract the denominator of the ratio in reduced form:
106 -- the numerator and denominator have no common factor and the denominator
107 -- is positive.
108 denominator :: Ratio a -> a
109
110
111 -- | 'reduce' is a subsidiary function used only in this module.
112 -- It normalises a ratio by dividing both numerator and denominator by
113 -- their greatest common divisor.
114 reduce :: (Integral a) => a -> a -> Ratio a
115 {-# SPECIALISE reduce :: Integer -> Integer -> Rational #-}
116 reduce _ 0 = ratioZeroDenominatorError
117 reduce x y = (x `quot` d) :% (y `quot` d)
118 where d = gcd x y
119
120 x % y = reduce (x * signum y) (abs y)
121
122 numerator (x :% _) = x
123 denominator (_ :% y) = y
124
125 --------------------------------------------------------------
126 -- Standard numeric classes
127 --------------------------------------------------------------
128
129 class (Num a, Ord a) => Real a where
130 -- | the rational equivalent of its real argument with full precision
131 toRational :: a -> Rational
132
133 -- | Integral numbers, supporting integer division.
134 --
135 -- The Haskell Report defines no laws for 'Integral'. However, 'Integral'
136 -- instances are customarily expected to define a Euclidean domain and have the
137 -- following properties for the `div`\/`mod` and `quot`\/`rem` pairs, given
138 -- suitable Euclidean functions @f@ and @g@:
139 --
140 -- * @x@ = @y * quot x y + rem x y@ with @rem x y@ = @fromInteger 0@ or
141 -- @g (rem x y)@ < @g y@
142 -- * @x@ = @y * div x y + mod x y@ with @mod x y@ = @fromInteger 0@ or
143 -- @f (mod x y)@ < @f y@
144 --
145 -- An example of a suitable Euclidean function, for `Integer`'s instance, is
146 -- 'abs'.
147 class (Real a, Enum a) => Integral a where
148 -- | integer division truncated toward zero
149 quot :: a -> a -> a
150 -- | integer remainder, satisfying
151 --
152 -- > (x `quot` y)*y + (x `rem` y) == x
153 rem :: a -> a -> a
154 -- | integer division truncated toward negative infinity
155 div :: a -> a -> a
156 -- | integer modulus, satisfying
157 --
158 -- > (x `div` y)*y + (x `mod` y) == x
159 mod :: a -> a -> a
160 -- | simultaneous 'quot' and 'rem'
161 quotRem :: a -> a -> (a,a)
162 -- | simultaneous 'div' and 'mod'
163 divMod :: a -> a -> (a,a)
164 -- | conversion to 'Integer'
165 toInteger :: a -> Integer
166
167 {-# INLINE quot #-}
168 {-# INLINE rem #-}
169 {-# INLINE div #-}
170 {-# INLINE mod #-}
171 n `quot` d = q where (q,_) = quotRem n d
172 n `rem` d = r where (_,r) = quotRem n d
173 n `div` d = q where (q,_) = divMod n d
174 n `mod` d = r where (_,r) = divMod n d
175
176 divMod n d = if signum r == negate (signum d) then (q-1, r+d) else qr
177 where qr@(q,r) = quotRem n d
178
179 -- | Fractional numbers, supporting real division.
180 --
181 -- The Haskell Report defines no laws for 'Fractional'. However, @('+')@ and
182 -- @('*')@ are customarily expected to define a division ring and have the
183 -- following properties:
184 --
185 -- [__'recip' gives the multiplicative inverse__]:
186 -- @x * recip x@ = @recip x * x@ = @fromInteger 1@
187 --
188 -- Note that it /isn't/ customarily expected that a type instance of
189 -- 'Fractional' implement a field. However, all instances in @base@ do.
190 class (Num a) => Fractional a where
191 {-# MINIMAL fromRational, (recip | (/)) #-}
192
193 -- | Fractional division.
194 (/) :: a -> a -> a
195 -- | Reciprocal fraction.
196 recip :: a -> a
197 -- | Conversion from a 'Rational' (that is @'Ratio' 'Integer'@).
198 -- A floating literal stands for an application of 'fromRational'
199 -- to a value of type 'Rational', so such literals have type
200 -- @('Fractional' a) => a@.
201 fromRational :: Rational -> a
202
203 {-# INLINE recip #-}
204 {-# INLINE (/) #-}
205 recip x = 1 / x
206 x / y = x * recip y
207
208 -- | Extracting components of fractions.
209 class (Real a, Fractional a) => RealFrac a where
210 -- | The function 'properFraction' takes a real fractional number @x@
211 -- and returns a pair @(n,f)@ such that @x = n+f@, and:
212 --
213 -- * @n@ is an integral number with the same sign as @x@; and
214 --
215 -- * @f@ is a fraction with the same type and sign as @x@,
216 -- and with absolute value less than @1@.
217 --
218 -- The default definitions of the 'ceiling', 'floor', 'truncate'
219 -- and 'round' functions are in terms of 'properFraction'.
220 properFraction :: (Integral b) => a -> (b,a)
221 -- | @'truncate' x@ returns the integer nearest @x@ between zero and @x@
222 truncate :: (Integral b) => a -> b
223 -- | @'round' x@ returns the nearest integer to @x@;
224 -- the even integer if @x@ is equidistant between two integers
225 round :: (Integral b) => a -> b
226 -- | @'ceiling' x@ returns the least integer not less than @x@
227 ceiling :: (Integral b) => a -> b
228 -- | @'floor' x@ returns the greatest integer not greater than @x@
229 floor :: (Integral b) => a -> b
230
231 {-# INLINE truncate #-}
232 truncate x = m where (m,_) = properFraction x
233
234 round x = let (n,r) = properFraction x
235 m = if r < 0 then n - 1 else n + 1
236 in case signum (abs r - 0.5) of
237 -1 -> n
238 0 -> if even n then n else m
239 1 -> m
240 _ -> errorWithoutStackTrace "round default defn: Bad value"
241
242 ceiling x = if r > 0 then n + 1 else n
243 where (n,r) = properFraction x
244
245 floor x = if r < 0 then n - 1 else n
246 where (n,r) = properFraction x
247
248 -- These 'numeric' enumerations come straight from the Report
249
250 numericEnumFrom :: (Fractional a) => a -> [a]
251 numericEnumFrom n = go 0
252 where
253 -- See Note [Numeric Stability of Enumerating Floating Numbers]
254 go !k = let !n' = n + k
255 in n' : go (k + 1)
256
257 numericEnumFromThen :: (Fractional a) => a -> a -> [a]
258 numericEnumFromThen n m = go 0
259 where
260 step = m - n
261 -- See Note [Numeric Stability of Enumerating Floating Numbers]
262 go !k = let !n' = n + k * step
263 in n' : go (k + 1)
264
265 numericEnumFromTo :: (Ord a, Fractional a) => a -> a -> [a]
266 numericEnumFromTo n m = takeWhile (<= m + 1/2) (numericEnumFrom n)
267
268 numericEnumFromThenTo :: (Ord a, Fractional a) => a -> a -> a -> [a]
269 numericEnumFromThenTo e1 e2 e3
270 = takeWhile predicate (numericEnumFromThen e1 e2)
271 where
272 mid = (e2 - e1) / 2
273 predicate | e2 >= e1 = (<= e3 + mid)
274 | otherwise = (>= e3 + mid)
275
276 {- Note [Numeric Stability of Enumerating Floating Numbers]
277 -----------------------------------------------------------
278 When enumerate floating numbers, we could add the increment to the last number
279 at every run (as what we did previously):
280
281 numericEnumFrom n = n `seq` (n : numericEnumFrom (n + 1))
282
283 This approach is concise and really fast, only needs an addition operation.
284 However when a floating number is large enough, for `n`, `n` and `n+1` will
285 have the same binary representation. For example (all number has type
286 `Double`):
287
288 9007199254740990 is: 0x433ffffffffffffe
289 9007199254740990 + 1 is: 0x433fffffffffffff
290 (9007199254740990 + 1) + 1 is: 0x4340000000000000
291 ((9007199254740990 + 1) + 1) + 1 is: 0x4340000000000000
292
293 When we evaluate ([9007199254740990..9007199254740991] :: Double), we would
294 never reach the condition in `numericEnumFromTo`
295
296 9007199254740990 + 1 + 1 + ... > 9007199254740991 + 1/2
297
298 We would fall into infinite loop (as reported in #15081).
299
300 To remedy the situation, we record the number of `1` that needed to be added
301 to the start number, rather than increasing `1` at every time. This approach
302 can improvement the numeric stability greatly at the cost of a multiplication.
303
304 Furthermore, we use the type of the enumerated number, `Fractional a => a`,
305 as the type of multiplier. In rare situations, the multiplier could be very
306 large and will lead to the enumeration to infinite loop, too, which should
307 be very rare. Consider the following example:
308
309 [1..9007199254740994]
310
311 We could fix that by using an Integer as multiplier but we don't do that.
312 The benchmark on T7954.hs shows that this approach leads to significant
313 degeneration on performance (33% increase allocation and 300% increase on
314 elapsed time).
315
316 See #15081 and Phab:D4650 for the related discussion about this problem.
317 -}
318
319 --------------------------------------------------------------
320 -- Instances for Int
321 --------------------------------------------------------------
322
323 -- | @since 2.0.1
324 instance Real Int where
325 toRational x = toInteger x :% 1
326
327 -- | @since 2.0.1
328 instance Integral Int where
329 toInteger (I# i) = smallInteger i
330
331 a `quot` b
332 | b == 0 = divZeroError
333 | b == (-1) && a == minBound = overflowError -- Note [Order of tests]
334 -- in GHC.Int
335 | otherwise = a `quotInt` b
336
337 a `rem` b
338 | b == 0 = divZeroError
339 -- The quotRem CPU instruction fails for minBound `quotRem` -1,
340 -- but minBound `rem` -1 is well-defined (0). We therefore
341 -- special-case it.
342 | b == (-1) = 0
343 | otherwise = a `remInt` b
344
345 a `div` b
346 | b == 0 = divZeroError
347 | b == (-1) && a == minBound = overflowError -- Note [Order of tests]
348 -- in GHC.Int
349 | otherwise = a `divInt` b
350
351 a `mod` b
352 | b == 0 = divZeroError
353 -- The divMod CPU instruction fails for minBound `divMod` -1,
354 -- but minBound `mod` -1 is well-defined (0). We therefore
355 -- special-case it.
356 | b == (-1) = 0
357 | otherwise = a `modInt` b
358
359 a `quotRem` b
360 | b == 0 = divZeroError
361 -- Note [Order of tests] in GHC.Int
362 | b == (-1) && a == minBound = (overflowError, 0)
363 | otherwise = a `quotRemInt` b
364
365 a `divMod` b
366 | b == 0 = divZeroError
367 -- Note [Order of tests] in GHC.Int
368 | b == (-1) && a == minBound = (overflowError, 0)
369 | otherwise = a `divModInt` b
370
371 --------------------------------------------------------------
372 -- Instances for @Word@
373 --------------------------------------------------------------
374
375 -- | @since 2.01
376 instance Real Word where
377 toRational x = toInteger x % 1
378
379 -- | @since 2.01
380 instance Integral Word where
381 quot (W# x#) y@(W# y#)
382 | y /= 0 = W# (x# `quotWord#` y#)
383 | otherwise = divZeroError
384 rem (W# x#) y@(W# y#)
385 | y /= 0 = W# (x# `remWord#` y#)
386 | otherwise = divZeroError
387 div (W# x#) y@(W# y#)
388 | y /= 0 = W# (x# `quotWord#` y#)
389 | otherwise = divZeroError
390 mod (W# x#) y@(W# y#)
391 | y /= 0 = W# (x# `remWord#` y#)
392 | otherwise = divZeroError
393 quotRem (W# x#) y@(W# y#)
394 | y /= 0 = case x# `quotRemWord#` y# of
395 (# q, r #) ->
396 (W# q, W# r)
397 | otherwise = divZeroError
398 divMod (W# x#) y@(W# y#)
399 | y /= 0 = (W# (x# `quotWord#` y#), W# (x# `remWord#` y#))
400 | otherwise = divZeroError
401 toInteger (W# x#) = wordToInteger x#
402
403 --------------------------------------------------------------
404 -- Instances for Integer
405 --------------------------------------------------------------
406
407 -- | @since 2.0.1
408 instance Real Integer where
409 toRational x = x :% 1
410
411 -- | @since 4.8.0.0
412 instance Real Natural where
413 toRational n = naturalToInteger n :% 1
414
415 -- Note [Integer division constant folding]
416 -- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
417 --
418 -- Constant folding of quot, rem, div, mod, divMod and quotRem for
419 -- Integer arguments depends crucially on inlining. Constant folding
420 -- rules defined in compiler/prelude/PrelRules.hs trigger for
421 -- quotInteger, remInteger and so on. So if calls to quot, rem and so on
422 -- were not inlined the rules would not fire. The rules would also not
423 -- fire if calls to quotInteger and so on were inlined, but this does not
424 -- happen because they are all marked with NOINLINE pragma - see documentation
425 -- of integer-gmp or integer-simple.
426
427 -- | @since 2.0.1
428 instance Integral Integer where
429 toInteger n = n
430
431 {-# INLINE quot #-}
432 _ `quot` 0 = divZeroError
433 n `quot` d = n `quotInteger` d
434
435 {-# INLINE rem #-}
436 _ `rem` 0 = divZeroError
437 n `rem` d = n `remInteger` d
438
439 {-# INLINE div #-}
440 _ `div` 0 = divZeroError
441 n `div` d = n `divInteger` d
442
443 {-# INLINE mod #-}
444 _ `mod` 0 = divZeroError
445 n `mod` d = n `modInteger` d
446
447 {-# INLINE divMod #-}
448 _ `divMod` 0 = divZeroError
449 n `divMod` d = case n `divModInteger` d of
450 (# x, y #) -> (x, y)
451
452 {-# INLINE quotRem #-}
453 _ `quotRem` 0 = divZeroError
454 n `quotRem` d = case n `quotRemInteger` d of
455 (# q, r #) -> (q, r)
456
457 -- | @since 4.8.0.0
458 instance Integral Natural where
459 toInteger = naturalToInteger
460
461 divMod = quotRemNatural
462 div = quotNatural
463 mod = remNatural
464
465 quotRem = quotRemNatural
466 quot = quotNatural
467 rem = remNatural
468
469 --------------------------------------------------------------
470 -- Instances for @Ratio@
471 --------------------------------------------------------------
472
473 -- | @since 2.0.1
474 instance (Integral a) => Ord (Ratio a) where
475 {-# SPECIALIZE instance Ord Rational #-}
476 (x:%y) <= (x':%y') = x * y' <= x' * y
477 (x:%y) < (x':%y') = x * y' < x' * y
478
479 -- | @since 2.0.1
480 instance (Integral a) => Num (Ratio a) where
481 {-# SPECIALIZE instance Num Rational #-}
482 (x:%y) + (x':%y') = reduce (x*y' + x'*y) (y*y')
483 (x:%y) - (x':%y') = reduce (x*y' - x'*y) (y*y')
484 (x:%y) * (x':%y') = reduce (x * x') (y * y')
485 negate (x:%y) = (-x) :% y
486 abs (x:%y) = abs x :% y
487 signum (x:%_) = signum x :% 1
488 fromInteger x = fromInteger x :% 1
489
490 -- | @since 2.0.1
491 {-# RULES "fromRational/id" fromRational = id :: Rational -> Rational #-}
492 instance (Integral a) => Fractional (Ratio a) where
493 {-# SPECIALIZE instance Fractional Rational #-}
494 (x:%y) / (x':%y') = (x*y') % (y*x')
495 recip (0:%_) = ratioZeroDenominatorError
496 recip (x:%y)
497 | x < 0 = negate y :% negate x
498 | otherwise = y :% x
499 fromRational (x:%y) = fromInteger x % fromInteger y
500
501 -- | @since 2.0.1
502 instance (Integral a) => Real (Ratio a) where
503 {-# SPECIALIZE instance Real Rational #-}
504 toRational (x:%y) = toInteger x :% toInteger y
505
506 -- | @since 2.0.1
507 instance (Integral a) => RealFrac (Ratio a) where
508 {-# SPECIALIZE instance RealFrac Rational #-}
509 properFraction (x:%y) = (fromInteger (toInteger q), r:%y)
510 where (q,r) = quotRem x y
511
512 -- | @since 2.0.1
513 instance (Show a) => Show (Ratio a) where
514 {-# SPECIALIZE instance Show Rational #-}
515 showsPrec p (x:%y) = showParen (p > ratioPrec) $
516 showsPrec ratioPrec1 x .
517 showString " % " .
518 -- H98 report has spaces round the %
519 -- but we removed them [May 04]
520 -- and added them again for consistency with
521 -- Haskell 98 [Sep 08, #1920]
522 showsPrec ratioPrec1 y
523
524 -- | @since 2.0.1
525 instance (Integral a) => Enum (Ratio a) where
526 {-# SPECIALIZE instance Enum Rational #-}
527 succ x = x + 1
528 pred x = x - 1
529
530 toEnum n = fromIntegral n :% 1
531 fromEnum = fromInteger . truncate
532
533 enumFrom = numericEnumFrom
534 enumFromThen = numericEnumFromThen
535 enumFromTo = numericEnumFromTo
536 enumFromThenTo = numericEnumFromThenTo
537
538 --------------------------------------------------------------
539 -- Coercions
540 --------------------------------------------------------------
541
542 -- | general coercion from integral types
543 {-# NOINLINE [1] fromIntegral #-}
544 fromIntegral :: (Integral a, Num b) => a -> b
545 fromIntegral = fromInteger . toInteger
546
547 {-# RULES
548 "fromIntegral/Int->Int" fromIntegral = id :: Int -> Int
549 #-}
550
551 {-# RULES
552 "fromIntegral/Int->Word" fromIntegral = \(I# x#) -> W# (int2Word# x#)
553 "fromIntegral/Word->Int" fromIntegral = \(W# x#) -> I# (word2Int# x#)
554 "fromIntegral/Word->Word" fromIntegral = id :: Word -> Word
555 #-}
556
557 {-# RULES
558 "fromIntegral/Natural->Natural" fromIntegral = id :: Natural -> Natural
559 "fromIntegral/Natural->Integer" fromIntegral = toInteger :: Natural->Integer
560 "fromIntegral/Natural->Word" fromIntegral = naturalToWord
561 #-}
562
563 {-# RULES
564 "fromIntegral/Word->Natural" fromIntegral = wordToNatural
565 "fromIntegral/Int->Natural" fromIntegral = intToNatural
566 #-}
567
568 -- | general coercion to fractional types
569 realToFrac :: (Real a, Fractional b) => a -> b
570 {-# NOINLINE [1] realToFrac #-}
571 realToFrac = fromRational . toRational
572
573 --------------------------------------------------------------
574 -- Overloaded numeric functions
575 --------------------------------------------------------------
576
577 -- | Converts a possibly-negative 'Real' value to a string.
578 showSigned :: (Real a)
579 => (a -> ShowS) -- ^ a function that can show unsigned values
580 -> Int -- ^ the precedence of the enclosing context
581 -> a -- ^ the value to show
582 -> ShowS
583 showSigned showPos p x
584 | x < 0 = showParen (p > 6) (showChar '-' . showPos (-x))
585 | otherwise = showPos x
586
587 even, odd :: (Integral a) => a -> Bool
588 even n = n `rem` 2 == 0
589 odd = not . even
590 {-# INLINABLE even #-}
591 {-# INLINABLE odd #-}
592
593 -------------------------------------------------------
594 -- | raise a number to a non-negative integral power
595 {-# SPECIALISE [1] (^) ::
596 Integer -> Integer -> Integer,
597 Integer -> Int -> Integer,
598 Int -> Int -> Int #-}
599 {-# INLINABLE [1] (^) #-} -- See Note [Inlining (^)]
600 (^) :: (Num a, Integral b) => a -> b -> a
601 x0 ^ y0 | y0 < 0 = errorWithoutStackTrace "Negative exponent"
602 | y0 == 0 = 1
603 | otherwise = f x0 y0
604 where -- f : x0 ^ y0 = x ^ y
605 f x y | even y = f (x * x) (y `quot` 2)
606 | y == 1 = x
607 | otherwise = g (x * x) (y `quot` 2) x -- See Note [Half of y - 1]
608 -- g : x0 ^ y0 = (x ^ y) * z
609 g x y z | even y = g (x * x) (y `quot` 2) z
610 | y == 1 = x * z
611 | otherwise = g (x * x) (y `quot` 2) (x * z) -- See Note [Half of y - 1]
612
613 -- | raise a number to an integral power
614 (^^) :: (Fractional a, Integral b) => a -> b -> a
615 {-# INLINABLE [1] (^^) #-} -- See Note [Inlining (^)
616 x ^^ n = if n >= 0 then x^n else recip (x^(negate n))
617
618 {- Note [Half of y - 1]
619 ~~~~~~~~~~~~~~~~~~~~~
620 Since y is guaranteed to be odd and positive here,
621 half of y - 1 can be computed as y `quot` 2, optimising subtraction away.
622 -}
623
624 {- Note [Inlining (^)
625 ~~~~~~~~~~~~~~~~~~~~~
626 The INLINABLE pragma allows (^) to be specialised at its call sites.
627 If it is called repeatedly at the same type, that can make a huge
628 difference, because of those constants which can be repeatedly
629 calculated.
630
631 Currently the fromInteger calls are not floated because we get
632 \d1 d2 x y -> blah
633 after the gentle round of simplification. -}
634
635 {- Rules for powers with known small exponent
636 see #5237
637 For small exponents, (^) is inefficient compared to manually
638 expanding the multiplication tree.
639 Here, rules for the most common exponent types are given.
640 The range of exponents for which rules are given is quite
641 arbitrary and kept small to not unduly increase the number of rules.
642 0 and 1 are excluded based on the assumption that nobody would
643 write x^0 or x^1 in code and the cases where an exponent could
644 be statically resolved to 0 or 1 are rare.
645
646 It might be desirable to have corresponding rules also for
647 exponents of other types (e. g., Word), but it's doubtful they
648 would fire, since the exponents of other types tend to get
649 floated out before the rule has a chance to fire.
650
651 Also desirable would be rules for (^^), but I haven't managed
652 to get those to fire.
653
654 Note: Trying to save multiplications by sharing the square for
655 exponents 4 and 5 does not save time, indeed, for Double, it is
656 up to twice slower, so the rules contain flat sequences of
657 multiplications.
658 -}
659
660 {-# RULES
661 "^2/Int" forall x. x ^ (2 :: Int) = let u = x in u*u
662 "^3/Int" forall x. x ^ (3 :: Int) = let u = x in u*u*u
663 "^4/Int" forall x. x ^ (4 :: Int) = let u = x in u*u*u*u
664 "^5/Int" forall x. x ^ (5 :: Int) = let u = x in u*u*u*u*u
665 "^2/Integer" forall x. x ^ (2 :: Integer) = let u = x in u*u
666 "^3/Integer" forall x. x ^ (3 :: Integer) = let u = x in u*u*u
667 "^4/Integer" forall x. x ^ (4 :: Integer) = let u = x in u*u*u*u
668 "^5/Integer" forall x. x ^ (5 :: Integer) = let u = x in u*u*u*u*u
669 #-}
670
671 -------------------------------------------------------
672 -- Special power functions for Rational
673 --
674 -- see #4337
675 --
676 -- Rationale:
677 -- For a legitimate Rational (n :% d), the numerator and denominator are
678 -- coprime, i.e. they have no common prime factor.
679 -- Therefore all powers (n ^ a) and (d ^ b) are also coprime, so it is
680 -- not necessary to compute the greatest common divisor, which would be
681 -- done in the default implementation at each multiplication step.
682 -- Since exponentiation quickly leads to very large numbers and
683 -- calculation of gcds is generally very slow for large numbers,
684 -- avoiding the gcd leads to an order of magnitude speedup relatively
685 -- soon (and an asymptotic improvement overall).
686 --
687 -- Note:
688 -- We cannot use these functions for general Ratio a because that would
689 -- change results in a multitude of cases.
690 -- The cause is that if a and b are coprime, their remainders by any
691 -- positive modulus generally aren't, so in the default implementation
692 -- reduction occurs.
693 --
694 -- Example:
695 -- (17 % 3) ^ 3 :: Ratio Word8
696 -- Default:
697 -- (17 % 3) ^ 3 = ((17 % 3) ^ 2) * (17 % 3)
698 -- = ((289 `mod` 256) % 9) * (17 % 3)
699 -- = (33 % 9) * (17 % 3)
700 -- = (11 % 3) * (17 % 3)
701 -- = (187 % 9)
702 -- But:
703 -- ((17^3) `mod` 256) % (3^3) = (4913 `mod` 256) % 27
704 -- = 49 % 27
705 --
706 -- TODO:
707 -- Find out whether special-casing for numerator, denominator or
708 -- exponent = 1 (or -1, where that may apply) gains something.
709
710 -- Special version of (^) for Rational base
711 {-# RULES "(^)/Rational" (^) = (^%^) #-}
712 (^%^) :: Integral a => Rational -> a -> Rational
713 (n :% d) ^%^ e
714 | e < 0 = errorWithoutStackTrace "Negative exponent"
715 | e == 0 = 1 :% 1
716 | otherwise = (n ^ e) :% (d ^ e)
717
718 -- Special version of (^^) for Rational base
719 {-# RULES "(^^)/Rational" (^^) = (^^%^^) #-}
720 (^^%^^) :: Integral a => Rational -> a -> Rational
721 (n :% d) ^^%^^ e
722 | e > 0 = (n ^ e) :% (d ^ e)
723 | e == 0 = 1 :% 1
724 | n > 0 = (d ^ (negate e)) :% (n ^ (negate e))
725 | n == 0 = ratioZeroDenominatorError
726 | otherwise = let nn = d ^ (negate e)
727 dd = (negate n) ^ (negate e)
728 in if even e then (nn :% dd) else (negate nn :% dd)
729
730 -------------------------------------------------------
731 -- | @'gcd' x y@ is the non-negative factor of both @x@ and @y@ of which
732 -- every common factor of @x@ and @y@ is also a factor; for example
733 -- @'gcd' 4 2 = 2@, @'gcd' (-4) 6 = 2@, @'gcd' 0 4@ = @4@. @'gcd' 0 0@ = @0@.
734 -- (That is, the common divisor that is \"greatest\" in the divisibility
735 -- preordering.)
736 --
737 -- Note: Since for signed fixed-width integer types, @'abs' 'minBound' < 0@,
738 -- the result may be negative if one of the arguments is @'minBound'@ (and
739 -- necessarily is if the other is @0@ or @'minBound'@) for such types.
740 gcd :: (Integral a) => a -> a -> a
741 {-# NOINLINE [1] gcd #-}
742 gcd x y = gcd' (abs x) (abs y)
743 where gcd' a 0 = a
744 gcd' a b = gcd' b (a `rem` b)
745
746 -- | @'lcm' x y@ is the smallest positive integer that both @x@ and @y@ divide.
747 lcm :: (Integral a) => a -> a -> a
748 {-# SPECIALISE lcm :: Int -> Int -> Int #-}
749 {-# SPECIALISE lcm :: Word -> Word -> Word #-}
750 {-# NOINLINE [1] lcm #-}
751 lcm _ 0 = 0
752 lcm 0 _ = 0
753 lcm x y = abs ((x `quot` (gcd x y)) * y)
754
755 {-# RULES
756 "gcd/Integer->Integer->Integer" gcd = gcdInteger
757 "lcm/Integer->Integer->Integer" lcm = lcmInteger
758 "gcd/Natural->Natural->Natural" gcd = gcdNatural
759 "lcm/Natural->Natural->Natural" lcm = lcmNatural
760 #-}
761
762 #if defined(MIN_VERSION_integer_gmp)
763 -- GMP defines a more efficient Int# and Word# GCD
764
765 gcdInt' :: Int -> Int -> Int
766 gcdInt' (I# x) (I# y) = I# (gcdInt x y)
767
768 gcdWord' :: Word -> Word -> Word
769 gcdWord' (W# x) (W# y) = W# (gcdWord x y)
770
771 {-# RULES
772 "gcd/Int->Int->Int" gcd = gcdInt'
773 "gcd/Word->Word->Word" gcd = gcdWord'
774 #-}
775
776 #endif
777
778 integralEnumFrom :: (Integral a, Bounded a) => a -> [a]
779 integralEnumFrom n = map fromInteger [toInteger n .. toInteger (maxBound `asTypeOf` n)]
780
781 integralEnumFromThen :: (Integral a, Bounded a) => a -> a -> [a]
782 integralEnumFromThen n1 n2
783 | i_n2 >= i_n1 = map fromInteger [i_n1, i_n2 .. toInteger (maxBound `asTypeOf` n1)]
784 | otherwise = map fromInteger [i_n1, i_n2 .. toInteger (minBound `asTypeOf` n1)]
785 where
786 i_n1 = toInteger n1
787 i_n2 = toInteger n2
788
789 integralEnumFromTo :: Integral a => a -> a -> [a]
790 integralEnumFromTo n m = map fromInteger [toInteger n .. toInteger m]
791
792 integralEnumFromThenTo :: Integral a => a -> a -> a -> [a]
793 integralEnumFromThenTo n1 n2 m
794 = map fromInteger [toInteger n1, toInteger n2 .. toInteger m]