%遺傳算法
% function [out]=GA1(in)

clear all;
close all;
tic
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                % 參數設置 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
M1 = 4;  
M2 = 4;
M3 = 4;    % 我方各載機所攜帶的武器數

total_W = M1+M2+M3;        % 總的武器數12
total_T = 10;                 % 總的目標數10

% 迎擊武器單一命中概率,行為目標，列為武器 
PP = [0.10  0.05  0.36  0.73  0.21  0.70  0.02  0.85  0.90  0.10  0.08  0.50 ;
      0.56  0.56  0.75  0.31  0.38  0.06  0.33  0.32  0.95  0.37  0.50  0.87;
      0.80  0.25  0.21  0.67  0.40  0.20  0.12  0.50  0.45  0.60  0.30  0.60;
      0.40  0.17  0.10  0.52  0.23  0.80  0.75  0.73  0.50  0.80  0.08  0.10;
      0.32  0.20  0.43  0.26  0.45  0.15  0.62  0.85  0.95  0.01  0.90  0.70;
      0.20  0.10  0.82  0.32  0.30  0.40  0.03  0.34  0.26  0.69  0.70  0.87;
      0.96  0.58  0.75  0.80  0.01  0.63  0.27  0.10  0.72  0.23  0.01  0.07;
      0.40  0.67  0.50  0.91  0.34  0.54  0.90  0.22  0.65  0.23  0.90  0.92;
      0.03  0.42  0.15  0.32  0.40  0.30  0.74  0.58  0.85  0.90  0.30  0.55;
      0.15  0.75  0.62  0.45  0.56  0.80  0.61  0.76  0.10  0.06  0.35  0.50];    
% 目標對飛機的威脅評估因子，行為飛機，列為目標  
TH = [0.5485  0.4987  0.4680  0.3867  0.4161  0.4529  0.2943  0.3466  0.2708  0.3320 ;
      0.5390  0.4881  0.4521  0.3418  0.3829  0.4652  0.2562  0.3463  0.2935  0.2548 ;
      0.5541  0.5402  0.4976  0.4765  0.4921  0.5654  0.3178  0.3735  0.3012  0.3030 ];   
  
% 求每個目標對我機群的總威脅，即TH每列求和     
TH_Z  = zeros(1,total_T);
for i = 1:total_T
    for j = 1:3
        TH_Z(i) = TH_Z(i)+TH(j,i);
    end
end                                    
% 對目標進行排序，按照其對我機群的威脅從大到小，輸出為敵機序號
mb_px = mbpx(TH_Z,total_T);          
%[order,index]=sort(TH_Z);
%index
weapon = [1:total_W];          % 武器列表

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                % 初始化種群 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
popsize=80;                        %種群規模
len=total_W;                       %染色體長度
maxgen=100;                         %最大迭代次數
pcross=0.72;                        %交叉概率
pmutation=0.08;                    %變異概率

%初始化種群
for popsiz=1:popsize
    t_num = total_W;     % 目標數
    target = [1:total_T,mb_px(1:(total_W-total_T))];   % 目標列表；每機一個多余導彈加分給威脅最大的目標
    for l=1:total_W
        m = rand_function(1,t_num);                    % 產生一個1-t_num之間的隨機數
        pop(popsiz,l)=target(m);                       % 列數即12個武器，依次排列
        target = [target(1:m-1),target(m+1:t_num)];    % 將目標列表中的第m個排除
        t_num = t_num-1;                               % 目標個數減去1
    end
end        

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                %遺傳結束條件 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
mingen=15;              %最小迭代次數
gen_min_tatio=0.02;     %子群體最小進化率
gen_die=5;              %子群體進化率小于最小進化率的可容忍進化代數
t=0;                    %子群體進化率小于最小進化率的當前代數
        
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                % 迭代尋優 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for gen=1:maxgen                   %執行至設定的代數終止
    [mpop,npop]=size(pop);         %取種群矩陣的行數和列數，即種群規模和染色體長度
    
cost = zeros(1,mpop); %失敗威脅值初始化
 %%%%%%%%%%%% 計算每個方案的失敗威脅概率和  %%%%%%%%%%%%%%
  for m=1:mpop  
    for kk = 1:npop             
        cost(m) = cost(m)+TH_Z(pop(m,kk))*(1-PP(pop(m,kk),kk)); %每個目標分一個導彈
    end 
         
    for i = 1:npop-1                 
        for j = i+1:npop
            if pop(m,i)==pop(m,j)       %當有一個目標分配兩個導彈時
               cost(m) = cost(m)-TH_Z(pop(m,i))*(1-PP(pop(m,i),i))-TH_Z(pop(m,j))*(1-PP(pop(m,j),j))...
                        +TH_Z(pop(m,i))*(1-PP(pop(m,i),i))*(1-PP(pop(m,j),j));
            end
        end
    end
 
 c=10;                      %適配值函數f(x)中的常數c
 f(m)=c/(0.001+cost(m));           %第m個方案的適配值函數值   
 
  end
  
%保存每一代最好的結果
 [cost_order,cost_index]=sort(cost);
 cost_min(gen)=cost_order(1);
 result(gen,:)=pop(cost_index(1),:);  
 
%結束的條件判斷 
 
      if (gen>mingen)&((cost_min(gen-1)-cost_min(gen))<gen_min_tatio*cost_min(gen-1))
          t=t+1;
      else
          t=0;
      end

  
     if t>=gen_die
        break;
     end

    %輪盤選擇
    selectpop=select2(f,pop);
    %交叉操作策略B
    crosspop=crossover2(pcross,selectpop); 
    %逆轉變異操作
    pop=mutation2(pmutation,crosspop);

end  

%%%%%  輸出最優分配  %%%%%
time1=toc;   %測CPU時間
gen
cost_min(gen)  %遺傳結束時的最好評估值結果

MM = 10;    % 螞蟻數m 
for k=1:MM     %得到前MM個最好值
pop_choose(k,:)=pop(cost_index(k),:);
end
cost_min_2=ant3_10_GAAA(pop_choose);  %轉向蟻群尋優

% final_result = result(gen,:);  %排序結果
load parameter
time=time1+time2  %CPU總耗時
die=gen+nc        %總代數
cost=[cost_min,cost_min_2];   %最優解序列
min=cost(die+1)

 T=[1:(die+1)];   %畫圖
plot(T,cost);
xlabel('迭代次數');
ylabel('最佳評估值');
%plot(T,result(:,1));