{-# LANGUAGE FlexibleInstances #-}
{-# LANGUAGE OverloadedLists #-}

{-| 
Module      : Internal.Gates
Description : Gate library internals
Stability   : experimental

Internal matrix operations
-}
module Lib.Internal.Gates where

import Lib.QM ( getState, get, put, QM, QState(QState), QBit(..), Ix)
import Numeric.LinearAlgebra
    ( Complex(..),
      C,
      (#>),
      (><),
      dispcf,
      ident,
      kronecker,
      Matrix,
      Linear(scale) )

instance {-# OVERLAPS#-} Show (Matrix C) where
  show :: Matrix C -> String
show mx :: Matrix C
mx = Int -> Matrix C -> String
dispcf 3 Matrix C
mx

-- | The imaginary unit
i :: Complex Double
i :: C
i = 0 Double -> Double -> C
forall a. a -> a -> Complex a
:+ 1

applyParallel :: Matrix C -> Matrix C -> Matrix C
applyParallel :: Matrix C -> Matrix C -> Matrix C
applyParallel = Matrix C -> Matrix C -> Matrix C
forall t. Product t => Matrix t -> Matrix t -> Matrix t
kronecker

-- | Apply gate to the current quantum state
applyGate :: Matrix C -> QM ()
applyGate :: Matrix C -> QM ()
applyGate g :: Matrix C
g = do
    (QState v :: Vector C
v) <- QM QState
get
    QState -> QM ()
put (QState -> QM ()) -> QState -> QM ()
forall a b. (a -> b) -> a -> b
$ Vector C -> QState
QState (Vector C -> QState) -> Vector C -> QState
forall a b. (a -> b) -> a -> b
$ Matrix C
g Matrix C -> Vector C -> Vector C
forall t. Numeric t => Matrix t -> Vector t -> Vector t
#> Vector C
v

-- | Changes an element at an index in a list.
--  index 0 will change the first element.
changeAt :: a -> Int -> [a] -> [a]
changeAt :: a -> Int -> [a] -> [a]
changeAt x :: a
x index :: Int
index [] = String -> [a]
forall a. HasCallStack => String -> a
error "changeAt: Can't change an element in an empty list"
changeAt x :: a
x 0 (_:ys :: [a]
ys) = a
xa -> [a] -> [a]
forall a. a -> [a] -> [a]
:[a]
ys
changeAt x :: a
x index :: Int
index (y :: a
y:ys :: [a]
ys) = a
y a -> [a] -> [a]
forall a. a -> [a] -> [a]
: a -> Int -> [a] -> [a]
forall a. a -> Int -> [a] -> [a]
changeAt a
x (Int
index Int -> Int -> Int
forall a. Num a => a -> a -> a
- 1)  [a]
ys

-- | Apply a 2x2 gate, to a specific qubit.
--
-- It will update the qstate. 
runGate :: Matrix C -> (QBit -> QM QBit)
runGate :: Matrix C -> QBit -> QM QBit
runGate g :: Matrix C
g x :: QBit
x = do
    (state :: QState
state, size :: Int
size) <- QM (QState, Int)
getState
    let ids :: [Matrix C]
ids = Int -> Matrix C -> [Matrix C]
forall a. Int -> a -> [a]
replicate Int
size (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident 2)
    let list :: [Matrix C]
list = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
g (QBit -> Int
link QBit
x) [Matrix C]
ids
    let m :: Matrix C
m = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
list
    Matrix C -> QM ()
applyGate Matrix C
m
    QBit -> QM QBit
forall (m :: * -> *) a. Monad m => a -> m a
return QBit
x

-- run specified gates in parallel
-- TODO: should work, but not fully tested. does NOT work with controlled gates
parallel :: Int -> [(Matrix C, QBit)] -> Matrix C
parallel :: Int -> [(Matrix C, QBit)] -> Matrix C
parallel size :: Int
size as :: [(Matrix C, QBit)]
as = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
list
  where list :: [Matrix C]
list = ((Matrix C, QBit) -> [Matrix C] -> [Matrix C])
-> [Matrix C] -> [(Matrix C, QBit)] -> [Matrix C]
forall (t :: * -> *) a b.
Foldable t =>
(a -> b -> b) -> b -> t a -> b
foldr (Matrix C, QBit) -> [Matrix C] -> [Matrix C]
apply (Int -> Matrix C -> [Matrix C]
forall a. Int -> a -> [a]
replicate Int
size (Matrix C -> [Matrix C]) -> Matrix C -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident 2) [(Matrix C, QBit)]
as
        apply :: (Matrix C, QBit) -> [Matrix C] -> [Matrix C]
        apply :: (Matrix C, QBit) -> [Matrix C] -> [Matrix C]
apply (mx :: Matrix C
mx, Ptr q :: Int
q) nx :: [Matrix C]
nx = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
mx Int
q [Matrix C]
nx

-- | Produce matrix running a gate controlled by another bit
controlMatrix :: Int -> QBit -> QBit -> Matrix C -> Matrix C
controlMatrix :: Int -> QBit -> QBit -> Matrix C -> Matrix C
controlMatrix size :: Int
size (Ptr c :: Int
c) (Ptr t :: Int
t) g :: Matrix C
g = Matrix C
fl Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C
fr
  where idsl :: [Matrix C]
idsl = Int -> Matrix C -> [Matrix C]
forall a. Int -> a -> [a]
replicate Int
size (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident 2)
        idsr :: [Matrix C]
idsr = Int -> Matrix C -> [Matrix C]
forall a. Int -> a -> [a]
replicate Int
size (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident 2)
        l :: [Matrix C]
l = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj0 Int
c [Matrix C]
idsl
        rc :: [Matrix C]
rc = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj1 Int
c [Matrix C]
idsr
        r :: [Matrix C]
r = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
g Int
t [Matrix C]
rc
        fl :: Matrix C
fl = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
l
        fr :: Matrix C
fr = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
r
-- | Produce a matrix running a gate controlled by two other bits
ccontrolMatrix :: Int -> QBit -> QBit -> QBit -> Matrix C -> Matrix C
ccontrolMatrix :: Int -> QBit -> QBit -> QBit -> Matrix C -> Matrix C
ccontrolMatrix size :: Int
size (Ptr c1 :: Int
c1) (Ptr c2 :: Int
c2) (Ptr t :: Int
t) g :: Matrix C
g = Matrix C
f00 Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C
f01 Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C
f10 Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C
f11
  where ids :: [Matrix C]
ids = Int -> Matrix C -> [Matrix C]
forall a. Int -> a -> [a]
replicate Int
size (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident 2)
        m00c :: [Matrix C]
m00c = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj0 Int
c2 ([Matrix C] -> [Matrix C]) -> [Matrix C] -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj0 Int
c1 [Matrix C]
ids
        m01c :: [Matrix C]
m01c = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj1 Int
c2 ([Matrix C] -> [Matrix C]) -> [Matrix C] -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj0 Int
c1 [Matrix C]
ids
        m10c :: [Matrix C]
m10c = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj0 Int
c2 ([Matrix C] -> [Matrix C]) -> [Matrix C] -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj1 Int
c1 [Matrix C]
ids
        m11c :: [Matrix C]
m11c = Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj1 Int
c2 ([Matrix C] -> [Matrix C]) -> [Matrix C] -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
proj1 Int
c1 ([Matrix C] -> [Matrix C]) -> [Matrix C] -> [Matrix C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Int -> [Matrix C] -> [Matrix C]
forall a. a -> Int -> [a] -> [a]
changeAt Matrix C
g Int
t [Matrix C]
ids
        f00 :: Matrix C
f00  = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
m00c
        f01 :: Matrix C
f01  = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
m01c
        f10 :: Matrix C
f10  = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
m10c
        f11 :: Matrix C
f11  = (Matrix C -> Matrix C -> Matrix C) -> [Matrix C] -> Matrix C
forall (t :: * -> *) a. Foldable t => (a -> a -> a) -> t a -> a
foldr1 Matrix C -> Matrix C -> Matrix C
applyParallel [Matrix C]
m11c

-- | Quantum fourier transform matrix
qftMatrix :: Int -> Matrix C
qftMatrix :: Int -> Matrix C
qftMatrix n :: Int
n = (1 Matrix C -> Matrix C -> Matrix C
forall a. Fractional a => a -> a -> a
/ Matrix C -> Matrix C
forall a. Floating a => a -> a
sqrt (Int -> Matrix C
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n)) Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
* (Int
n Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< Int
n)
  [ C
ωC -> Int -> C
forall a b. (Num a, Integral b) => a -> b -> a
^(Int
jInt -> Int -> Int
forall a. Num a => a -> a -> a
*Int
k) | Int
j <- [0..Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-1], Int
k <- [0..Int
nInt -> Int -> Int
forall a. Num a => a -> a -> a
-1] ]
  where ω :: C
ω = C -> C
forall a. Floating a => a -> a
exp ((2C -> C -> C
forall a. Num a => a -> a -> a
*C
forall a. Floating a => a
piC -> C -> C
forall a. Num a => a -> a -> a
*C
i) C -> C -> C
forall a. Fractional a => a -> a -> a
/ Int -> C
forall a b. (Integral a, Num b) => a -> b
fromIntegral Int
n)

notAdjacent :: [Ix] -> Bool
notAdjacent :: [Int] -> Bool
notAdjacent [a :: Item [Int]
a]      = Bool
False
notAdjacent [a :: Item [Int]
a, b :: Item [Int]
b]   = Int
Item [Int]
bInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
Item [Int]
a Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= 1
notAdjacent (a :: Int
a:b :: Int
b:as :: [Int]
as) = Int
bInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
a Int -> Int -> Bool
forall a. Eq a => a -> a -> Bool
/= 1 Bool -> Bool -> Bool
|| [Int] -> Bool
notAdjacent (Int
bInt -> [Int] -> [Int]
forall a. a -> [a] -> [a]
:[Int]
as)

-- | Projection of the zero basis vector
proj0 :: Matrix C
proj0 :: Matrix C
proj0 = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 1 , 0
  , 0 , 0 ]

-- | Projection of the one basis vector
proj1 :: Matrix C
proj1 :: Matrix C
proj1 = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 0 , 0
  , 0 , 1 ]

-- | Hadamard matrix
hmat :: Matrix C
hmat :: Matrix C
hmat = C -> Matrix C -> Matrix C
forall t (c :: * -> *). Linear t c => t -> c t -> c t
scale (C -> C
forall a. Floating a => a -> a
sqrt 0.5) (Matrix C -> Matrix C) -> Matrix C -> Matrix C
forall a b. (a -> b) -> a -> b
$ (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
    [ 1 ,  1
    , 1 , -1 ]

-- | CNOT matrix
cmat :: Matrix C
cmat :: Matrix C
cmat = (4 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 4)
  [ 1, 0, 0, 0
  , 0, 1, 0, 0
  , 0, 0, 0, 1
  , 0, 0, 1, 0 ]

-- | Generic phase matrix, takes in phase change as radians
phasemat :: Double -> Matrix C
phasemat :: Double -> Matrix C
phasemat r :: Double
r = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 1, 0
  , 0, C
Item [C]
p ]
  where p :: C
p = C -> C
forall a. Floating a => a -> a
exp (C
iC -> C -> C
forall a. Num a => a -> a -> a
*(Double
r Double -> Double -> C
forall a. a -> a -> Complex a
:+ 0))

-- | PauliX matrix
pXmat :: Matrix C
pXmat :: Matrix C
pXmat = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 0, 1
  , 1, 0 ]

-- | PauliY matrix
pYmat :: Matrix C
pYmat :: Matrix C
pYmat = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 0, -C
i
  , C
Item [C]
i,  0 ]

-- | PauliZ matrix
pZmat :: Matrix C
pZmat :: Matrix C
pZmat = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 1 ,  0
  , 0 , -1 ]

-- | Identity matrix
idmat :: Matrix C
idmat :: Matrix C
idmat = (2 Int -> Int -> [C] -> Matrix C
forall a. Storable a => Int -> Int -> [a] -> Matrix a
>< 2)
  [ 1 , 0
  , 0 , 1 ]