```{-# OPTIONS_GHC -fglasgow-exts -fallow-undecidable-instances #-}

module QIO.Shor where

import Data.Monoid as Monoid
import QIO.QIORandom
import QIO.QioSyn
import QIO.Qdata
import QIO.Qio
import QIO.QExamples
import QIO.QArith
import QIO.Qft
import System.Time

qftI :: QInt -> U
qftI (QInt i) = urev (qft i)

shorU :: QInt -> QInt -> Int -> Int -> U
shorU k i1 x n = hadamardsI k
`mappend`
modExp n x k i1
`mappend`
qftI k

shor :: Int -> Int -> QIO Int
shor x n = do i0 <- mkQ 0
i1 <- mkQ 1
applyU (shorU i0 i1 x n)
p <- measQ i0
return p

period :: Int -> Int -> Int
period m q = r where (_,r) = reduce (m,q)

factor :: Int -> QIO (Int,Int)
factor n | even n = return (2,2)
| otherwise = do x <- rand_coprime n
a <- shor x n
let xa = x^(half a)
in if   odd a || xa == (n-1) `mod` n || a == 0
then factor n
else  return (gcd (xa+1) n,gcd (xa-1) n)
--this function can only be run too, for similar reasons to the rand_co'
--function below

runTime :: QIO a -> IO a
runTime a = do start <- getClockTime
result <- run a
stop <- getClockTime
putStr ("The total time taken was " ++
(timeDiffToString (diffClockTimes stop start) ++ "\n"))
return result

factorV' :: Int -> IO (Int,Int)
factorV' n | even n = return (2,2)
| otherwise = do start <- getClockTime
putStr ("Started at " ++ (show start) ++ "\n")
x <- run (rand_coprime n)
putStr ("Calling \"shor " ++ show x ++ " " ++ show n ++ "\"\n")
a <- run (shor x n)
stop <- getClockTime
putStr ("Shor took " ++ (timeDiffToString (diffClockTimes stop start)) ++ "\n")
putStr ("period a = " ++ show a)
let xa = x^(half a)
in do putStr (", giving xa = " ++ show xa ++ "\n")
if odd a || xa == (n-1) `mod` n || (gcd (xa+1) n,gcd (xa-1) n) == (1,n) || (gcd (xa+1) n,gcd (xa-1) n) == (n,1) || (gcd (xa+1) n,gcd (xa-1) n) == (1,1)
then do putStr "Recalling factorV\n"
factorV' n
else do putStr "Result: "
return (gcd (xa+1) n,gcd (xa-1) n)

factorV :: Int -> IO ()
factorV n = do start <- getClockTime
(a,b) <- factorV' n
stop <- getClockTime
putStr ( "Factors of "
++ (show n)
++ " include "
++ (show a)
++ " and "
++ (show b)
++ ".\n The total time taken was "
++ (timeDiffToString (diffClockTimes stop start) ++ "\n"))

rand_co' :: Int -> QIO Int
rand_co' n = do x <- randomQIO (2,n)
if gcd x n == 1 then return x else rand_co' n
--simulating this (with the sim function) gives rise to infinite paths in
--the computation, e.g. each path where gcd x n /= 1. However, this function
--can still be run (with the run function) always returning a single value.

rand_coprime :: Int -> QIO Int
rand_coprime n = do x <- randomQIO (0,(length cps)-1)
return (cps!!x)
where cps = [x | x <- [0..n], gcd x n == 1]

half :: Int -> Int
half x = floor (fromIntegral x/2.0)

reduce :: (Int,Int) -> (Int,Int)
reduce (x,y) = if g == 1 then (x,y) else (floor ((fromIntegral x)/(fromIntegral g)),floor ((fromIntegral y)/(fromIntegral g)))
where g = gcd x y

```