{-# LANGUAGE ViewPatterns #-}
-- | Chinese Remainder for Rational numbers.
--
--   Since 0.4.0.0
module Algebra.Algorithms.ChineseRemainder
       ( recoverRat
       , rationalChineseRemainder
       ) where
import           Algebra.Instances        ()
import           AlgebraicPrelude
import           Data.List                (findIndices)
import           Numeric.Domain.Euclidean (euclid)
import           Numeric.Domain.Euclidean (chineseRemainder)
import qualified Prelude                  as P

-- | Recovers rational number from Z/pZ.
recoverRat :: Integer                    -- ^ Bound for numerator
           -> Integer                    -- ^ modulus
           -> Integer                    -- ^ integer corresponds to the rational number.
           -> Maybe (Fraction Integer)   -- ^ recovered rational number
recoverRat :: Integer -> Integer -> Integer -> Maybe (Fraction Integer)
recoverRat (Integer -> Integer
forall a. Num a => a -> a
abs -> Integer
k) Integer
m Integer
g
  | Integer
g Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
0 = Fraction Integer -> Maybe (Fraction Integer)
forall a. a -> Maybe a
Just Fraction Integer
0
  | Bool
otherwise =
  let ps :: [(Integer, Integer, Integer)]
ps = Integer -> Integer -> [(Integer, Integer, Integer)]
forall d. Euclidean d => d -> d -> [(d, d, d)]
euclid Integer
m Integer
g
      ixs :: [Int]
ixs = ((Integer, Integer, Integer) -> Bool)
-> [(Integer, Integer, Integer)] -> [Int]
forall a. (a -> Bool) -> [a] -> [Int]
findIndices (\(Integer
rj, Integer
_, Integer
_) -> Integer -> Integer
forall a. Num a => a -> a
P.abs Integer
rj Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
k) [(Integer, Integer, Integer)]
ps
  in if [Int] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [Int]
ixs
     then Maybe (Fraction Integer)
forall a. Maybe a
Nothing
     else
       let j :: Int
j = [Int] -> Int
forall a. [a] -> a
last [Int]
ixs
           (Integer
r, Integer
_, Integer
t)   = [(Integer, Integer, Integer)]
ps [(Integer, Integer, Integer)] -> Int -> (Integer, Integer, Integer)
forall a. [a] -> Int -> a
!! Int
j
           (Integer
r0,Integer
_ , Integer
t0) = [(Integer, Integer, Integer)]
ps [(Integer, Integer, Integer)] -> Int -> (Integer, Integer, Integer)
forall a. [a] -> Int -> a
!! (Int
j Int -> Int -> Int
forall r. Additive r => r -> r -> r
+ Int
1)
           q :: Integer
q | Int
j Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
== Int
0  = Integer
0
             | Bool
otherwise = [Integer] -> Integer
forall a. [a] -> a
head ([Integer] -> Integer) -> [Integer] -> Integer
forall a b. (a -> b) -> a -> b
$ (Integer -> Bool) -> [Integer] -> [Integer]
forall a. (a -> Bool) -> [a] -> [a]
filter (\Integer
v -> Integer
r0 Integer -> Integer -> Integer
forall r. Group r => r -> r -> r
- Integer
vInteger -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
*Integer
r Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
< Integer
k Bool -> Bool -> Bool
&& Integer
k Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
r0 Integer -> Integer -> Integer
forall r. Group r => r -> r -> r
- (Integer
vInteger -> Integer -> Integer
forall r. Group r => r -> r -> r
-Integer
1)Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
*Integer
r) [Integer
1..]
           (Integer
r', Integer
t') = (Integer
r0 Integer -> Integer -> Integer
forall r. Group r => r -> r -> r
- Integer
q Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Integer
r, Integer
t0 Integer -> Integer -> Integer
forall r. Group r => r -> r -> r
- Integer
q Integer -> Integer -> Integer
forall r. Multiplicative r => r -> r -> r
* Integer
t)
       in if Integer -> Integer -> Integer
forall d. GCDDomain d => d -> d -> d
gcd Integer
r Integer
t Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
1
          then Fraction Integer -> Maybe (Fraction Integer)
forall a. a -> Maybe a
Just (Integer
r Integer -> Integer -> Fraction Integer
forall d. GCDDomain d => d -> d -> Fraction d
% Integer
t)
          else if Integer -> Integer -> Integer
forall d. GCDDomain d => d -> d -> d
gcd Integer
r' Integer
t' Integer -> Integer -> Bool
forall a. Eq a => a -> a -> Bool
== Integer
1 Bool -> Bool -> Bool
&& Integer -> Integer
forall a. Num a => a -> a
abs Integer
t' Integer -> Integer -> Bool
forall a. Ord a => a -> a -> Bool
<= Integer
m Integer -> Integer -> Integer
forall d. Euclidean d => d -> d -> d
`quot` Integer
k
               then Fraction Integer -> Maybe (Fraction Integer)
forall a. a -> Maybe a
Just (Integer
r' Integer -> Integer -> Fraction Integer
forall d. GCDDomain d => d -> d -> Fraction d
% Integer
t')
               else Maybe (Fraction Integer)
forall a. Maybe a
Nothing

-- | Chinese Remainder for raional numbers.
rationalChineseRemainder :: Integer
                         -> [(Integer, Integer)]
                         -> Maybe (Fraction Integer)
rationalChineseRemainder :: Integer -> [(Integer, Integer)] -> Maybe (Fraction Integer)
rationalChineseRemainder Integer
k [(Integer, Integer)]
mvs =
  let m :: Integer
m = [Integer] -> Integer
forall (f :: * -> *) r. (Foldable f, Unital r) => f r -> r
product ([Integer] -> Integer) -> [Integer] -> Integer
forall a b. (a -> b) -> a -> b
$ ((Integer, Integer) -> Integer)
-> [(Integer, Integer)] -> [Integer]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
map (Integer, Integer) -> Integer
forall a b. (a, b) -> a
fst [(Integer, Integer)]
mvs
      g :: Integer
g = [(Integer, Integer)] -> Integer
forall r. Euclidean r => [(r, r)] -> r
chineseRemainder [(Integer, Integer)]
mvs
  in Integer -> Integer -> Integer -> Maybe (Fraction Integer)
recoverRat Integer
k Integer
m Integer
g