Shank Tonelli Algorithm

Shanks tonelli algorithm is for finding square roots modulo prime number.
It solves where n is quadratic residue (mod p ) and p is odd prime. Implementation in Haskell.

import Data.List
--Calculates x^2=n (mod p) where p is prime n is quadratic residue. 

legendrE::Integer->Integer->Integer
legendrE a p = modpoW a (div (p-1) 2) p

modpoW::Integer->Integer->Integer->Integer
modpoW a d n= fuN a d where
        fuN a 1 = mod a n
        fuN a 2 = mod (a^2) n
        fuN a d = let
                        p=fuN (mod (a^2) n) (div d 2)
                        q=if odd d then mod (a*p) n else p
                  in mod q n

shankS::Integer->Integer->Integer
shankS n p
    | legendrE n p /= 1 = error "non quadratic residude"
    | mod p 4 == 3 = n^(div (p+1) 4)
    |otherwise = 
      let 
	 q = until (\x->mod x 2 /= 0) (\x->div x 2) $ pred p
	 s = until (\x->q*2^x == pred p) (+1) 0
	 z = until (\x->legendrE x p == pred p) (+1) 2
	 c = modpoW z q p
	 r = modpoW n (div (q+1) 2) p
	 t = modpoW n q p
	 m = s
      in tonellI p c r t m
	
tonellI::Integer->Integer->Integer->Integer->Integer->Integer
tonellI p c r t m = 
    if mod t p == 1 then r 
    else let 
	   i =repsqR (t^2) p 1 --until (\x->modpoW t (2^x) p == 1) (+1) 1
	   b = modpoW c (2^(m-i-1)) p
	   r' = mod (r*b) p
	   t' = mod (t*b^2) p
	   c' = mod (b^2) p
	   m' = i
	  in tonellI p c' r' t' m'
           

repsqR t p cnt| mod t p == 1 =cnt
	      | otherwise = repsqR (t^2) p (cnt+1)

main::IO()
main=do
	n<-fmap read getLine
	p<-fmap read getLine;
	print $ shankS n p

About these ads