Graphene Nanoribbon Band Structure Plot in Matlab: Bug Fixes Included

[anahita]

Hi
I have plot band structure zigzag graphene nanoribbons with Matlab.
I do not know if it is properly written program anyone can help Bug fixes.

NU=10; % Number of atoms
Nbnd=4*NU; % number of bands
aa=2.26;
a=sqrt(3)*aa;
Csoc=0.0;

X(1)=0;
Y(1)=0;
for ixy=2:NU
  if mod(ixy,4)==2
  X(ixy)=X(ixy-1)-aa*cosd(30);
  Y(ixy)=Y(ixy-1)+aa*sind(30);
  end
  if mod(ixy,4)==3
  X(ixy)=X(ixy-1);
  Y(ixy)=Y(ixy-1)+aa;
  end
  if mod(ixy,4)==0
  X(ixy)=X(ixy-1)+aa*cosd(30);
  Y(ixy)=Y(ixy-1)+aa*sind(30);
  end
  if mod(ixy,4)==1
  X(ixy)=X(ixy-1);
  Y(ixy)=Y(ixy-1)+aa;
  end
end

for iz=1:NU
if mod(iz,2)==1
Z(iz)=0.45;
else
Z(iz)=-0.45;
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);

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.26)<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;
  h=hamiltonian(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)+h*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)
The functions h and H0 attached.
function [h] = hamiltonian(l,m,n)

  h=zeros(4,4);
  tsss=-2.08;
  tsps=2.48;
  tpps=2.72;
  tppp=-0.72;
  
  % gharar dad ----->% S, Px, Py, Pz
  %  S, Px, Py, Pz
  %1  2  3  4  
  h(1,1)=tsss;
  h(1,2)=l*tsps;
  h(1,3)=m*tsps;
  h(1,4)=n*tsps;
  
  h(2,1)=-(l*tsps)';
  h(2,2)=l*l*tpps+(1-l*l)*tppp;
  h(2,3)=l*m*tpps-l*m*tppp;
  h(2,4)=l*n*tpps-l*n*tppp;
  
  h(3,1)=-(m*tsps)';
  h(3,2)=(l*m*tpps-l*m*tppp)';
  h(3,3)=m*m*tpps+(1-m*m)*tppp;
  h(3,4)=m*n*tpps-m*n*tppp;
  
  h(4,1)=-(n*tsps)';
  h(4,2)=(l*n*tpps-l*n*tppp)';
  h(4,3)=(m*n*tpps-m*n*tppp)';
  h(4,4)=n*n*tpps+(1-n*n)*tppp;
end
and
function [H0] = H0(Nbnd)
Es=-4.2;
Ep=1.715;
H0=zeros(Nbnd);
for ih0=1:Nbnd
if mod(ih0,4)==1
  H0(ih0,ih0)=Es;
else
  H0(ih0,ih0)=Ep;
end
end

 

[DrClaude]

Does it do what you expect it to do? Without an explanation of what you are trying to achieve, it is very hard to judge if something is correct or not.

 

[anahita]

Does it do what you expect it to do? Without an explanation of what you are trying to achieve, it is very hard to judge if something is correct or not.

I want to calculate structure graphene nanoribbons under strain, but the density of states thinks is wrong.

 

[hokhani]

I want to calculate structure graphene nanoribbons under strain, but the density of states thinks is wrong.

Your program is not clear. You have to at least clarify which method you have used.

 

[anahita]

Your program is not clear. You have to at least clarify which method you have used.

The program is written using the tight binding model.

 

[hokhani]

The program is written using the tight binding model.

Ok, right. But the commands are not clear and it seems that you take 4 orbitals s and p (in gharadrad) and also some other details are not clear.To interpret a program you should bring the algorithm or the exact physical method you have used. For example I can not realize which bases you have taken into account for your Tight Binding calculation and what is the approximation (nearest neighbor, second,...). If your problem is with programming you should ask it in "Math software and Latex".If you want to check the accuracy of the result you can compare with the prepared results in many sites such as http://demonstrations.wolfram.com/ElectronicBandStructureOfArmchairAndZigzagGrapheneNanoribbon/.

 

[anahita]

Ok, right. But the commands are not clear and it seems that you take 4 orbitals s and p (in gharadrad) and also some other details are not clear.To interpret a program you should bring the algorithm or the exact physical method you have used. For example I can not realize which bases you have taken into account for your Tight Binding calculation and what is the approximation (nearest neighbor, second,...). If your problem is with programming you should ask it in "Math software and Latex".If you want to check the accuracy of the result you can compare with the prepared results in many sites such as http://demonstrations.wolfram.com/ElectronicBandStructureOfArmchairAndZigzagGrapheneNanoribbon/.

What is meant by the following statement? " I can not realize which bases you have taken into account for your Tight Binding calculation"
I have used nearest neighbour approximation.

 

[hokhani]

What is meant by the following statement? " I can not realize which bases you have taken into account for your Tight Binding calculation"
I have used nearest neighbour approximation.

I meant which bases you have taken for the Hamiltonian. As I said earlier, this problem requires that details to be explained entirely.

 

[anahita]

I meant which bases you have taken for the Hamiltonian. As I said earlier, this problem requires that details to be explained entirely.

The unit cell contains 2N0 carbon atoms which leads to a 2N0x2N0 integral matrix H. matrices developed for the unit cell based on the first nearest neighbour assumption.