|
| 1 | +from scipy.linalg import expm |
| 2 | +import numpy as np |
| 3 | + |
| 4 | +#solve expm error |
| 5 | +A=[[1, 0], [0, 1]]; A = np.array(A); expm(-1j * A) |
| 6 | + |
| 7 | +import inspect |
| 8 | +import tensorcircuit as tc |
| 9 | +import random |
| 10 | +import math |
| 11 | +import matplotlib.pyplot as plt |
| 12 | +import time |
| 13 | + |
| 14 | +tc.set_backend("tensorflow") |
| 15 | + |
| 16 | +#calculate the matirx of kth qubit exert matrix[[a, b], [c, d]] |
| 17 | +def up_to_matrixx(k, a, b, c, d): |
| 18 | + I2 = np.array([[1,0],[0,1]])*(1+0j); K=np.array([[a,b],[c,d]])*(1+0j); um=I2; |
| 19 | + if k == 0: |
| 20 | + um = K; |
| 21 | + for i in range(1, N): |
| 22 | + if i == k: |
| 23 | + um = np.kron(um, K) |
| 24 | + else: |
| 25 | + um = np.kron(um, I2) |
| 26 | + return um |
| 27 | + |
| 28 | +#realize R gates in paper |
| 29 | +def R_gate(k): |
| 30 | + if door[k][0] == 0: |
| 31 | + c.rx(door[k][1]+1,theta=ODE_theta[k]) |
| 32 | + if door[k][0] == 1: |
| 33 | + c.ry(door[k][1]+1,theta=ODE_theta[k]) |
| 34 | + if door[k][0] == 2: |
| 35 | + c.rz(door[k][1]+1,theta=ODE_theta[k]) |
| 36 | + if door[k][0] == 3: |
| 37 | + c.rxx(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 38 | + if door[k][0] == 4: |
| 39 | + c.ryy(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 40 | + if door[k][0] == 5: |
| 41 | + c.rzz(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 42 | + if door[k][0] == 6: |
| 43 | + c.crx(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 44 | + if door[k][0] == 7: |
| 45 | + c.cry(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 46 | + if door[k][0] == 8: |
| 47 | + c.crz(door[k][1]+1,door[k][2]+1,theta=ODE_theta[k]) |
| 48 | + |
| 49 | +#realize U gates in paper |
| 50 | +def U_gate(k): |
| 51 | + if door[k][0] == 0: |
| 52 | + c.cx(0,door[k][1]+1) |
| 53 | + if door[k][0] == 1: |
| 54 | + c.cy(0,door[k][1]+1) |
| 55 | + if door[k][0] == 2: |
| 56 | + c.cz(0,door[k][1]+1) |
| 57 | + if door[k][0] == 3: |
| 58 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0],unitary=tc.gates._xx_matrix) |
| 59 | + if door[k][0] == 4: |
| 60 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0],unitary=tc.gates._yy_matrix) |
| 61 | + if door[k][0] == 5: |
| 62 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0],unitary=tc.gates._zz_matrix) |
| 63 | + if door[k][0] == 6: |
| 64 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._xx_matrix) |
| 65 | + if door[k][0] == 7: |
| 66 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._yy_matrix) |
| 67 | + if door[k][0] == 8: |
| 68 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._zz_matrix) |
| 69 | + |
| 70 | +#realize Hamilton gates in ancillary circuit |
| 71 | +def H_gate(q): |
| 72 | + if h_door[q][0] == 0: |
| 73 | + c.cx(0,h_door[q][1]+1) |
| 74 | + if h_door[q][0] == 1: |
| 75 | + c.cy(0,h_door[q][1]+1) |
| 76 | + if h_door[q][0] == 2: |
| 77 | + c.cz(0,h_door[q][1]+1) |
| 78 | + if h_door[q][0] == 3: |
| 79 | + c.multicontrol(0,h_door[q][1]+1,h_door[q][2]+1,ctrl=[0],unitary=tc.gates._xx_matrix) |
| 80 | + if h_door[q][0] == 4: |
| 81 | + c.multicontrol(0,h_door[q][1]+1,h_door[q][2]+1,ctrl=[0],unitary=tc.gates._yy_matrix) |
| 82 | + if h_door[q][0] == 5: |
| 83 | + c.multicontrol(0,h_door[q][1]+1,h_door[q][2]+1,ctrl=[0],unitary=tc.gates._zz_matrix) |
| 84 | + if h_door[q][0] == 6: |
| 85 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._xx_matrix) |
| 86 | + if h_door[q][0] == 7: |
| 87 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._yy_matrix) |
| 88 | + if h_door[q][0] == 8: |
| 89 | + c.multicontrol(0,door[k][1]+1,door[k][2]+1,ctrl=[0,door[k][1]+1],unitary=tc.gates._zz_matrix) |
| 90 | + |
| 91 | +#use quantum circuit to calculate coefficient of variation A and C in paper |
| 92 | +def find_ACkq(mod, theta_x, k, q, whi): |
| 93 | + #mod: a in paper; theta_x: theta in paper; k, q: A[k, q] or C[k] qth term; whi: whi=0 A whi=1 C |
| 94 | + global c |
| 95 | + ancilla = np.array([1, np.exp(1j * theta_x)]) / np.sqrt(2) |
| 96 | + c = tc.Circuit(N+1,inputs = np.kron(ancilla, state)) |
| 97 | + for i in range(len(door)): |
| 98 | + if i == k: |
| 99 | + c.x(0) |
| 100 | + U_gate(i) |
| 101 | + c.x(0) |
| 102 | + if whi == 0 and i == q: |
| 103 | + U_gate(i) |
| 104 | + R_gate(i) |
| 105 | + break |
| 106 | + R_gate(i) |
| 107 | + if whi == 1: |
| 108 | + H_gate(q) |
| 109 | + pstar = np.real(np.array(c.expectation([np.array([[1, 1], [1, 1]]) / 2, [0]]))) |
| 110 | + return mod * (2 * pstar - 1) |
| 111 | + |
| 112 | +#use original quantum circuit simulate with c |
| 113 | +def simulation(): |
| 114 | + global c |
| 115 | + c=tc.Circuit(N,inputs=state) |
| 116 | + for k in range(len(door)): |
| 117 | + if door[k][0]==0: |
| 118 | + c.rx(door[k][1],theta=ODE_theta[k]) |
| 119 | + if door[k][0]==1: |
| 120 | + c.ry(door[k][1],theta=ODE_theta[k]) |
| 121 | + if door[k][0]==2: |
| 122 | + c.rz(door[k][1],theta=ODE_theta[k]) |
| 123 | + if door[k][0]==3: |
| 124 | + c.rxx(door[k][1],door[k][2],theta=ODE_theta[k]) |
| 125 | + if door[k][0]==4: |
| 126 | + c.ryy(door[k][1],door[k][2],theta=ODE_theta[k]) |
| 127 | + if door[k][0]==5: |
| 128 | + c.rzz(door[k][1],door[k][2],theta=ODE_theta[k]) |
| 129 | + |
| 130 | +if __name__ == '__main__': |
| 131 | + |
| 132 | + #l: layers; h and J: coefficient of Hamilton; L_var and L_num: results of variation method and numerical method |
| 133 | + N=3; l=2; J=1/4; dt=0.05; t=1; h=[]; L_var=[]; L_num=[]; x_value=[]; |
| 134 | + |
| 135 | + how_variation = 0 #0 McLachlan 1 time-dependent |
| 136 | + |
| 137 | + #the priciple correspond with all gates |
| 138 | + #the first term: 0rx,1ry,2rz,3rxx,4ryy,5rzz,6crx,7cry,8crz; |
| 139 | + #the second and the third term: num/ctrl+num |
| 140 | + #f: coefficient with simulation gates in paper |
| 141 | + door = []; h_door = []; f = [] |
| 142 | + for k in range(l): |
| 143 | + for i in range(N): |
| 144 | + f.append(-0.5j) |
| 145 | + door.append([0, i]) |
| 146 | + for i in range(N - 1): |
| 147 | + f.append(-1j) |
| 148 | + door.append([5, i, i + 1]) |
| 149 | + for i in range(N - 1): |
| 150 | + f.append(-1j) |
| 151 | + door.append([3, i, i + 1]) |
| 152 | + for i in range(N): |
| 153 | + h.append(1) |
| 154 | + h_door.append([0, i]) |
| 155 | + for i in range(N-1): |
| 156 | + h.append(J); h_door.append([5, i, i + 1]) |
| 157 | + |
| 158 | + #initial state |
| 159 | + state = np.zeros(1 << N); state[0]=1 |
| 160 | + |
| 161 | + #numerical realize H |
| 162 | + H = np.zeros((1<<N, 1<<N)) * 1j |
| 163 | + for i in range(N-1): |
| 164 | + H += J*up_to_matrixx(i, 1, 0, 0, -1) @ up_to_matrixx(i + 1, 1, 0, 0, -1) |
| 165 | + for i in range(N): |
| 166 | + H += h[i] * up_to_matrixx(i, 0, 1, 1, 0) |
| 167 | + |
| 168 | + #variation realize |
| 169 | + ODE_theta = np.zeros(len(door)) |
| 170 | + for T in range(int(t / dt)): |
| 171 | + #calculate coefficient in paper |
| 172 | + A = np.zeros((len(door), len(door))); C = np.zeros(len(door)) |
| 173 | + for k in range(len(door)): |
| 174 | + for q in range(len(door)): |
| 175 | + if k > q: |
| 176 | + A[k, q] = A[q, k] |
| 177 | + continue |
| 178 | + if how_variation == 0: |
| 179 | + A[k, q] = find_ACkq(abs(f[k] * f[q]), np.angle(f[q]) - np.angle(f[k]), k, q, 0) |
| 180 | + if how_variation == 1: |
| 181 | + A[k, q] = find_ACkq(abs(f[k] * f[q]), np.angle(f[q]) - np.angle(f[k]) - math.pi / 2, k, q, 0) |
| 182 | + for k in range(len(door)): |
| 183 | + for q in range(len(h)): |
| 184 | + if how_variation == 0: |
| 185 | + C[k] += find_ACkq(abs(f[k] * h[q]), np.angle(h[q]) - np.angle(f[k]) - math.pi / 2, k, q, 1) |
| 186 | + if how_variation == 1: |
| 187 | + C[k] += find_ACkq(-abs(f[k] * h[q]), np.angle(h[q]) - np.angle(f[k]), k, q, 1) |
| 188 | + |
| 189 | + #calculate parameter and its derivative |
| 190 | + A += np.eye(len(door)) * 1e-5 |
| 191 | + ODE_dtheta = np.linalg.solve(A, C) |
| 192 | + print(ODE_dtheta) |
| 193 | + for i in range(len(door)): |
| 194 | + ODE_theta[i] += ODE_dtheta[i] * dt |
| 195 | + |
| 196 | + #numerical results |
| 197 | + simulation() |
| 198 | + ep = expm(-1j * H * (T + 1) * dt) @ state |
| 199 | + L_num.append(np.real(np.array(ep.conj().T @ up_to_matrixx(1, 0, 1, 1, 0) @ ep)).tolist()) |
| 200 | + |
| 201 | + #variation results |
| 202 | + L_var.append(np.real(np.array(c.expectation([tc.gates.x(), [1]]))).tolist()) |
| 203 | + |
| 204 | + x_value.append((T + 1) * dt) |
| 205 | + print([(T + 1) * dt, L_num[T] - L_var[T]]) |
| 206 | + plt.plot(x_value, L_var, color = 'green') |
| 207 | + plt.plot(x_value, L_num, color = 'red') |
| 208 | + plt.show() |
0 commit comments