-
Notifications
You must be signed in to change notification settings - Fork 0
/
MainSimulationNqbits.py
114 lines (89 loc) · 2.53 KB
/
MainSimulationNqbits.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
import math
import numpy as np
from QuantumPackage import getSci
from QuantumPackage import getInitialState
from QuantumPackage import getHamiltonian
import matplotlib.pyplot as plt
import scipy.linalg as scl
import time
start_time = time.time()
#Variables
Nspins = 5
j = 13
Bz = 5
hbar = 1
# Matrices
Sx = np.matrix([[0,1],[1,0]])
Sz = np.matrix([[1,0],[0,-1]])
I = np.matrix([[1,0],[0,1]])
# Hamiltonianio
H = getHamiltonian(j,Bz,Nspins,Sx,Sz)
# Condicion inicial
X = np.matrix([[1],[-1]])/math.sqrt(2)
Psi_0 = getInitialState(X,Nspins)
# Tiempo
ti = 0
tf = 15/Bz
N = 1000
dt =(tf-ti)/(N)
t = np.arange(ti,tf,dt)
# Funcion de onda inicial
psi = Psi_0
# Magnetizacion promedio en x y z
Mx = 0
Mz = 0
# Inicialización operadores de magnetizacion
for i in range(1,Nspins+1):
Sxi = getSci(Sx,i,Nspins)
Szi = getSci(Sz,i,Nspins)
Mx = Mx + Sxi/Nspins
Mz = Mz + Szi/Nspins
# Coeficientes function de onda
c1 = np.zeros([1,N], dtype=np.complex_)
c2 = np.zeros([1,N], dtype=np.complex_)
c3 = np.zeros([1,N], dtype=np.complex_)
c4 = np.zeros([1,N], dtype=np.complex_)
# Magnetizacion promedio
PromMx = np.zeros([1,N], dtype=np.complex_)
PromMZ = np.zeros([1,N], dtype=np.complex_)
# Operador evolucion temporal
U = scl.expm(-1j*H/hbar*dt)
# Simulacion de la ecuacion de Schrodinger
for n in range(0,N):
if n == 0:
psi = Psi_0
PromMx[0][n] = psi.getH() * Mx * psi
PromMZ[0][n] = psi.getH() * Mz * psi
else:
psi = U*psi
PromMx[0][n] = psi.getH() * Mx * psi
PromMZ[0][n] = psi.getH() * Mz * psi
# Funciones temporales para la magnetizacion
MeanMz = PromMZ[0][:]
MeanMx = PromMx[0][:]
lineW = 2 # Line thickness
plt.figure(figsize=(10,6), tight_layout=True)
plt.plot(t*Bz, MeanMz,'-r', linewidth=3);
plt.ylabel(r'$\langle M_z \rangle$');plt.xlabel(r'$B_z t$')
axes = plt.gca()
axes.xaxis.label.set_size(22)
axes.yaxis.label.set_size(22)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim(0, tf*Bz)
plt.savefig('Mz.png')
lineW = 2 # Line thickness
plt.figure(figsize=(10,6), tight_layout=True)
plt.plot(t*Bz, MeanMx,'-b', linewidth=3);
plt.ylabel(r'$\langle M_x \rangle$');plt.xlabel(r'$B_z t$')
axes = plt.gca()
axes.xaxis.label.set_size(22)
axes.yaxis.label.set_size(22)
plt.xticks(fontsize=20)
plt.yticks(fontsize=20)
plt.xlim(0, tf*Bz)
plt.savefig('Mx.png')
print('Bz : ', Bz)
print('j : ', j)
print('t : ', t)
print("--- %s seconds ---" % (time.time() - start_time))