1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
|
-- |
-- Module : Unicode.Internal.Division
-- Copyright : (c) 2020 Andrew Lelechenko
--
-- License : BSD-3-Clause
-- Maintainer : streamly@composewell.com
-- Stability : experimental
-- Portability : GHC
--
-- Fast division by known constants.
--
-- Division by a constant can be replaced by a double-word multiplication.
-- Roughly speaking, instead of dividing by x, multiply by 2^64/x,
-- obtaining 128-bit-long product, and take upper 64 bits. The peculiar
-- details can be found in Hacker's Delight, Ch. 10.
--
-- Even GHC 8.10 does not provide a primitive for a signed double-word
-- multiplication, but since our applications does not involve negative
-- integers, we convert 'Int' to 'Word' and use 'GHC.Exts.timesWord#'.
--
-- Textbook unsigned division by 21 or 28 becomes involved, when an argument
-- is allowed to take the full range of 'Word' up to 2^64. Luckily, in our
-- case the argument was casted from 'Int', so we can guarantee that it is
-- below 2^63.
module Unicode.Internal.Division
( quotRem21
, quotRem28
) where
import Data.Bits (Bits(..), FiniteBits(..))
import GHC.Exts (Word(..), timesWord2#)
highMul :: Word -> Word -> Word
highMul (W# x#) (W# y#) = W# high#
where
!(# high#, _ #) = timesWord2# x# y#
-- | Input must be non-negative.
--
-- Instead of division by 21, we compute
-- floor(floor((2^68+17)/21 * n) / 2^68) = floor((2^68+17)/21 * n/2^68) =
-- floor(n/21 + (n/2^63 * 17/32)/21) = floor(n/21),
-- because n/2^63 * 17/32 < 1.
{-# INLINE quotRem21 #-}
quotRem21 :: Int -> (Int, Int)
quotRem21 n
| finiteBitSize (0 :: Word) /= 64
= n `quotRem` 21
| otherwise
= (fromIntegral q, fromIntegral (w - 21 * q))
where
w = fromIntegral n
high = highMul w 14054662151397753613 -- (2^68+17)/21
q = high `shiftR` 4
-- | Input must be non-negative.
--
-- Instead of division by 28, we compute
-- floor(floor((2^65+3)/7 * n) / 2^67) = floor((2^65+3)/7 * n/2^67) =
-- floor(n/28 + (n/2^63 * 3/4)/28) = floor(n/28),
-- because n/2^63 * 3/4 < 1.
{-# INLINE quotRem28 #-}
quotRem28 :: Int -> (Int, Int)
quotRem28 n
| finiteBitSize (0 :: Word) /= 64
= n `quotRem` 28
| otherwise
= (fromIntegral q, fromIntegral r)
where
w = fromIntegral n
high = highMul w 5270498306774157605 -- (2^65+3)/7
q = high `shiftR` 3
prod = (q `shiftL` 3 - q) `shiftL` 2
r = w - prod
|
