{-# OPTIONS_GHC -Wall #-}
{-# LANGUAGE Trustworthy #-}
{-# LANGUAGE CPP #-}
module Physics.Learn.BeamStack
(
BeamStack()
, randomBeam
, split
, recombine
, applyBField
, dropBeam
, flipBeams
, numBeams
, detect
, splitX
, splitY
, splitZ
, applyBFieldX
, applyBFieldY
, applyBFieldZ
, recombineX
, recombineY
, recombineZ
, xpFilter
, xmFilter
, zpFilter
, zmFilter
)
where
import Physics.Learn.QuantumMat
( zp
, zm
, nm
, np
, couter
, oneQubitMixed
)
import Numeric.LinearAlgebra
( C
, Vector
, Matrix
, iC
, (<>)
, kronecker
, fromLists
, toList
, toLists
, scale
, size
, takeDiag
, ident
, tr
)
import Data.Complex
( Complex(..)
, realPart
, imagPart
)
import Data.List
( intercalate
)
#if MIN_VERSION_base(4,11,0)
import Prelude hiding ((<>))
#endif
data BeamStack = BeamStack (Matrix C)
showOneBeam :: Double -> String
showOneBeam :: Double -> String
showOneBeam Double
r = String
"Beam of intensity " String -> String -> String
forall a. [a] -> [a] -> [a]
++ Double -> String
forall a. Show a => a -> String
show Double
r
instance Show BeamStack where
show :: BeamStack -> String
show BeamStack
b = String -> [String] -> String
forall a. [a] -> [[a]] -> [a]
intercalate String
"\n" ([String] -> String) -> [String] -> String
forall a b. (a -> b) -> a -> b
$ (Double -> String) -> [Double] -> [String]
forall a b. (a -> b) -> [a] -> [b]
map Double -> String
showOneBeam (BeamStack -> [Double]
detect BeamStack
b)
randomBeam :: BeamStack
randomBeam :: BeamStack
randomBeam = Matrix C -> BeamStack
BeamStack Matrix C
oneQubitMixed
extendWithZeros :: Matrix C -> Matrix C
extendWithZeros :: Matrix C -> Matrix C
extendWithZeros Matrix C
m
= let (Int
_,Int
q) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
ml :: [[C]]
ml = Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists Matrix C
m
in [[C]] -> Matrix C
forall t. Element t => [[t]] -> Matrix t
fromLists ([[C]] -> Matrix C) -> [[C]] -> Matrix C
forall a b. (a -> b) -> a -> b
$ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map ([C] -> [C] -> [C]
forall a. [a] -> [a] -> [a]
++ [C
0,C
0]) [[C]]
ml
[[C]] -> [[C]] -> [[C]]
forall a. [a] -> [a] -> [a]
++ [Int -> C -> [C]
forall a. Int -> a -> [a]
replicate (Int
qInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
2) C
0, Int -> C -> [C]
forall a. Int -> a -> [a]
replicate (Int
qInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
2) C
0]
reduceMat :: Matrix C -> Matrix C
reduceMat :: Matrix C -> Matrix C
reduceMat Matrix C
m
= let (Int
p,Int
q) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
ml :: [[C]]
ml = Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists Matrix C
m
in [[C]] -> Matrix C
forall t. Element t => [[t]] -> Matrix t
fromLists ([[C]] -> Matrix C) -> [[C]] -> Matrix C
forall a b. (a -> b) -> a -> b
$ Int -> [[C]] -> [[C]]
forall a. Int -> [a] -> [a]
take (Int
pInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
2) ([[C]] -> [[C]]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> a -> b
$ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> [C] -> [C]
forall a. Int -> [a] -> [a]
take (Int
qInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
2)) [[C]]
ml
checkedRealPart :: C -> Double
checkedRealPart :: C -> Double
checkedRealPart C
c
= let eps :: Double
eps = Double
1e-14
in if C -> Double
forall a. Complex a -> a
imagPart C
c Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
< Double
eps
then C -> Double
forall a. Complex a -> a
realPart C
c
else String -> Double
forall a. HasCallStack => String -> a
error (String -> Double) -> String -> Double
forall a b. (a -> b) -> a -> b
$ String
"checkRealPart: imagPart = " String -> String -> String
forall a. [a] -> [a] -> [a]
++ Double -> String
forall a. Show a => a -> String
show (C -> Double
forall a. Complex a -> a
imagPart C
c)
detect :: BeamStack -> [Double]
detect :: BeamStack -> [Double]
detect (BeamStack Matrix C
m)
= [C] -> [Double]
addAlternate ([C] -> [Double]) -> [C] -> [Double]
forall a b. (a -> b) -> a -> b
$ Vector C -> [C]
forall a. Storable a => Vector a -> [a]
toList (Vector C -> [C]) -> Vector C -> [C]
forall a b. (a -> b) -> a -> b
$ Matrix C -> Vector C
forall t. Element t => Matrix t -> Vector t
takeDiag Matrix C
m
addAlternate :: [C] -> [Double]
addAlternate :: [C] -> [Double]
addAlternate [] = []
addAlternate [C
_] = String -> [Double]
forall a. HasCallStack => String -> a
error String
"addAlternate needs even number of elements"
addAlternate (C
x:C
y:[C]
xs) = C -> Double
checkedRealPart (C
xC -> C -> C
forall a. Num a => a -> a -> a
+C
y) Double -> [Double] -> [Double]
forall a. a -> [a] -> [a]
: [C] -> [Double]
addAlternate [C]
xs
dropBeam :: BeamStack -> BeamStack
dropBeam :: BeamStack -> BeamStack
dropBeam (BeamStack Matrix C
m) = Matrix C -> BeamStack
BeamStack (Matrix C -> Matrix C
reduceMat Matrix C
m)
numBeams :: BeamStack -> Int
numBeams :: BeamStack -> Int
numBeams (BeamStack Matrix C
m)
= let (Int
p,Int
_) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
in Int
p Int -> Int -> Int
forall a. Integral a => a -> a -> a
`div` Int
2
flipBeams :: BeamStack -> BeamStack
flipBeams :: BeamStack -> BeamStack
flipBeams (BeamStack Matrix C
m)
= let (Int
d,Int
_) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
fl :: Matrix C
fl = Int -> Matrix C
flipMat Int
d
in Matrix C -> BeamStack
BeamStack (Matrix C -> BeamStack) -> Matrix C -> BeamStack
forall a b. (a -> b) -> a -> b
$ Matrix C
fl Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C
m Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C -> Matrix C
forall m mt. Transposable m mt => m -> mt
tr Matrix C
fl
flipMat :: Int -> Matrix C
flipMat :: Int -> Matrix C
flipMat Int
d = Int -> Matrix C -> Matrix C
bigM Int
d ([[C]] -> Matrix C
forall t. Element t => [[t]] -> Matrix t
fromLists [[C
0,C
0,C
1,C
0]
,[C
0,C
0,C
0,C
1]
,[C
1,C
0,C
0,C
0]
,[C
0,C
1,C
0,C
0]])
bigM2 :: Int -> Matrix C -> Matrix C
bigM2 :: Int -> Matrix C -> Matrix C
bigM2 Int
d Matrix C
m
| Int
d Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
2 = String -> Matrix C
forall a. HasCallStack => String -> a
error String
"bigM2 requires d >= 2"
| Int -> Bool
forall a. Integral a => a -> Bool
odd Int
d = String -> Matrix C
forall a. HasCallStack => String -> a
error String
"bigM2 requires even d"
| Bool
otherwise = [[C]] -> Matrix C
forall t. Element t => [[t]] -> Matrix t
fromLists ([[C]] -> Matrix C) -> [[C]] -> Matrix C
forall a b. (a -> b) -> a -> b
$ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map ([C] -> [C] -> [C]
forall a. [a] -> [a] -> [a]
++ [C
0,C
0]) (Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident (Int
dInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
2)))
[[C]] -> [[C]] -> [[C]]
forall a. [a] -> [a] -> [a]
++ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> C -> [C]
forall a. Int -> a -> [a]
replicate (Int
dInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
2) C
0 [C] -> [C] -> [C]
forall a. [a] -> [a] -> [a]
++) (Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists Matrix C
m)
bigM :: Int -> Matrix C -> Matrix C
bigM :: Int -> Matrix C -> Matrix C
bigM Int
d Matrix C
m
| Int
d Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
4 = String -> Matrix C
forall a. HasCallStack => String -> a
error String
"bigM requires d >= 4"
| Int -> Bool
forall a. Integral a => a -> Bool
odd Int
d = String -> Matrix C
forall a. HasCallStack => String -> a
error String
"bigM requires even d"
| Bool
otherwise = [[C]] -> Matrix C
forall t. Element t => [[t]] -> Matrix t
fromLists ([[C]] -> Matrix C) -> [[C]] -> Matrix C
forall a b. (a -> b) -> a -> b
$ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map ([C] -> [C] -> [C]
forall a. [a] -> [a] -> [a]
++ [C
0,C
0,C
0,C
0]) (Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists (Int -> Matrix C
forall a. (Num a, Element a) => Int -> Matrix a
ident (Int
dInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
4)))
[[C]] -> [[C]] -> [[C]]
forall a. [a] -> [a] -> [a]
++ ([C] -> [C]) -> [[C]] -> [[C]]
forall a b. (a -> b) -> [a] -> [b]
map (Int -> C -> [C]
forall a. Int -> a -> [a]
replicate (Int
dInt -> Int -> Int
forall a. Num a => a -> a -> a
-Int
4) C
0 [C] -> [C] -> [C]
forall a. [a] -> [a] -> [a]
++) (Matrix C -> [[C]]
forall t. Element t => Matrix t -> [[t]]
toLists Matrix C
m)
s :: Double -> Double -> Matrix C
s :: Double -> Double -> Matrix C
s Double
theta Double
phi = Matrix C -> Matrix C -> Matrix C
forall t. Product t => Matrix t -> Matrix t -> Matrix t
kronecker (Vector C
u Vector C -> Vector C -> Matrix C
`couter` Vector C
u) (Double -> Double -> Vector C
np Double
theta Double
phi Vector C -> Vector C -> Matrix C
`couter` Double -> Double -> Vector C
np Double
theta Double
phi)
Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C -> Matrix C -> Matrix C
forall t. Product t => Matrix t -> Matrix t -> Matrix t
kronecker (Vector C
l Vector C -> Vector C -> Matrix C
`couter` Vector C
u) (Double -> Double -> Vector C
nm Double
theta Double
phi Vector C -> Vector C -> Matrix C
`couter` Double -> Double -> Vector C
nm Double
theta Double
phi)
Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C -> Matrix C -> Matrix C
forall t. Product t => Matrix t -> Matrix t -> Matrix t
kronecker (Vector C
u Vector C -> Vector C -> Matrix C
`couter` Vector C
l) (Double -> Double -> Vector C
nm Double
theta Double
phi Vector C -> Vector C -> Matrix C
`couter` Double -> Double -> Vector C
nm Double
theta Double
phi)
Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ Matrix C -> Matrix C -> Matrix C
forall t. Product t => Matrix t -> Matrix t -> Matrix t
kronecker (Vector C
l Vector C -> Vector C -> Matrix C
`couter` Vector C
l) (Double -> Double -> Vector C
np Double
theta Double
phi Vector C -> Vector C -> Matrix C
`couter` Double -> Double -> Vector C
np Double
theta Double
phi)
u :: Vector C
u :: Vector C
u = Vector C
zp
l :: Vector C
l :: Vector C
l = Vector C
zm
split :: Double -> Double -> BeamStack -> BeamStack
split :: Double -> Double -> BeamStack -> BeamStack
split Double
theta Double
phi (BeamStack Matrix C
m)
= let m' :: Matrix C
m' = Matrix C -> Matrix C
extendWithZeros Matrix C
m
(Int
p,Int
_) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m'
ss :: Matrix C
ss = Int -> Matrix C -> Matrix C
bigM Int
p (Double -> Double -> Matrix C
s Double
theta Double
phi)
in Matrix C -> BeamStack
BeamStack (Matrix C -> BeamStack) -> Matrix C -> BeamStack
forall a b. (a -> b) -> a -> b
$ Matrix C
ss Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C
m' Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C -> Matrix C
forall m mt. Transposable m mt => m -> mt
tr Matrix C
ss
recombine :: Double -> Double -> BeamStack -> BeamStack
recombine :: Double -> Double -> BeamStack -> BeamStack
recombine Double
theta Double
phi (BeamStack Matrix C
m)
= let (Int
d,Int
_) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
ss :: Matrix C
ss = Int -> Matrix C -> Matrix C
bigM Int
d (Double -> Double -> Matrix C
s Double
theta Double
phi)
in BeamStack -> BeamStack
dropBeam (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall a b. (a -> b) -> a -> b
$ Matrix C -> BeamStack
BeamStack (Matrix C -> BeamStack) -> Matrix C -> BeamStack
forall a b. (a -> b) -> a -> b
$ Matrix C
ss Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C
m Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C -> Matrix C
forall m mt. Transposable m mt => m -> mt
tr Matrix C
ss
mag2x2 :: Double -> Double -> Double -> Matrix C
mag2x2 :: Double -> Double -> Double -> Matrix C
mag2x2 Double
theta Double
phi Double
omegaT
= let z :: C
z = C
iC C -> C -> C
forall a. Num a => a -> a -> a
* (Double
omegaT Double -> Double -> C
forall a. a -> a -> Complex a
:+ Double
0) C -> C -> C
forall a. Fractional a => a -> a -> a
/ C
2
np' :: Vector C
np' = Double -> Double -> Vector C
np Double
theta Double
phi
nm' :: Vector C
nm' = Double -> Double -> Vector C
nm Double
theta Double
phi
in 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
exp C
z ) (Vector C
np' Vector C -> Vector C -> Matrix C
`couter` Vector C
np')
Matrix C -> Matrix C -> Matrix C
forall a. Num a => a -> a -> a
+ 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
exp (-C
z)) (Vector C
nm' Vector C -> Vector C -> Matrix C
`couter` Vector C
nm')
applyBField :: Double -> Double -> Double -> BeamStack -> BeamStack
applyBField :: Double -> Double -> Double -> BeamStack -> BeamStack
applyBField Double
theta Double
phi Double
omegaT (BeamStack Matrix C
m)
= let (Int
d,Int
_) = Matrix C -> IndexOf Matrix
forall (c :: * -> *) t. Container c t => c t -> IndexOf c
size Matrix C
m
uu :: Matrix C
uu = Int -> Matrix C -> Matrix C
bigM2 Int
d (Double -> Double -> Double -> Matrix C
mag2x2 Double
theta Double
phi Double
omegaT)
in Matrix C -> BeamStack
BeamStack (Matrix C -> BeamStack) -> Matrix C -> BeamStack
forall a b. (a -> b) -> a -> b
$ Matrix C
uu Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C
m Matrix C -> Matrix C -> Matrix C
forall t. Numeric t => Matrix t -> Matrix t -> Matrix t
<> Matrix C -> Matrix C
forall m mt. Transposable m mt => m -> mt
tr Matrix C
uu
splitX :: BeamStack -> BeamStack
splitX :: BeamStack -> BeamStack
splitX = Double -> Double -> BeamStack -> BeamStack
split (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) Double
0
splitY :: BeamStack -> BeamStack
splitY :: BeamStack -> BeamStack
splitY = Double -> Double -> BeamStack -> BeamStack
split (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
splitZ :: BeamStack -> BeamStack
splitZ :: BeamStack -> BeamStack
splitZ = Double -> Double -> BeamStack -> BeamStack
split Double
0 Double
0
applyBFieldX :: Double -> BeamStack -> BeamStack
applyBFieldX :: Double -> BeamStack -> BeamStack
applyBFieldX = Double -> Double -> Double -> BeamStack -> BeamStack
applyBField (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) Double
0
applyBFieldY :: Double -> BeamStack -> BeamStack
applyBFieldY :: Double -> BeamStack -> BeamStack
applyBFieldY = Double -> Double -> Double -> BeamStack -> BeamStack
applyBField (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
applyBFieldZ :: Double -> BeamStack -> BeamStack
applyBFieldZ :: Double -> BeamStack -> BeamStack
applyBFieldZ = Double -> Double -> Double -> BeamStack -> BeamStack
applyBField Double
0 Double
0
recombineX :: BeamStack -> BeamStack
recombineX :: BeamStack -> BeamStack
recombineX = Double -> Double -> BeamStack -> BeamStack
recombine (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) Double
0
recombineY :: BeamStack -> BeamStack
recombineY :: BeamStack -> BeamStack
recombineY = Double -> Double -> BeamStack -> BeamStack
recombine (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2) (Double
forall a. Floating a => a
piDouble -> Double -> Double
forall a. Fractional a => a -> a -> a
/Double
2)
recombineZ :: BeamStack -> BeamStack
recombineZ :: BeamStack -> BeamStack
recombineZ = Double -> Double -> BeamStack -> BeamStack
recombine Double
0 Double
0
xpFilter :: BeamStack -> BeamStack
xpFilter :: BeamStack -> BeamStack
xpFilter = BeamStack -> BeamStack
dropBeam (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
splitX
xmFilter :: BeamStack -> BeamStack
xmFilter :: BeamStack -> BeamStack
xmFilter = BeamStack -> BeamStack
dropBeam (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
flipBeams (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
splitX
zpFilter :: BeamStack -> BeamStack
zpFilter :: BeamStack -> BeamStack
zpFilter = BeamStack -> BeamStack
dropBeam (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
splitZ
zmFilter :: BeamStack -> BeamStack
zmFilter :: BeamStack -> BeamStack
zmFilter = BeamStack -> BeamStack
dropBeam (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
flipBeams (BeamStack -> BeamStack)
-> (BeamStack -> BeamStack) -> BeamStack -> BeamStack
forall b c a. (b -> c) -> (a -> b) -> a -> c
. BeamStack -> BeamStack
splitZ