%**********************%主函數(shù)*****************************************
function main()
global chrom lchrom oldpop newpop varible fitness  popsize  sumfitness %定義全局變量
global pcross pmutation temp bestfit maxfit gen bestgen
global maxgen po pp mp np
lchrom=12;         %染色體長度
popsize=30;        %種群大小
pcross=0.6;        %交叉概率
pmutation=0.1;     %變異概率
maxgen=100;        %最大代數(shù)
po=0.1;            %淘汰概率
pp=0.1;            %保護概率
mp=floor(pp*popsize);  %保護的個數(shù)
np=floor(po*popsize);  %淘汰的個數(shù)
%
initpop;              % 初始種群
%
for gen=1:maxgen
    generation;
end
%
best;
bestfit       % 最佳個體適應度值輸出
bestgen       % 最佳個體所在代數(shù)輸出
gen=1:maxgen;
plot(gen,maxfit(1,gen));  % 進化曲線
%
%**********************  產(chǎn)生初始種群  ************************************ 
%
function initpop()
global lchrom oldpop popsize chrom
for i=1:popsize
    chrom=rand(1,lchrom);    %  lchrom=12  染色體長度
    for j=1:lchrom
        if chrom(1,j)<0.5
            chrom(1,j)=0;
       else 
           chrom(1,j)=1;
        end
    end
    oldpop(i,1:lchrom)=chrom;    
end
%
%*************************%產(chǎn)生新一代個體**********************************
%
function generation()
objfun;                  %計算適應度值
pp_po;                   %保護淘汰操作
select;                  %選擇操作
crossover;
mutation;
%
%************************%計算適應度值************************************
%
function objfun()
global lchrom oldpop fitness popsize chrom maxfit gen varible
a=0;
%b=100;
b=30;
for i=1:popsize
    chrom=oldpop(i,:);
    c=decimal(chrom);
    varible(1,i)=a+c*(b-a)/(2.^lchrom-1);  %對應變量值
    fitness(1,i)=varible(1,i)*varible(1,i);
%    fitness(1,i)=20*cos(0.25*varible(1,i))-12*sin(0.33*varible(1,i))+40;  %個體適應度函數(shù)值
end
%
lsort;                          % 個體排序
%
maxfit(1,gen)=max(fitness);     %求最大適應度值
%
%************************%二進制轉(zhuǎn)十進制**********************************
%
function c=decimal(chrom)
global lchrom popsize
c=0;
for j=1:lchrom
    c=c+chrom(1,j)*2.^(lchrom-j);
end
%
%*************************%個體排序****************************************
%
function lsort()
global popsize fitness oldpop
for i=1:popsize
    j=i+1;
    while j<=popsize
        if fitness(1,i)>fitness(1,j)
            tf=fitness(1,i);
            tc=oldpop(i,:);
            fitness(1,i)=fitness(1,j);
            oldpop(i,:)=oldpop(j,:);
            fitness(1,j)=tf;
            oldpop(j,:)=tc;
        end
        j=j+1;
    end
end
%**************************************************************************
%
%*************************%保護/淘汰操作************************************
function pp_po()
global popsize oldpop np
i=np+1;                      % np=floor(po*popsize);  %淘汰的個數(shù)
j=1;
while i<=popsize
    toldpop(j,:)=oldpop(i,:);
    j=j+1;
    i=i+1;
end
for i=1:(popsize-np)
    oldpop(i,:)=toldpop(i,:);
end
%**************************************************************************
%*************************%轉(zhuǎn)輪法選擇操作**********************************
function select()
global fitness popsize sumfitness oldpop temp mp np
sumfitness=0;
for i=1:(popsize-np-mp)
    sumfitness=sumfitness+fitness(1,i);
end
for i=1:(popsize-mp-np)
    p(1,i)=fitness(1,i)/sumfitness;   % 個體染色體的選擇概率
end
q=cumsum(p);                          % 個體染色體的累積概率
b=sort(rand(1,(popsize-mp)));         % mp保護個體數(shù)
j=1;
k=1;
while j<=(popsize-mp)
    if b(1,j)<q(1,k)
       temp(j,:)=oldpop(k,:);
       j=j+1;
   else
       k=k+1;
   end
end
for i=(popsize-mp+1):popsize
    j=popsize-np-mp+1;
    temp(i,:)=oldpop(j,:);
    j=j+1;
end
%**************************************************************************
%**************************%交叉操作***************************************
function crossover()
global temp popsize pcross lchrom mp
n=floor(pcross*(popsize-mp));           %交叉發(fā)生的次數(shù)
if rem(n,2)~=0
    n=n+1;
end
j=1;
m=0;
for i=1:(popsize-mp)
    p=rand;
    if p<pcross
        parent(j,:)=temp(i,:);      %產(chǎn)生兩個父代
        k(1,j)=i;
        j=j+1;
        m=m+1;
        if (j==3)&(m<=n)
            pos=round(rand*(lchrom-1))+1;
            for i=1:pos
                child1(1,i)=parent(1,i);
                child2(1,i)=parent(2,i);
            end
            for i=(pos+1):lchrom
                child1(1,i)=parent(2,i);
                child2(1,i)=parent(1,i);
            end
            i=k(1,1);
            j=k(1,2);
            temp(i,:)=child1(1,:);
            temp(j,:)=child2(1,:);
        end
        j=1;
    end
end
%**************************************************************************
%****************************%變異操作*************************************
function mutation()
global popsize lchrom pmutation temp newpop oldpop mp 
m=lchrom*(popsize-mp);               %總的基因數(shù)
n=round(pmutation*m);           %變異發(fā)生的次數(shù)
for i=1:n
    k=round(rand*(m-1))+1;
    j=ceil(k/lchrom);
    l=rem(k,lchrom);
    if l==0
       temp(j,lchrom)=~temp(j,lchrom);
   else
       temp(j,l)=~temp(j,l);
   end
end
for i=1:popsize
    newpop(i,:)=temp(i,:);     %產(chǎn)生新的個體
    oldpop(i,:)=newpop(i,:);
end
%**************************************************************************
%*********************%最佳個體********************************************
function best()
global maxfit bestfit gen maxgen bestgen
bestfit=maxfit(1,1);
gen=2;
while gen<=maxgen
      if maxfit(1,1)<maxfit(1,gen)
         bestfit=maxfit(1,gen);
         bestgen=gen;
     end
     gen=gen+1;
 end
%**************************************************************************
function[sol,eval]=fitness(sol,options) 
x=sol(1); 
eval=x+10*sin(5*x)+7*cos(4*x);
%eval為個體適應度
%遺傳算法求最大值
figure(1)
fplot('x+10*sin(5*x)+7*cos(4*x)',[0,9])
hold on
initPop=initializega(10,[0 9],'fitness'); 
figure(1)
plot(initPop(:,1),initPop(:,2),'b*')
[x endPop]=ga([0 9],'fitness',[],initPop,[1e-6 1 1],'maxGenTerm',25,'normGeomSelect',... 
     [0.08],['arithXover'],[20],'nonUnifMutation',[2 25 3]) 
figure
plot(endPop(:,1),endPop(:,2),'bo')
[x endPop,bPop,trace]=ga([0 9],'fitness',[],initPop,[1e-6 1 1],'maxGenTerm',25,'normGeomSelect',... 
     [0.08],['arithXover'],[2],'nonUnifMutation',[2 25 3]) 
figure
plot(endPop(:,1),endPop(:,2),'y*')
figure
plot(trace(:,1),trace(:,3),'y-')
hold on
plot(trace(:,1),trace(:,2),'r-')
xlabel('Generation');ylabel('Fittness');
legend('解的變化','種群平均值的變化');