二维光子晶体禁带的遗传优化算法MATLAB源码

1、二维光子晶体禁带的遗传优化算法MATLAB源码 光子晶体中因周期性结构而存在的频率禁带称为光子禁带，光子禁带的存在是光子晶体具有广泛应用前景的重要原因。 禁带越大，可控光的频带也越宽，因此如何设计合适的晶体结构以获得大的光子禁带一直是研究和应用的重点与热点。本程序实现了参考文献快速遗传算法优化计算二维光子晶体，是GreenSim团队原创作品，转载请注明。function Zp,Yp,Xp,LC1,LC2=PBGGA(M,N,Pm)%-% PBGGA.m% Photonic Band Gap Genetic Algorithm% 二维光子晶体禁带的遗传优化算法% GreenSim团队原创作品，转

2、载请注明% 欢迎访问GreenSim算法仿真团队%-% 函数功能% 在平面波展开法的基础上，使用遗传算法寻找具有最大光子禁带的二维光子晶体结构%-% 主要参考文献% 龚春娟,胡雄伟.快速遗传算法优化计算二维光子晶体J.% 半导体学报.2006,27(6):1098-1102%-% 输入参数列表% M 遗传进化迭代次数% N 种群规模(取偶数)% Pm 变异概率%-% Zp 最大禁带% Yp 最大禁带的区间上下限% Xp 光子晶体结构决策变量（10×10的01矩阵）% LC1 收敛曲线1，各代最优个体适应值的记录% LC2 收敛曲线2，各代群体平均适应值的记录%-%第一步：变量初始化e

3、psa=11.56;%介质介电常数epsb=1;%空气介电常数Xp=zeros(10,10);%最优结构Yp=zeros(1,2);%禁带上下限Zp=0;%禁带宽度LC1=zeros(1,M);%收敛曲线1LC2=LC1;%收敛曲线2%第二步：随机产生初始钟群farm=cell(1,N);%采用细胞结构存储种群for i=1:NA=unidrnd(2,5,10)-ones(5,10);%由于基本单元结构是中心对称的，所以只需要对矩阵的上半部分进行初始化B=flipud(fliplr(A);%剩余部分先做上下倒置，再做左右旋转 C=A;B;%把两步分合并farmi=reshape(C,100,1

4、);%个体被存储在细胞结构farm中endcounter=0;%设置迭代计数器while counter<M%停止条件为达到最大迭代次数%第三步：交叉newfarm=cell(1,N);%交叉产生的新种群存在这里面Ser=randperm(N);%用这个函数保证随机配对for i=1:2:(N-1)AA=reshape(farmSer(i),10,10);% 待交叉的两个父代个体 BB=reshape(farmSer(i+1),10,10);A=AA(1:5,:);%为保证对称性，需先提取其中一半的基因B=BB(1:5,:);Manner=unidrnd(2);%随机选择交叉方式if M

5、anner=1%横交叉cp=unidrnd(4);%随机选择交叉点a=A(1:cp,:);B(cp+1):5,:);%双亲双子单点交叉b=B(1:cp,:);A(cp+1):5,:);else%列交叉cp=unidrnd(9);%随机选择交叉点a=A(:,1:cp),B(:,(cp+1):10);%双亲双子单点交叉b=B(:,1:cp),A(:,(cp+1):10);enda=a;flipud(fliplr(a);%反演对称b=b;flipud(fliplr(b);newfarmi=reshape(a,100,1);%交叉后的子代存入newfarm newfarmi+1=reshape(b,1

6、00,1);end%新旧种群合并FARM=farm,newfarm;%第四步：选择复制FITNESS=zeros(1,2*N);fitness=zeros(1,N);for i=1:(2*N)X=FARMi;Band=pixelTM(X,10,epsa,epsb,0);%调用平面波展开法的程序，中间结果不画图B1=Band(:,1);B2=Band(:,2);for j=1:length(B1)if B1(j)>1|B2(j)>1B1(j)=0;B2(j)=0;endendBB=abs(B1-B2);pos=find(BB=max(BB);Z=BB(pos(1);%最宽的禁带作为适

7、应值FITNESS(i)=Z;endSer=randperm(2*N);%选择复制采取两两随机配对竞争的方式，具有保留最优个体的能力for i=1:Nf1=FITNESS(Ser(2*i-1);f2=FITNESS(Ser(2*i);if f1>=f2farmi=FARMSer(2*i-1);fitness(i)=FITNESS(Ser(2*i-1);elsefarmi=FARMSer(2*i);fitness(i)=FITNESS(Ser(2*i);endend%记录最佳个体和收敛曲线maxfitness=max(fitness)meanfitness=mean(fitness)LC1

8、(counter+1)=maxfitness;%收敛曲线1，各代最优个体适应值的记录 LC2(counter+1)=meanfitness;%收敛曲线2，各代群体平均适应值的记录 pos=find(fitness=maxfitness);xp=farmpos(1);Zp=maxfitness;Xp=reshape(xp,10,10);%第五步：变异for i=1:Nif Pm>rand;%变异概率为PmX=farmi;XX=reshape(X,10,10);%为保证对称性，需先将个体基因变形 XA=XX(1:5,:);posX=unidrnd(5,1,1);%随机选取变异位posY=un

9、idrnd(10,1,1);XA(posX(1),posY(1)=mod(XA(posX(1),posY(1)+1,2);%变异，0变成1，1编程0%XA(posX(2),posY(2)=mod(XA(posX(2),posY(2)+1,2); XB=flipud(fliplr(XA);XX=XA;XB;farmi=reshape(XX,100,1);%还原个体形态endendfarmpos(1)=reshape(Xp,100,1);counter=counter+1end%输出结果并绘图figure(1);xp=reshape(Xp,100,1);Band=pixelTM(xp,10,epsa,epsb,1);%画最优个体的禁带图B1=Band