1{-# LANGUAGE BangPatterns #-}
2-- |
3-- Module      : Crypto.Number.ModArithmetic
4-- License     : BSD-style
5-- Maintainer  : Vincent Hanquez <vincent@snarc.org>
6-- Stability   : experimental
7-- Portability : Good
8
9module Crypto.Number.ModArithmetic
10    (
11    -- * Exponentiation
12      expSafe
13    , expFast
14    -- * Inverse computing
15    , inverse
16    , inverseCoprimes
17    , inverseFermat
18    -- * Squares
19    , jacobi
20    , squareRoot
21    ) where
22
23import Control.Exception (throw, Exception)
24import Crypto.Number.Basic
25import Crypto.Number.Compat
26
27-- | Raised when two numbers are supposed to be coprimes but are not.
28data CoprimesAssertionError = CoprimesAssertionError
29    deriving (Show)
30
31instance Exception CoprimesAssertionError
32
33-- | Compute the modular exponentiation of base^exponent using
34-- algorithms design to avoid side channels and timing measurement
35--
36-- Modulo need to be odd otherwise the normal fast modular exponentiation
37-- is used.
38--
39-- When used with integer-simple, this function is not different
40-- from expFast, and thus provide the same unstudied and dubious
41-- timing and side channels claims.
42--
43-- Before GHC 8.4.2, powModSecInteger is missing from integer-gmp,
44-- so expSafe has the same security as expFast.
45expSafe :: Integer -- ^ base
46        -> Integer -- ^ exponent
47        -> Integer -- ^ modulo
48        -> Integer -- ^ result
49expSafe b e m
50    | odd m     = gmpPowModSecInteger b e m `onGmpUnsupported`
51                  (gmpPowModInteger b e m   `onGmpUnsupported`
52                  exponentiation b e m)
53    | otherwise = gmpPowModInteger b e m    `onGmpUnsupported`
54                  exponentiation b e m
55
56-- | Compute the modular exponentiation of base^exponent using
57-- the fastest algorithm without any consideration for
58-- hiding parameters.
59--
60-- Use this function when all the parameters are public,
61-- otherwise 'expSafe' should be preferred.
62expFast :: Integer -- ^ base
63        -> Integer -- ^ exponent
64        -> Integer -- ^ modulo
65        -> Integer -- ^ result
66expFast b e m = gmpPowModInteger b e m `onGmpUnsupported` exponentiation b e m
67
68-- | @exponentiation@ computes modular exponentiation as /b^e mod m/
69-- using repetitive squaring.
70exponentiation :: Integer -> Integer -> Integer -> Integer
71exponentiation b e m
72    | b == 1    = b
73    | e == 0    = 1
74    | e == 1    = b `mod` m
75    | even e    = let p = exponentiation b (e `div` 2) m `mod` m
76                   in (p^(2::Integer)) `mod` m
77    | otherwise = (b * exponentiation b (e-1) m) `mod` m
78
79-- | @inverse@ computes the modular inverse as in /g^(-1) mod m/.
80inverse :: Integer -> Integer -> Maybe Integer
81inverse g m = gmpInverse g m `onGmpUnsupported` v
82  where
83    v
84        | d > 1     = Nothing
85        | otherwise = Just (x `mod` m)
86    (x,_,d) = gcde g m
87
88-- | Compute the modular inverse of two coprime numbers.
89-- This is equivalent to inverse except that the result
90-- is known to exists.
91--
92-- If the numbers are not defined as coprime, this function
93-- will raise a 'CoprimesAssertionError'.
94inverseCoprimes :: Integer -> Integer -> Integer
95inverseCoprimes g m =
96    case inverse g m of
97        Nothing -> throw CoprimesAssertionError
98        Just i  -> i
99
100-- | Computes the Jacobi symbol (a/n).
101-- 0 ≤ a < n; n ≥ 3 and odd.
102--
103-- The Legendre and Jacobi symbols are indistinguishable exactly when the
104-- lower argument is an odd prime, in which case they have the same value.
105--
106-- See algorithm 2.149 in "Handbook of Applied Cryptography" by Alfred J. Menezes et al.
107jacobi :: Integer -> Integer -> Maybe Integer
108jacobi a n
109    | n < 3 || even n  = Nothing
110    | a == 0 || a == 1 = Just a
111    | n <= a           = jacobi (a `mod` n) n
112    | a < 0            =
113      let b = if n `mod` 4 == 1 then 1 else -1
114       in fmap (*b) (jacobi (-a) n)
115    | otherwise        =
116      let (e, a1) = asPowerOf2AndOdd a
117          nMod8   = n `mod` 8
118          nMod4   = n `mod` 4
119          a1Mod4  = a1 `mod` 4
120          s'      = if even e || nMod8 == 1 || nMod8 == 7 then 1 else -1
121          s       = if nMod4 == 3 && a1Mod4 == 3 then -s' else s'
122          n1      = n `mod` a1
123       in if a1 == 1 then Just s
124          else fmap (*s) (jacobi n1 a1)
125
126-- | Modular inverse using Fermat's little theorem.  This works only when
127-- the modulus is prime but avoids side channels like in 'expSafe'.
128inverseFermat :: Integer -> Integer -> Integer
129inverseFermat g p = expSafe g (p - 2) p
130
131-- | Raised when the assumption about the modulus is invalid.
132data ModulusAssertionError = ModulusAssertionError
133    deriving (Show)
134
135instance Exception ModulusAssertionError
136
137-- | Modular square root of @g@ modulo a prime @p@.
138--
139-- If the modulus is found not to be prime, the function will raise a
140-- 'ModulusAssertionError'.
141--
142-- This implementation is variable time and should be used with public
143-- parameters only.
144squareRoot :: Integer -> Integer -> Maybe Integer
145squareRoot p
146    | p < 2     = throw ModulusAssertionError
147    | otherwise =
148        case p `divMod` 8 of
149           (v, 3) -> method1 (2 * v + 1)
150           (v, 7) -> method1 (2 * v + 2)
151           (u, 5) -> method2 u
152           (_, 1) -> tonelliShanks p
153           (0, 2) -> \a -> Just (if even a then 0 else 1)
154           _      -> throw ModulusAssertionError
155
156  where
157    x `eqMod` y = (x - y) `mod` p == 0
158
159    validate g y | (y * y) `eqMod` g = Just y
160                 | otherwise         = Nothing
161
162    -- p == 4u + 3 and u' == u + 1
163    method1 u' g =
164        let y = expFast g u' p
165         in validate g y
166
167    -- p == 8u + 5
168    method2 u g =
169        let gamma = expFast (2 * g) u p
170            g_gamma = g * gamma
171            i = (2 * g_gamma * gamma) `mod` p
172            y = (g_gamma * (i - 1)) `mod` p
173         in validate g y
174
175tonelliShanks :: Integer -> Integer -> Maybe Integer
176tonelliShanks p a
177    | aa == 0   = Just 0
178    | otherwise =
179        case expFast aa p2 p of
180            b | b == p1   -> Nothing
181              | b == 1    -> Just $ go (expFast aa ((s + 1) `div` 2) p)
182                                       (expFast aa s p)
183                                       (expFast n  s p)
184                                       e
185              | otherwise -> throw ModulusAssertionError
186  where
187    aa = a `mod` p
188    p1 = p - 1
189    p2 = p1 `div` 2
190    n  = findN 2
191
192    x `mul` y = (x * y) `mod` p
193
194    pow2m 0 x = x
195    pow2m i x = pow2m (i - 1) (x `mul` x)
196
197    (e, s) = asPowerOf2AndOdd p1
198
199    -- find a quadratic non-residue
200    findN i
201        | expFast i p2 p == p1 = i
202        | otherwise            = findN (i + 1)
203
204    -- find m such that b^(2^m) == 1 (mod p)
205    findM b i
206        | b == 1    = i
207        | otherwise = findM (b `mul` b) (i + 1)
208
209    go !x b g !r
210        | b == 1    = x
211        | otherwise =
212            let r' = findM b 0
213                z = pow2m (r - r' - 1) g
214                x' = x `mul` z
215                b' = b `mul` g'
216                g' = z `mul` z
217             in go x' b' g' r'
218