import numpy as np
import matplotlib.pyplot as plt
import math
import cmath
from scipy import linalg
d = 1.0e-6
a = 0.82*d
epsia = 1.0
epsib = 1.0
e_p0p0 = epsib+(a/d)*(a/d)*(epsia-epsib)
e_p1p0 = (a/d)*(a/d)*(epsia-epsib)*(math.sin(+math.pi*(a/d))/(+math.pi*(a/d)))
e_m1p0 = (a/d)*(a/d)*(epsia-epsib)*(math.sin(-math.pi*(a/d))/(-math.pi*(a/d)))
Nq = 100
qi = 0.0
qf = 0.5
dq = (qf-qi)/Nq
qx = np.zeros(Nq+1)
qy = np.zeros(Nq+1)
Om1_M = np.zeros(Nq+1,dtype=complex)
Om2_M = np.zeros(Nq+1,dtype=complex)
Om1_X = np.zeros(Nq+1,dtype=complex)
Om2_X = np.zeros(Nq+1,dtype=complex)
# MT
for iq in range(0,Nq+1):
qx[iq] = qi+iq*dq
qy[iq] = qi+iq*dq
A = [ [(qx[iq]-1)**2+(qy[iq]-1)**2,0],
[0,qx[iq]**2+qy[iq]**2]]
B = [[e_p0p0,e_p1p0],[e_m1p0,e_p0p0]]
x,E = linalg.eig(A,B)
x.sort()
Om1_M[iq] = cmath.sqrt(x[0])
Om2_M[iq] = cmath.sqrt(x[1])
qM = np.sqrt(qx*qx+qy*qy)
# TX
for iq in range(0,Nq+1):
qx[iq] = qi+iq*dq
qy[iq] = 0.0
A = [ [(qx[iq]-1)**2+qy[iq]**2,0],
[0,qx[iq]**2+qy[iq]**2]]
B = [[e_p0p0,e_p1p0],[e_m1p0,e_p0p0]]
x,E = linalg.eig(A,B)
x.sort()
Om1_X[iq] = cmath.sqrt(x[0])
Om2_X[iq] = cmath.sqrt(x[1])
qX = np.sqrt(qx*qx+qy*qy)
plt.plot(-qM,Om2_M.real,'-r',-qM,Om1_M.real,'-b')
plt.plot(+qX,Om2_X.real,'-r',+qX,Om1_X.real,'-b')
plt.grid()
plt.show()
import numpy as np
import matplotlib.pyplot as plt
import math
N = 10
# Trayectoria X
x = np.zeros(N+1)
Om_p0p0_x = np.zeros(N+1)
Om_m1p0_x = np.zeros(N+1)
Om_p1p0_x = np.zeros(N+1)
Om_p0m1_x = np.zeros(N+1)
Om_p0p1_x = np.zeros(N+1)
Om_p1p1_x = np.zeros(N+1)
Om_m1m1_x = np.zeros(N+1)
for iz in range(0,N+1):
x_x = 0.5*(iz/N)
x_y = 0.0*(iz/N)
x[iz] = math.sqrt(x_x*x_x+x_y+x_y)
Om_p0p0_x[iz] = math.sqrt((x_x+0)**2+(x_y+0)**2)
Om_m1p0_x[iz] = math.sqrt((x_x-1)**2+(x_y+0)**2)
Om_p1p0_x[iz] = math.sqrt((x_x+1)**2+(x_y+0)**2)
Om_p0m1_x[iz] = math.sqrt((x_x+0)**2+(x_y-1)**2)
Om_p0p1_x[iz] = math.sqrt((x_x+0)**2+(x_y+1)**2)
Om_p1p1_x[iz] = math.sqrt((x_x+1)**2+(x_y+1)**2)
Om_m1m1_x[iz] = math.sqrt((x_x-1)**2+(x_y-1)**2)
# Trayectoria M
m = np.zeros(N+1)
Om_p0p0_m = np.zeros(N+1)
Om_m1p0_m = np.zeros(N+1)
Om_p1p0_m = np.zeros(N+1)
Om_p0m1_m = np.zeros(N+1)
Om_p0p1_m = np.zeros(N+1)
Om_p1p1_m = np.zeros(N+1)
Om_m1m1_m = np.zeros(N+1)
for iz in range(0,N+1):
m_x = 0.5*(iz/N)
m_y = 0.5*(iz/N)
m[iz] = -math.sqrt(m_x*m_x+m_y*m_y)
Om_p0p0_m[iz] = math.sqrt((m_x+0)**2+(m_y+0)**2)
Om_m1p0_m[iz] = math.sqrt((m_x-1)**2+(m_y+0)**2)
Om_p1p0_m[iz] = math.sqrt((m_x+1)**2+(m_y+0)**2)
Om_p0m1_m[iz] = math.sqrt((m_x+0)**2+(m_y-1)**2)
Om_p0p1_m[iz] = math.sqrt((m_x+0)**2+(m_y+1)**2)
Om_p1p1_m[iz] = math.sqrt((m_x+1)**2+(m_y+1)**2)
Om_m1m1_m[iz] = math.sqrt((m_x-1)**2+(m_y-1)**2)
plt.plot(x,Om_p0p0_x,'-b',label='(0,0)')
plt.plot(x,Om_m1p0_x,'og',label='(-1,0)')
plt.plot(x,Om_p1p0_x,'-r',label='(+1,0)')
plt.plot(x,Om_p0m1_x,'-c',label='(0,-1)')
plt.plot(x,Om_p0p1_x,'oc',label='(0,+1)')
plt.plot(x,Om_p1p1_x,'-m',label='(+1,+1)')
plt.plot(x,Om_m1m1_x,'-y',label='(-1,-1)')
plt.plot(m,Om_p0p0_m,'-b',label='(0,0)')
plt.plot(m,Om_m1p0_m,'og',label='(-1,0)')
plt.plot(m,Om_p1p0_m,'-r',label='(+1,0)')
plt.plot(m,Om_p0m1_m,'-c',label='(0,-1)')
plt.plot(m,Om_p0p1_m,'oc',label='(0,+1)')
plt.plot(m,Om_p1p1_m,'-m',label='(+1,+1)')
plt.plot(m,Om_m1m1_m,'-y',label='(-1,-1)')
plt.legend()
plt.show()
plt.plot(x,Om_p0p0_x,'.b',label='(0,0)')
plt.plot(x,Om_m1p0_x,'.g',label='(-1,0)')
plt.plot(x,Om_p1p0_x,'.r',label='(+1,0)')
plt.plot(x,Om_p0m1_x,'.c',label='(0,-1)')
plt.plot(x,Om_p0p1_x,'.y',label='(0,+1)')
plt.plot(x,Om_p1p1_x,'.k',label='(+1,+1)')
plt.plot(x,Om_m1m1_x,'ob',label='(-1,-1)')
plt.plot(m,Om_p0p0_m,'.b',label='(0,0)')
plt.plot(m,Om_m1p0_m,'.g',label='(-1,0)')
plt.plot(m,Om_p1p0_m,'.r',label='(+1,0)')
plt.plot(m,Om_p0m1_m,'.c',label='(0,-1)')
plt.plot(m,Om_p0p1_m,'.y',label='(0,+1)')
plt.plot(m,Om_p1p1_m,'.k',label='(+1,+1)')
plt.plot(m,Om_m1m1_m,'ob',label='(-1,-1)')
plt.plot(-qM,Om2_M.real,'-r',-qM,Om1_M.real,'-b')
plt.plot(+qX,Om2_X.real,'-r',+qX,Om1_X.real,'-b')
plt.legend()
plt.show()
