-- | Triangular numbers and various helper functions.
--
-- Main use is for linearization of triangular array indexing.
-- 
-- Triangular numbers:
-- @
-- T_n = Σ_{k=1)^n k = 1 + 2 + 3 + ... + n =
--
-- n * (n+1) / 2 = (n+1) `choose` 2
-- @
--
--

module Math.TriangularNumbers where



-- | Triangular numbers.
--
-- https://oeis.org/A000217

triangularNumber :: Int -> Int
triangularNumber :: Int -> Int
triangularNumber Int
x = (Int
x Int -> Int -> Int
forall a. Num a => a -> a -> a
* (Int
xInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)) Int -> Int -> Int
forall a. Integral a => a -> a -> a
`quot` Int
2
{-# INLINE triangularNumber #-}

-- | Size of an upper triangle starting at 'i' and ending at 'j'. "(0,N)" what
-- be the normal thing to use.

linearizeUppertri :: (Int,Int) -> Int
linearizeUppertri :: (Int, Int) -> Int
linearizeUppertri (Int
i,Int
j) = Int -> Int
triangularNumber (Int -> Int) -> Int -> Int
forall a b. (a -> b) -> a -> b
$ Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
iInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1
{-# INLINE linearizeUppertri #-}

-- | Subword indexing. Given the longest subword and the current subword,
-- calculate a linear index "[0,..]". "(l,n)" in this case means "l"ower bound,
-- length "n". And "(i,j)" is the normal index.
--
-- @
-- 0 1 2 3    <- j = ...
--
-- 0 1 2 3    i=0
-- _ 4 5 6    i=1
-- _ _ 7 8    i=2
--       9    i=3
--
-- i=2, j=3  -> (4+1) * i - tri i + j
--
-- _
-- _ _  the triangular number to subtract.
-- @

toLinear :: Int -> (Int,Int) -> Int
toLinear :: Int -> (Int, Int) -> Int
toLinear Int
n (Int
i,Int
j) = Int -> (Int, Int) -> Int
adr Int
n (Int
i,Int
j)
  where
    adr :: Int -> (Int, Int) -> Int
adr Int
n (Int
i,Int
j) = (Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1)Int -> Int -> Int
forall a. Num a => a -> a -> a
*Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int -> Int
triangularNumber Int
i Int -> Int -> Int
forall a. Num a => a -> a -> a
+ Int
j
    {-# Inline adr #-}
{-# INLINE toLinear #-}



-- | Linear index to paired.
--
-- We have indices in @[0,N]@, and linear index @k@.
--
-- @
-- (N+1)*i - (i*(i+1)/2) + j == K
-- @

fromLinear :: Int -> Int -> (Int,Int)
fromLinear :: Int -> Int -> (Int, Int)
fromLinear Int
n' Int
k' = (Int
i,Int
j)
  where ll :: Double
ll = (Double
2Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
nDouble -> Double -> Double
forall a. Num a => a -> a -> a
+Double
1) Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
2
        rr :: Double
rr = Double -> Double
forall a. Floating a => a -> a
sqrt (Double -> Double) -> Double -> Double
forall a b. (a -> b) -> a -> b
$ ((Double
2Double -> Double -> Double
forall a. Num a => a -> a -> a
*(Double
nDouble -> Double -> Double
forall a. Num a => a -> a -> a
+Double
1)Double -> Double -> Double
forall a. Num a => a -> a -> a
+Double
1) Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ Double
2)Double -> Integer -> Double
forall a b. (Num a, Integral b) => a -> b -> a
^Integer
2 Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
2Double -> Double -> Double
forall a. Num a => a -> a -> a
*Double
k
        n :: Double
n  = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n'
        k :: Double
k  = Int -> Double
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
k'
        i :: Int
i  = Double -> Int
forall a b. (RealFrac a, Integral b) => a -> b
floor (Double -> Int) -> Double -> Int
forall a b. (a -> b) -> a -> b
$ Double
ll Double -> Double -> Double
forall a. Num a => a -> a -> a
- Double
rr Double -> Double -> Double
forall a. Num a => a -> a -> a
+ Double
1
        j :: Int
j  = Int
k' Int -> Int -> Int
forall a. Num a => a -> a -> a
- Int -> (Int, Int) -> Int
toLinear Int
n' (Int
i,Int
0)
{-# Inline fromLinear #-}