✅作者简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，matlab项目合作可私信。
🍎个人主页：Matlab科研工作室
🍊个人信条：格物致知。
⛄ 内容介绍
光子晶体具有光子带隙特性的电介质结构,是一种能够操控光的人造物质.具有简便,廉价,高速,大容量和强抗干扰性能的更多光子晶体产品将在各个领域发挥重要作用.
⛄ 部分代码
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% this program calculates and plots the wave-vector diagram (i.e. photonic bands at constant frequency)
%%% for a 2D photonic crystal consisting of cylinders with circular cross-section and
%%% infinite height, arranged in a triangular lattice; oblique propagation is implicit, so
%%% the polarization states cannot be separated in E-pol and H-pol; 'omega'is taken as input;
%%% Fourier coefficients for the expansion of dielectric constant are calculated analytically;
%%% the materials considered here are dielectric and dispersionless;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% the package contains the following programs:
%%% pwem2Db.m - main program
%%% epsgg.m - routine for calculating the matrix of Fourier coefficients
%%% of dielectric function
%%% bz_irr2.m - routine for 2D discretization of irreducible Brillouin zone polygon;
%%% kvect2.m - routine for calculating diagonal matrices with elements
%%% (kx+Gx) and (ky+Gy), where G=(Gx,Gy) is a reciprocal
%%% lattice vector
%%% oblic_eigs.m - routine for solving the eigenvalue problem for
%%% H-field
close all
clear all
tic
omega=0.45; % normalized frequency "a/lambda"
r=0.43; % radius of cylindrical holes (normalized w.r.t. lattice constant "a")
na=1; nb=3.45; % refractive indices (cylinders-atoms, background)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
No1=7; No2=No1;
N1=2*No1+1; N2=2*No2+1;
N=N1*N2; % total number of plane waves used in Fourier expansions
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% primitive vectors of direct lattice (normalized w.r.t. lattice constant "a")
a1=[sqrt(3)/2, -1/2, 0]; a2=[sqrt(3)/2, 1/2, 0];
%%% area of primitive cell
ac=norm(cross(a1,a2));
%%% primitive vectors of direct lattice (normalized w.r.t. lattice constant "2*pi/a"): b1=[1/sqrt(3),-1]; b2=[1/sqrt(3),1];
b1=(1/ac)*[a2(2),-a2(1)]; b2=(1/ac)*[-a1(2), a1(1)];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% matrix of Fourier coefficients
eps1 = feval（'epsgg',r,na,nb,b1,b2,N1,N2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% 2D discretization of irreducible Brillouin zone polygon
Nr=20; % even number
[BZx,BZy]=feval（'bz_irr2', Nr);
kx=[]; ky=[]; kz=[];
S=2.5; % point size for scatter plot
for j=1:length(BZx)
%%% diagonal matrices with elements (kx+Gx) si (ky+Gy)
[kGx, kGy] = feval（'kvect2',BZx(j),BZy(j),b1,b2,N1,N2);
[P, beta]=feval（'oblic_eigs',omega,kGx,kGy,eps1,N);
L=imag(beta)==0;
qp=sort(beta(L)); %%% keep only the propagative modes
display(sprintf('Calculation for k[%d] is finished',j));
for r=1:length(qp)
kx(j,r)=BZx(j); ky(j,r)=BZy(j); kz(j,r)=qp(r);
end
end
M=length(BZx)*length(qp);
scatter3(reshape(kx,1,M), reshape(ky,1,M), reshape(kz,1,M), S,'r','filled'), view(65,10)
title(sprintf('Wavevector diagram for omega=%0.5g',omega));
xlabel('kx'); ylabel('ky'); zlabel('kz');
toc
⛄ 运行结果
⛄ 参考文献
[1]许江勇. 基于Matlab研究光子晶体的特性[J]. 兴义民族师范学院学报, 2012(4):4.
[2]吴炳坚, 沈廷根. 基于Matlab的光子晶体波导仿真研究[J]. 微计算机信息, 2007(01S):3.
