Matlab code for calculate band structure graphene nanoribbons under strain is below:
clear;
clc;
close all;
NU=12; % Number of atoms in unit cell
Nbnd=4*NU; % number of bands
q=0.03;
w=0.41;
aa=2.232*(1+q);
a=3.866*(1+q);
X(1)=1.9330*(1+q);
Y(1)=0;
xswitch = 0;
for ixy=2:NU
if mod(ixy,2)==1
Y(ixy)=Y(ixy-1)+aa;
else
Y(ixy)=Y(ixy-1)+aa*sind(30);
X(ixy) = xswitch;
if (ixy+1)<=NU
X(ixy+1)=xswitch;
end
if xswitch == 0
xswitch =aa*cosd(30);
else
xswitch = 0;
end
end
end
for iz=1:NU
if mod(iz,2)==1
Z(iz)=0.46152;
else
Z(iz)=0;
end
end
sho=0;
for is=[0,-1,1]
for ks=1:NU
sho=sho+1;
XT(sho)=X(ks)+is*a;
YT(sho)=Y(ks);
ZT(sho)=Z(ks);
Ax(sho)=is*a;
No(sho)=ks;
end
end
figure(1)
plot(XT,YT,'*')
Ax=Ax/a;for ik=1:101
K(ik)=(-pi+(ik-1)*((2*pi)/100))*(1-q);
H=H0(Nbnd);
for is=1:NU
for js=1:sho
dis=sqrt(((XT(is)-XT(js))^2)+((YT(is)-YT(js))^2));
if abs(dis-2.232)<0.1 & abs(No(is)-No(js))>0
l=(XT(is)-XT(js))/dis;
m=(YT(is)-YT(js))/dis;
n=(ZT(is)-ZT(js))/dis;
h1=hamiltonian1(l,m,n);
h2=hamiltonian2(l,m,n);
H((No(is)-1)*4+1:No(is)*4,(No(js)-1)*4+1:No(js)*4)=H((No(is)-1)*4+1:No(is)*4,(No(js)-1)*4+1:No(js)*4)+h1*exp(i*K(ik)*Ax(js))+h2*exp(i*K(ik)*Ax(js));
end
end
end
E(ik,1:Nbnd)=sort(real(eig(H)));pl(ik)=(ik-1)/100;
end
figure(2)
plot(E)
[V,D]=eig(H);
f=diag(D);
g=f<0;
g1=f(g);
r=sum(g1)