-
Notifications
You must be signed in to change notification settings - Fork 2
/
utilities.py
219 lines (169 loc) · 4.51 KB
/
utilities.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
import numpy as np
import scipy as scp
import time
import cmath
from lazyarray import larray
"""
A file for storing various mathematical helper functions for each implementation.
"""
###################################################################
#tensor products
def tensor(b,a):
"""
Function that takes the tensor product of 2 matrices or vectors.
Used by Dense implementation (dense.py)
Behaves identically to np.kron()
"""
#Dimension of output
a0 = a.shape[0]
a1 = a.shape[1]
b0 = b.shape[0]
b1 = b.shape[1]
outdim = (a0*b0,a1*b1)
#Initialise output matrix with zeros
output = np.zeros(outdim,dtype=complex)
#Calculate output matrix
for x in range(outdim[0]):
for y in range(outdim[1]):
output[x,y] = a[x%a0,y%a1]*b[x//a0,y//a1]
return output
def tensor_sparse_qubit(A,B):
'''
Function that takes the tensor product of 2 Qubits.
Used by sparse implementation (sparse.py) but does not take
advantage of sparse form, since speed difference is minimal
'''
#convert to dense matrices
a = A.toarray()[0]
b = B.toarray()[0]
#Dimension of output
a0 = a.shape[0]
b0 = b.shape[0]
outdim = (1,a0*b0)
#Initialise output matrix with zeros
output = np.zeros(outdim,dtype=complex)
#Calculate output matrix
for x in range(outdim[0]):
output[0,x] = a[x%a0]*b[x//a0]
return(output)
def tensor_sparse_gate(a,b):
'''
Function that takes the tensor product of 2 Gates.
Used by sparse implementation (sparse.py)
Gives the same result as sp.kron
'''
#import sparse inside as it breaks lazy if outside
from scipy import sparse as sp
#output dimensions
a0 = a.shape[0]
a1 = a.shape[1]
b0 = b.shape[0]
b1 = b.shape[1]
outdim = (a0*b0,a1*b1)
#catch for multiplying empty matrices
if a.nnz == 0 or b.nnz == 0:
return sp.coo_matrix(outdim)
#convert a from bsr to more standard csc version of sparse
a = sp.csc_matrix(a)
#convert b to dense
b = b.toarray()
#reshape a into data and calculate
data = a.data.repeat(b.size).reshape(-1,b0,b1)
data = data * b
return sp.bsr_matrix((data,a.indices,a.indptr), shape=outdim)
def tensor_lazy(b,a):
"""
Function that takes the tensor product of 2 matrices or vectors.
Used by lazy implementation (lazy.py)
"""
#Dimension of output
a0 = a.shape[0]
a1 = a.shape[1]
b0 = b.shape[0]
b1 = b.shape[1]
outdim = (a0*b0,a1*b1)
#function to feed creation of lazy output matrix
#i and j refer to matrix index
def kron(i,j):
return a[i%a0,j%a1]*b[i//a0,j//a1]
#create output
output = larray(kron,shape=outdim)
return output
#########################################################################
#matrix multiplication
def lazy_mul_gate(b,a):
"""
Function for the matrix multiplication of 2 gates.
Used by lazy implementation (lazy.py)
"""
#Dimension of output
a0 = a.shape[0]
a1 = a.shape[1]
b0 = b.shape[0]
b1 = b.shape[1]
outdim = (b0,a1)
#function to feed creation of lazy output matrix
#i and j refer to matrix index
def mul(i,j):
elem = sum(map(lambda n: b[i][n]*a[n][j],range(b1)))
return elem
output = larray(mul,shape=outdim)
return(output)
def lazy_mul(b,a):
"""
Function for the matrix multiplication of 2 qubits.
Used by lazy implementation (lazy.py)
"""
#Dimension of output
a0 = a.shape[0]
a1 = a.shape[1]
b0 = b.shape[0]
b1 = b.shape[1]
outdim = (b0,a1)
#function which feeds creation of output matrix
#i and j refer to matrix index
def mul(i,j):
elem = sum(map(lambda n: b[i][n]*a[n],range(b1)))
return elem
#create output
output = larray(mul,shape = outdim)
return output
#######################################################################
#qubit swapping
def perm_matrix(n,index1,index2):
'''
generates a permutation matrix from a list of pairs of numbers
to swap
Used by all implementations
'''
#catch user errors
assert index1!=index2, "Cant swap qubit with itself"
assert (index1<n) and (index2<n), "Cant swap qubits beyond size of gate"
#shuffle index
index1 = n-index1-1
index2 = n-index2-1
b = 2**index1 + 2**index2
#work out which parts correspond to the index being swapped
swaps = []
for x in range(2**n):
for y in range(x):
if((x^y==b) and (count_bits(x)==count_bits(y))):
swaps.append((x,y))
#initialise output matrix with an identity
size = 2**n
i = np.identity(size)
#alter identity into swap matrix
for pairs in swaps:
temp = i[pairs[0]].copy()
i[pairs[0]] = i[pairs[1]]
i[pairs[1]] = temp
return i
def count_bits(n):
'''
helper function for perm_matrix
counts bits recursively
'''
if n==0:
return 0
else:
return (n&1)+count_bits(n>>1)