temp2.m

用传输矩阵的方法来计算光子晶体带隙

n=1.4+0.24i;
o=[1:25];
Csca=[];
krange=2*pi./[2:0.02:4];

for k=krange
% k=3;

% p1=besselh(o,1,k)+o./(2*o+1)*k.*besselh(o-1,1,k)-(o+1)./(2*o+1)*k.*besselh(o+1,1,k);
% p2=besselj(o,k)+o./(2*o+1)*k.*besselj(o-1,k)-(o+1)./(2*o+1)*k.*besselj(o+1,k);
% p3=besselj(o,n*k)+(n)*k*(o./(2*o+1).*besselj(o-1,n*k)-(o+1)./(2*o+1).*besselj(o+1,n*k));
% p4=besselj(o,k);
% p5=besselh(o,1,k);
% p6=besselj(o,n*k);

p1=k*besselh(o-1+1/2,1,k)-o.*besselh(o+1/2,1,k);
p2=k*besselj(o-1+1/2,k)-o.*besselj(o+1/2,k);
p3=n*k*besselj(o-1+1/2,n*k)-o.*besselj(o+1/2,n*k);
p4=besselj(o+1/2,k);
p5=besselh(o+1/2,1,k);
p6=besselj(o+1/2,n*k);

a=(n^2*p6.*p2-p4.*p3)./(n^2*p6.*p1-p5.*p3);
b=(p6.*p2-p4.*p3)./(p6.*p1-p5.*p3);
c=(p4.*p1-p5.*p2)./(p6.*p1-p5.*p3);
d=(n*p4.*p1-n*p5.*p2)./(n^2*p6.*p1-p5.*p3);

Csca=[Csca sum(real(2*pi/k^2*(2*o+1).*(real(a+b)-a*conj(a)-b*conj(b))))];
% -(a.*conj(a)+b.*conj(b))
end

plot(2*pi./krange,Csca/pi,'r')