    -- |
    -- Module:      Math.NumberTheory.Roots.Cubes
    -- Copyright:   (c) 2011 Daniel Fischer, 2016-2020 Andrew Lelechenko
    -- Licence:     MIT
    -- Maintainer:  Daniel Fischer <daniel.is.fischer@googlemail.com>
    --
    -- Functions dealing with cubes. Moderately efficient calculation of integer
    -- cube roots and testing for cubeness.

    {-# LANGUAGE BangPatterns #-}
    {-# LANGUAGE CPP          #-}
    {-# LANGUAGE MagicHash    #-}
    {- HLINT ignore "Use fewer imports" -}

    module Math.NumberTheory.Roots.Cubes
        ( integerCubeRoot
        , integerCubeRoot'
        , exactCubeRoot
        , isCube
        , isCube'
        , isPossibleCube
        ) where

    import Data.Bits (finiteBitSize, (.&.))
    import GHC.Exts (Int#, Ptr(..), int2Double#, double2Int#, isTrue#, (/##), (**##), (<#))
    import Numeric.Natural (Natural)

    #ifdef MIN_VERSION_integer_gmp
    import GHC.Exts (quotInt#, (*#), (-#))
    import GHC.Integer.GMP.Internals (Integer(..), shiftLInteger, shiftRInteger, sizeofBigNat#)
    import GHC.Integer.Logarithms (integerLog2#)
    #define IS S#
    #define IP Jp#
    #define bigNatSize sizeofBigNat
    #else
    import GHC.Exts (minusWord#, timesWord#, quotWord#)
    import GHC.Num.BigNat (bigNatSize#)
    import GHC.Num.Integer (Integer(..), integerLog2#, integerShiftR#, integerShiftL#)
    #endif

    import Math.NumberTheory.Utils.BitMask (indexBitSet)

    -- | For a given \( n \)
    --   calculate its integer cube root \( \lfloor \sqrt[3]{n} \rfloor \).
    --   Note that this is not symmetric about 0.
    --
    -- >>> map integerCubeRoot [7, 8, 9]
    -- [1,2,2]
    -- >>> map integerCubeRoot [-7, -8, -9]
    -- [-2,-2,-3]
    {-# SPECIALISE integerCubeRoot :: Int -> Int #-}
    {-# SPECIALISE integerCubeRoot :: Word -> Word #-}
    {-# SPECIALISE integerCubeRoot :: Integer -> Integer #-}
    {-# SPECIALISE integerCubeRoot :: Natural -> Natural #-}
    integerCubeRoot :: Integral a => a -> a
    integerCubeRoot :: forall a. Integral a => a -> a
    integerCubeRoot a
    0 = a
    0
    integerCubeRoot a
    n
        | a
    n a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    > a
    0     = a -> a
    forall a. Integral a => a -> a
    integerCubeRoot' a
    n
        | Bool
    otherwise =
          let m :: a
    m = a -> a
    forall a. Num a => a -> a
    negate a
    n
              r :: a
    r = if a
    m a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    < a
    0
                    then a -> a
    forall a. Num a => a -> a
    negate (a -> a) -> (Integer -> a) -> Integer -> a
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Integer -> a
    forall a. Num a => Integer -> a
    fromInteger (Integer -> a) -> Integer -> a
    forall a b. (a -> b) -> a -> b
    $ Integer -> Integer
    forall a. Integral a => a -> a
    integerCubeRoot' (Integer -> Integer
    forall a. Num a => a -> a
    negate (Integer -> Integer) -> Integer -> Integer
    forall a b. (a -> b) -> a -> b
    $ a -> Integer
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral a
    n)
                    else a -> a
    forall a. Num a => a -> a
    negate (a -> a
    forall a. Integral a => a -> a
    integerCubeRoot' a
    m)
          in if a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    r a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == a
    n then a
    r else a
    r a -> a -> a
    forall a. Num a => a -> a -> a
    - a
    1

    -- | Calculate the integer cube root of a nonnegative integer @n@,
    --   that is, the largest integer @r@ such that @r^3 <= n@.
    --   The precondition @n >= 0@ is not checked.
    {-# RULES
    "integerCubeRoot'/Int"     integerCubeRoot' = cubeRootInt'
    "integerCubeRoot'/Word"    integerCubeRoot' = cubeRootWord
    "integerCubeRoot'/Integer" integerCubeRoot' = cubeRootIgr
      #-}
    {-# INLINE [1] integerCubeRoot' #-}
    integerCubeRoot' :: Integral a => a -> a
    integerCubeRoot' :: forall a. Integral a => a -> a
    integerCubeRoot' a
    0 = a
    0
    integerCubeRoot' a
    n = a -> a -> a
    forall a. Integral a => a -> a -> a
    newton3 a
    n (a -> a
    forall a. Integral a => a -> a
    approxCuRt a
    n)

    -- | Calculate the exact integer cube root if it exists,
    -- otherwise return 'Nothing'.
    --
    -- >>> map exactCubeRoot [-9, -8, -7, 7, 8, 9]
    -- [Nothing,Just (-2),Nothing,Nothing,Just 2,Nothing]
    {-# SPECIALISE exactCubeRoot :: Int -> Maybe Int #-}
    {-# SPECIALISE exactCubeRoot :: Word -> Maybe Word #-}
    {-# SPECIALISE exactCubeRoot :: Integer -> Maybe Integer #-}
    {-# SPECIALISE exactCubeRoot :: Natural -> Maybe Natural #-}
    exactCubeRoot :: Integral a => a -> Maybe a
    exactCubeRoot :: forall a. Integral a => a -> Maybe a
    exactCubeRoot a
    0 = a -> Maybe a
    forall a. a -> Maybe a
    Just a
    0
    exactCubeRoot a
    n
        | a
    n a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    < a
    0     =
          if a
    m a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    < a
    0
            then (Integer -> a) -> Maybe Integer -> Maybe a
    forall a b. (a -> b) -> Maybe a -> Maybe b
    forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
    fmap (a -> a
    forall a. Num a => a -> a
    negate (a -> a) -> (Integer -> a) -> Integer -> a
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . Integer -> a
    forall a. Num a => Integer -> a
    fromInteger) (Maybe Integer -> Maybe a) -> Maybe Integer -> Maybe a
    forall a b. (a -> b) -> a -> b
    $ Integer -> Maybe Integer
    forall a. Integral a => a -> Maybe a
    exactCubeRoot (Integer -> Integer
    forall a. Num a => a -> a
    negate (Integer -> Integer) -> Integer -> Integer
    forall a b. (a -> b) -> a -> b
    $ a -> Integer
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral a
    n)
            else (a -> a) -> Maybe a -> Maybe a
    forall a b. (a -> b) -> Maybe a -> Maybe b
    forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
    fmap a -> a
    forall a. Num a => a -> a
    negate (a -> Maybe a
    forall a. Integral a => a -> Maybe a
    exactCubeRoot a
    m)
        | a -> Bool
    forall a. Integral a => a -> Bool
    isPossibleCube a
    n Bool -> Bool -> Bool
    && a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    r a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == a
    n    = a -> Maybe a
    forall a. a -> Maybe a
    Just a
    r
        | Bool
    otherwise = Maybe a
    forall a. Maybe a
    Nothing
          where
            m :: a
    m = a -> a
    forall a. Num a => a -> a
    negate a
    n
            r :: a
    r = a -> a
    forall a. Integral a => a -> a
    integerCubeRoot' a
    n

    -- | Test whether the argument is a perfect cube.
    --
    -- >>> map isCube [-9, -8, -7, 7, 8, 9]
    -- [False,True,False,False,True,False]
    {-# SPECIALISE isCube :: Int -> Bool #-}
    {-# SPECIALISE isCube :: Word -> Bool #-}
    {-# SPECIALISE isCube :: Integer -> Bool #-}
    {-# SPECIALISE isCube :: Natural -> Bool #-}
    isCube :: Integral a => a -> Bool
    isCube :: forall a. Integral a => a -> Bool
    isCube a
    0 = Bool
    True
    isCube a
    n
        | a
    n a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    > a
    0     = a -> Bool
    forall a. Integral a => a -> Bool
    isCube' a
    n
        | a
    m a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    > a
    0     = a -> Bool
    forall a. Integral a => a -> Bool
    isCube' a
    m
        | Bool
    otherwise = Integer -> Bool
    forall a. Integral a => a -> Bool
    isCube' (Integer -> Integer
    forall a. Num a => a -> a
    negate (a -> Integer
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral a
    n) :: Integer)
          where
            m :: a
    m = a -> a
    forall a. Num a => a -> a
    negate a
    n

    -- | Test whether a nonnegative integer is a cube.
    --   Before 'integerCubeRoot' is calculated, a few tests
    --   of remainders modulo small primes weed out most non-cubes.
    --   On average, assuming that the majority of inputs aren't cubes,
    --   this is much faster than @let r = cubeRoot n in r*r*r == n@.
    --   The condition @n >= 0@ is /not/ checked.
    {-# SPECIALISE isCube' :: Int -> Bool #-}
    {-# SPECIALISE isCube' :: Word -> Bool #-}
    {-# SPECIALISE isCube' :: Integer -> Bool #-}
    {-# SPECIALISE isCube' :: Natural -> Bool #-}
    isCube' :: Integral a => a -> Bool
    isCube' :: forall a. Integral a => a -> Bool
    isCube' !a
    n = a -> Bool
    forall a. Integral a => a -> Bool
    isPossibleCube a
    n
                 Bool -> Bool -> Bool
    && (a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    ra -> a -> a
    forall a. Num a => a -> a -> a
    *a
    r a -> a -> Bool
    forall a. Eq a => a -> a -> Bool
    == a
    n)
          where
            r :: a
    r    = a -> a
    forall a. Integral a => a -> a
    integerCubeRoot' a
    n

    -- | Test whether a nonnegative number is possibly a cube.
    --   Only about 0.08% of all numbers pass this test.
    --   The precondition @n >= 0@ is /not/ checked.
    {-# SPECIALISE isPossibleCube :: Int -> Bool #-}
    {-# SPECIALISE isPossibleCube :: Word -> Bool #-}
    {-# SPECIALISE isPossibleCube :: Integer -> Bool #-}
    {-# SPECIALISE isPossibleCube :: Natural -> Bool #-}
    isPossibleCube :: Integral a => a -> Bool
    isPossibleCube :: forall a. Integral a => a -> Bool
    isPossibleCube a
    n'
      =  Ptr Word -> Int -> Bool
    indexBitSet Ptr Word
    mask512 (Integer -> Int
    forall a. Num a => Integer -> a
    fromInteger (Integer
    n Integer -> Integer -> Integer
    forall a. Bits a => a -> a -> a
    .&. Integer
    511))
      Bool -> Bool -> Bool
    && Ptr Word -> Int -> Bool
    indexBitSet Ptr Word
    mask837 (Integer -> Int
    forall a. Num a => Integer -> a
    fromInteger (Integer
    n Integer -> Integer -> Integer
    forall a. Integral a => a -> a -> a
    `rem` Integer
    837))
      Bool -> Bool -> Bool
    && Ptr Word -> Int -> Bool
    indexBitSet Ptr Word
    mask637 (Integer -> Int
    forall a. Num a => Integer -> a
    fromInteger (Integer
    n Integer -> Integer -> Integer
    forall a. Integral a => a -> a -> a
    `rem` Integer
    637))
      Bool -> Bool -> Bool
    && Ptr Word -> Int -> Bool
    indexBitSet Ptr Word
    mask703 (Integer -> Int
    forall a. Num a => Integer -> a
    fromInteger (Integer
    n Integer -> Integer -> Integer
    forall a. Integral a => a -> a -> a
    `rem` Integer
    703))
      where
        n :: Integer
    n = a -> Integer
    forall a. Integral a => a -> Integer
    toInteger a
    n'

    ----------------------------------------------------------------------
    --                         Utility Functions                        --
    ----------------------------------------------------------------------

    -- Special case for 'Int', a little faster.
    -- For @n <= 2^64@, the truncated 'Double' is never
    -- more than one off. Things might overflow for @n@
    -- close to @maxBound@, so check for overflow.
    cubeRootInt' :: Int -> Int
    cubeRootInt' :: Int -> Int
    cubeRootInt' Int
    0 = Int
    0
    cubeRootInt' Int
    n
        | Int
    n Int -> Int -> Bool
    forall a. Ord a => a -> a -> Bool
    < Int
    c Bool -> Bool -> Bool
    || Int
    c Int -> Int -> Bool
    forall a. Ord a => a -> a -> Bool
    < Int
    0    = Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    -Int
    1
        | Int
    0 Int -> Int -> Bool
    forall a. Ord a => a -> a -> Bool
    < Int
    d Bool -> Bool -> Bool
    && Int
    d Int -> Int -> Bool
    forall a. Ord a => a -> a -> Bool
    < Int
    n    = Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +Int
    1
        | Bool
    otherwise         = Int
    r
          where
            x :: Double
    x = Int -> Double
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral Int
    n :: Double
            r :: Int
    r = Double -> Int
    forall b. Integral b => Double -> b
    forall a b. (RealFrac a, Integral b) => a -> b
    truncate (Double
    x Double -> Double -> Double
    forall a. Floating a => a -> a -> a
    ** (Double
    1Double -> Double -> Double
    forall a. Fractional a => a -> a -> a
    /Double
    3))
            c :: Int
    c = Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    *Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    *Int
    r
            d :: Int
    d = Int
    cInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +Int
    3Int -> Int -> Int
    forall a. Num a => a -> a -> a
    *Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    *(Int
    rInt -> Int -> Int
    forall a. Num a => a -> a -> a
    +Int
    1)

    cubeRootWordLimit :: Word
    cubeRootWordLimit :: Word
    cubeRootWordLimit = if Word -> Int
    forall b. FiniteBits b => b -> Int
    finiteBitSize (Word
    0 :: Word) Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    64 then Word
    2642245 else Word
    1625

    cubeRootWord :: Word -> Word
    cubeRootWord :: Word -> Word
    cubeRootWord Word
    0 = Word
    0
    cubeRootWord Word
    w
        | Word
    r Word -> Word -> Bool
    forall a. Ord a => a -> a -> Bool
    > Word
    cubeRootWordLimit = Word
    cubeRootWordLimit
        | Word
    w Word -> Word -> Bool
    forall a. Ord a => a -> a -> Bool
    < Word
    c             = Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    -Word
    1
        | Word
    c Word -> Word -> Bool
    forall a. Ord a => a -> a -> Bool
    < Word
    w Bool -> Bool -> Bool
    && Word
    e Word -> Word -> Bool
    forall a. Ord a => a -> a -> Bool
    < Word
    w Bool -> Bool -> Bool
    && Word
    c Word -> Word -> Bool
    forall a. Ord a => a -> a -> Bool
    < Word
    e  = Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    +Word
    1
        | Bool
    otherwise         = Word
    r
          where
            r :: Word
    r = Double -> Word
    forall b. Integral b => Double -> b
    forall a b. (RealFrac a, Integral b) => a -> b
    truncate (Word -> Double
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral Word
    w Double -> Double -> Double
    forall a. Floating a => a -> a -> a
    ** (Double
    1Double -> Double -> Double
    forall a. Fractional a => a -> a -> a
    /Double
    3) :: Double)
            c :: Word
    c = Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    *Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    *Word
    r
            d :: Word
    d = Word
    3Word -> Word -> Word
    forall a. Num a => a -> a -> a
    *Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    *(Word
    rWord -> Word -> Word
    forall a. Num a => a -> a -> a
    +Word
    1)
            e :: Word
    e = Word
    cWord -> Word -> Word
    forall a. Num a => a -> a -> a
    +Word
    d

    cubeRootIgr :: Integer -> Integer
    cubeRootIgr :: Integer -> Integer
    cubeRootIgr Integer
    0 = Integer
    0
    cubeRootIgr Integer
    n = Integer -> Integer -> Integer
    forall a. Integral a => a -> a -> a
    newton3 Integer
    n (Integer -> Integer
    forall a. Integral a => a -> a
    approxCuRt Integer
    n)

    {-# SPECIALISE newton3 :: Integer -> Integer -> Integer #-}
    newton3 :: Integral a => a -> a -> a
    newton3 :: forall a. Integral a => a -> a -> a
    newton3 a
    n a
    a = a -> a
    go (a -> a
    step a
    a)
          where
            step :: a -> a
    step a
    k = (a
    2a -> a -> a
    forall a. Num a => a -> a -> a
    *a
    k a -> a -> a
    forall a. Num a => a -> a -> a
    + a
    n a -> a -> a
    forall a. Integral a => a -> a -> a
    `quot` (a
    ka -> a -> a
    forall a. Num a => a -> a -> a
    *a
    k)) a -> a -> a
    forall a. Integral a => a -> a -> a
    `quot` a
    3
            go :: a -> a
    go a
    k
                | a
    m a -> a -> Bool
    forall a. Ord a => a -> a -> Bool
    < a
    k     = a -> a
    go a
    m
                | Bool
    otherwise = a
    k
                  where
                    m :: a
    m = a -> a
    step a
    k

    {-# SPECIALISE approxCuRt :: Integer -> Integer #-}
    approxCuRt :: Integral a => a -> a
    approxCuRt :: forall a. Integral a => a -> a
    approxCuRt a
    0 = a
    0
    approxCuRt a
    n = Integer -> a
    forall a. Num a => Integer -> a
    fromInteger (Integer -> a) -> Integer -> a
    forall a b. (a -> b) -> a -> b
    $ Integer -> Integer
    appCuRt (a -> Integer
    forall a b. (Integral a, Num b) => a -> b
    fromIntegral a
    n)

    -- | approximate cube root, about 50 bits should be correct for large numbers
    appCuRt :: Integer -> Integer
    appCuRt :: Integer -> Integer
    appCuRt (IS Int#
    i#) = case Double# -> Int#
    double2Int# (Int# -> Double#
    int2Double# Int#
    i# Double# -> Double# -> Double#
    **## (Double#
    1.0## Double# -> Double# -> Double#
    /## Double#
    3.0##)) of
                        Int#
    r# -> Int# -> Integer
    IS Int#
    r#
    appCuRt n :: Integer
    n@(IP ByteArray#
    bn#)
        | Int# -> Bool
    isTrue# (ByteArray# -> Int#
    bigNatSize# ByteArray#
    bn# Int# -> Int# -> Int#
    <# Int#
    thresh#) =
              Double -> Integer
    forall b. Integral b => Double -> b
    forall a b. (RealFrac a, Integral b) => a -> b
    floor (Integer -> Double
    forall a. Num a => Integer -> a
    fromInteger Integer
    n Double -> Double -> Double
    forall a. Floating a => a -> a -> a
    ** (Double
    1.0Double -> Double -> Double
    forall a. Fractional a => a -> a -> a
    /Double
    3.0) :: Double)
        | Bool
    otherwise = case Integer -> Word#
    integerLog2# Integer
    n of
    #ifdef MIN_VERSION_integer_gmp
                        l# -> case (l# `quotInt#` 3#) -# 51# of
                                h# -> case shiftRInteger n (3# *# h#) of
                                        m -> case floor (fromInteger m ** (1.0/3.0) :: Double) of
                                               r -> shiftLInteger r h#
    #else
                        Word#
    l# -> case (Word#
    l# Word# -> Word# -> Word#
    `quotWord#` Word#
    3##) Word# -> Word# -> Word#
    `minusWord#` Word#
    51## of
                                Word#
    h# -> case Integer -> Word# -> Integer
    integerShiftR# Integer
    n (Word#
    3## Word# -> Word# -> Word#
    `timesWord#` Word#
    h#) of
                                        Integer
    m -> case Double -> Integer
    forall b. Integral b => Double -> b
    forall a b. (RealFrac a, Integral b) => a -> b
    floor (Integer -> Double
    forall a. Num a => Integer -> a
    fromInteger Integer
    m Double -> Double -> Double
    forall a. Floating a => a -> a -> a
    ** (Double
    1.0Double -> Double -> Double
    forall a. Fractional a => a -> a -> a
    /Double
    3.0) :: Double) of
                                                Integer
    r -> Integer -> Word# -> Integer
    integerShiftL# Integer
    r Word#
    h#
    #endif
        where
            -- threshold for shifting vs. direct fromInteger
            -- we shift when we expect more than 256 bits
            thresh# :: Int#
            thresh# :: Int#
    thresh# = if Word -> Int
    forall b. FiniteBits b => b -> Int
    finiteBitSize (Word
    0 :: Word) Int -> Int -> Bool
    forall a. Eq a => a -> a -> Bool
    == Int
    64 then Int#
    5# else Int#
    9#

    -- There's already handling for negative in integerCubeRoot,
    -- but integerCubeRoot' is exported directly too.
    appCuRt Integer
    _ = [Char] -> Integer
    forall a. HasCallStack => [Char] -> a
    error [Char]
    "integerCubeRoot': negative argument"

    -----------------------------------------------------------------------------
    -- Generated by 'Math.NumberTheory.Utils.BitMask.vectorToAddrLiteral'

    mask512 :: Ptr Word
    mask512 :: Ptr Word
    mask512 = Addr# -> Ptr Word
    forall a. Addr# -> Ptr a
    Ptr Addr#
    "\171\171\170\171\170\171\170\171\171\171\170\171\170\171\170\171\170\171\170\171\170\171\170\171\171\171\170\171\170\171\170\171\170\171\170\171\170\171\170\171\171\171\170\171\170\171\170\171\170\171\170\171\170\171\170\171\171\171\170\171\170\171\170\171"#

    mask837 :: Ptr Word
    mask837 :: Ptr Word
    mask837 = Addr# -> Ptr Word
    forall a. Addr# -> Ptr a
    Ptr Addr#
    "\ETX\SOH\NUL\b\b@@@\SOH\NUL\NUL\n\NUL0\DLE \NUL\NUL\NUL\EOT\b\EOT\NULP\NUL\NUL\128\ETX\NUL\STX\DLE\NUL\NUL\128@\129\NUL\NUL\NUL\EOT \NUL\160\NUL\NUL\NUL\ENQ\NUL\DLE\b0\NUL\NUL\NUL\ETX\EOT\STX\NUL(\NUL\NUL@\SOH\NUL\SOH\b\NUL\NUL@\160@\NUL\NUL\NUL\STX\DLE\NULp\NUL\NUL\128\STX\NUL\b\EOT\b\NUL\NUL\NUL\SOH\STX\ETX\NUL\DC4\NUL\NUL\160\128\128\NUL\EOT\EOT\NUL \DLE"#

    mask637 :: Ptr Word
    mask637 :: Ptr Word
    mask637 = Addr# -> Ptr Word
    forall a. Addr# -> Ptr a
    Ptr Addr#
    "\ETX!\NUL\b\EOT\NUL\NUL\STX\SOH@\b\DC4\b\SOH@ \NUL\NUL\DLE\b\NULB\160@\CAN\NUL\STX\SOH\NUL\128@\NUL\DLE\STX\ENQB@\DLE\b\NUL\NUL\EOT\130\128\DLE(\DLE\STX\128@\NUL\NUL \DLE\NUL\134@\129\DLE\NUL\EOT\STX\NUL\NUL\129\NUL \EOT\n\132\NUL \DLE\NUL\NUL\b\EOT\NUL!\DLE"#

    mask703 :: Ptr Word
    mask703 :: Ptr Word
    mask703 = Addr# -> Ptr Word
    forall a. Addr# -> Ptr a
    Ptr Addr#
    "\ETX\t\NUL\140` \NUL\NUL\DC1\b\DLE\SOH\128\NUL\NUL&@\DLE\NUL\128\NUL\b\b\128\NUL\NUL\SOH!\DLE\DLE\NUL\SOH\f\n\STX`\NUL\ETX\SOH\NUL\f@\160\NUL\NUL\ENQ\STX0\NUL\128\192\NUL\ACK@P0\128\NUL\b\b\132\128\NUL\NUL\SOH\DLE\DLE\NUL\SOH\NUL\b\STXd\NUL\NUL\SOH\128\b\DLE\136\NUL\NUL\EOT\ACK1\NUL\144@"#

    -- -- not very discriminating, but cheap, so it's an overall gain
    -- cr512 :: V.Vector Bool
    -- cr512 = runST $ do
    --     ar <- MV.replicate 512 True
    --     let note s i
    --             | i < 512   = MV.unsafeWrite ar i False >> note s (i+s)
    --             | otherwise = return ()
    --     note 4 2
    --     note 8 4
    --     note 32 16
    --     note 64 32
    --     note 256 128
    --     MV.unsafeWrite ar 256 False
    --     V.unsafeFreeze ar

    -- -- Remainders modulo @3^3 * 31@
    -- cubeRes837 :: V.Vector Bool
    -- cubeRes837 = runST $ do
    --     ar <- MV.replicate 837 False
    --     let note 837 = return ()
    --         note k = MV.unsafeWrite ar ((k*k*k) `rem` 837) True >> note (k+1)
    --     note 0
    --     V.unsafeFreeze ar

    -- -- Remainders modulo @7^2 * 13@
    -- cubeRes637 :: V.Vector Bool
    -- cubeRes637 = runST $ do
    --     ar <- MV.replicate 637 False
    --     let note 637 = return ()
    --         note k = MV.unsafeWrite ar ((k*k*k) `rem` 637) True >> note (k+1)
    --     note 0
    --     V.unsafeFreeze ar

    -- -- Remainders modulo @19 * 37@
    -- cubeRes703 :: V.Vector Bool
    -- cubeRes703 = runST $ do
    --     ar <- MV.replicate 703 False
    --     let note 703 = return ()
    --         note k = MV.unsafeWrite ar ((k*k*k) `rem` 703) True >> note (k+1)
    --     note 0
    --     V.unsafeFreeze ar