Higher sigma operators: basic building blocks for quantum observables

Sigmas.hs

-- Sigmas Module ~ Higher Sigma operators and their operator algebra
-- Cliff Harvey
-- February 2014

module Sigmas where

import Control.Applicative

-- 4 basic complex phases: 1, -1, i, -i
data Phase = P1 | M1 | Pi | Mi deriving Eq

instance Show Phase where
  show P1 = "+1"
  show M1 = "-1"
  show Pi = "+i"
  show Mi = "-i"

-- Phase multiply
multP :: Phase -> Phase -> Phase
multP P1 x = x
multP M1 M1 = P1
multP M1 Pi = Mi
multP M1 Mi = Pi
multP Pi Pi = M1
multP Pi Mi = P1
multP Mi Mi = M1 -- End of commutatively-complete pattern list,
multP y z = multP z y -- reflect to complete definition.


-- 4 Pauli operators including identity 
data Pauli = PI | PX | PY | PZ deriving Eq

instance Show Pauli where 
  show PI = "I"
  show PX = "σx"
  show PY = "σy"
  show PZ = "σz"
  
-- Sigma operator consists of a Pauli with a complex phase  
data Sigma = Sigma Phase Pauli deriving Eq

instance Show Sigma where
  show (Sigma pa m) | pa == P1 = (show m)
  		            | pa == M1 = ("-" ++ show m)
  		            | pa == Pi = ("i" ++ show m)
                    | pa == Mi = ("-i" ++ show m)

-- Sigma multiply
multS :: Sigma -> Sigma -> Sigma
multS (Sigma pa a) (Sigma pb b) | a == b = Sigma (multP pa pb) PI -- When the Pauli operators multiply to identity
				                | a == PI = Sigma (multP pa pb) b -- When one of the two Pauli's is the identity
				                | (a==PY && b==PZ) = Sigma (pa `multP` pb `multP` Pi) PX 
				                | (a==PZ && b==PX) = Sigma (pa `multP` pb `multP` Pi) PY
				                | (a==PX && b==PY) = Sigma (pa `multP` pb `multP` Pi) PZ 
				                | otherwise = multS (Sigma pb b) (Sigma (pa `multP` M1) a) -- anti-reflect to complete definition

{- Physics note: The 3 main sigma operators' multiplication table is the 3-dimensional epsilon symbol times i:

σ_i σ_j = i ε^(ijk) σ_k 

So, in other words, σ_x σ_y = i σ_z, and σ_y σ_x = -i σ_x, with all other nonzero products given by cyclicly permuting the indices on that formula.

Consequently they represent the algebra of rotations in 3-dimensional space.
More importantly for our purposes, they span the space of possible observations that may be performed on the most basic
quantum system, a single quantum bit (or qubit) whose classic incarnation is as the two spin states of a spin-1/2 particle.
Any even-dimensional quantum system, and its operator algebra, can be built up out of these. -}


-- Higher tensor product space of simga operators. Always just one overall phase
data HigherSigma = Only Phase | More Pauli HigherSigma deriving Eq

-- HigherSigma represents a generator of the Lie algebra of quantum observables on the state space of multiple qubits.
-- Example: iσx ⊗ σy ⊗ σz ~ (More PX (More PY (More PZ (Only Pi))))

instance Show HigherSigma where
  show (Only p) = show p
  show (More sm (Only pa)) = show (Sigma pa sm)
  show (More sm hs) = (prefix $ getPhase hs) ++ (show sm) ++ " ⊗ " ++ showx hs
    where showx (More sn (Only pn)) = show sn -- showx ~ show without sign prefix
          showx (More sn hq) = (show sn) ++ " ⊗ " ++ (showx hq)
          prefix P1 = ""
          prefix Pi = "i"
          prefix M1 = "-"
          prefix Mi = "-i"


-- Some utility functions
getPhase (Only p) = p
getPhase (More ps hs) = getPhase hs

sigPhase (Sigma p s) = p

getDepth :: HigherSigma -> Integer
getDepth (Only p) = 0
getDepth (More s h) = 1 + getDepth h

lift :: Sigma -> HigherSigma
lift (Sigma p s) = More s (Only p)

liftPauli :: Pauli -> HigherSigma
liftPauli s = (More s (Only P1))

getPurePauli :: Sigma -> Pauli
getPurePauli (Sigma px sx) = sx


-- Sigma operator composition. If unequal tensor depth, fill in extra spaces with I
multH :: HigherSigma -> HigherSigma -> HigherSigma
multH (Only p) (Only q) = Only (p `multP` q)
multH (Only p) (More s h) = More s (Only p `multH` h)
multH (More s h) (Only p) = More s (Only p `multH` h)
multH (More sa ha) (More sb hb) = More (getPurePauli subProduct) (ha `multH` (hb `multH` (Only (sigPhase subProduct))))
  where subProduct = (Sigma P1 sa) `multS` (Sigma P1 sb)

-- Tensor multiplication 
otimes :: HigherSigma -> HigherSigma -> HigherSigma
otimes (Only pa) (Only pb) = Only (pa `multP` pb)
otimes (Only p) (More s h) = More s (Only p `multH` h)
otimes (More s h) (Only p) = More s (Only p `multH` h)
otimes (More sa (Only pa)) (More sb (Only pb)) = More sa (More sb (Only (pa `multP` pb)))
otimes (More sa ha) (More sb (Only pb)) = (More sa (ha `otimes` (More sb (Only pb))))
otimes ha (More sb hb) = (ha `otimes` (More sb (Only P1))) `otimes` hb


-- Some convenient shorthands:
sx = (More PX (Only P1))
sy = (More PY (Only P1))
sz = (More PZ (Only P1))
sI = (More PI (Only P1))
isx = (More PX (Only Pi))
isy = (More PY (Only Pi))
isz = (More PZ (Only Pi))
isI = (More PI (Only Pi))
i = (Only Pi)
p1 = (Only P1)
m1 = (Only M1)
mi = (Only Mi)
x = multH
ox = otimes

--Tensor exponent
opower :: HigherSigma -> Integer -> HigherSigma
opower h 0 = sI
opower h 1 = h
opower h n = h `otimes` (h `opower` (n-1))

-- Exclusive OR
xor :: Bool -> Bool -> Bool
xor True tf = not tf
xor False tf = tf

-- Binary commutation ~ True iff [σ1,σ2] = σ1σ2 - σ2σ1 = 0
comm :: HigherSigma -> HigherSigma -> Bool
comm (Only p) _ = True -- If either is just a phase, then it obviously commutes
comm _ (Only p) = True
comm (More PI (Only pa)) (More _ (Only pb)) = True -- For single Pauli's, if either is identity then yes it obviously commutes
comm (More _ (Only pa)) (More PI (Only pb)) = True
comm (More sa (Only pa)) (More sb (Only pb)) = (sa == sb) -- For single Pauli's other than identity, they commute iff sa == sb
comm (More sa ha) (More sb hb) = (comm (More sa (Only P1)) (More sb (Only P1))) `xor` not (comm ha hb) 
--  \_--> For higher sigma operators, they commute iff the subcommutator of a two bit subspace vanishes XOR
--          if the subcommutator from the removal of that particular one-qubit subspace vanishes


-- matrify - Create basic matrix representations for sigmas out of Maybe Phase's. See liftPhase utility function below
matrify :: HigherSigma -> [[Maybe Phase]]
matrify (Only p) = [[Just p]]
matrify (More PI (Only p)) = (map $ map $ fmap $ multP p) [[Just P1, Nothing], [Nothing, Just P1]]
matrify (More PX (Only p)) = (map $ map $ fmap $ multP p) [[Nothing, Just P1], [Just P1, Nothing]]
matrify (More PY (Only p)) = (map $ map $ fmap $ multP p) [[Nothing, Just Pi], [Just Mi, Nothing]]
matrify (More PZ (Only p)) = (map $ map $ fmap $ multP p) [[Just P1, Nothing], [Nothing, Just M1]]
matrify (More s h) = zipWith (++) ((liftPhase (sigPhaseIndex s 0 0) (matrify h)) ++ (liftPhase (sigPhaseIndex s 1 0) (matrify h)))
								  ((liftPhase (sigPhaseIndex s 0 1) (matrify h)) ++ (liftPhase (sigPhaseIndex s 1 1) (matrify h)))
  where sigPhaseIndex s n m = matrify (More s (Only P1)) !! n !! m 
        -- give a pauli and two numbers and get the phase at row n, col m

liftPhase :: Maybe Phase -> ([[Maybe Phase]] -> [[Maybe Phase]])
liftPhase p = (map (map ((multP <$> p) <*>))) -- lifts a phase p to a function [[Maybe Phase]] -> [[Maybe Phase]]

-- Matrix-Show, use in combination with matrify to print sigma matrices to terminal
mShow :: [[Maybe Phase]] -> IO ()
mShow k = putStr $ (foldl1 (++)) $ map (\x -> (printRow (getRow k x))) [0..(length k - 1)]
  where printRow m = "[" ++ foldr (++) "] \n" (map ((\x -> " " ++ x ++ " ") . thePrinter) m)
        thePrinter Nothing = " 0"
        thePrinter (Just p) = show p
        getRow w n = (map  (!! n)) w