covrp_ga_good9.dpr

OVRP问题的遗传算法研究实现

   {判断某一整数是否在整数数组中}
   Function judgein(const u4 : integer;a3: array of integer):boolean;
   var
   u5: integer;
   begin
      result:=true;
      for u5:=low(a3) to high(a3) do
        begin
          if a3[u5]=u4 then
           begin
             Result:=false;
             exit;
           end
            else Result:=true ;
        end;
   end;


   {*PMX*}
   procedure PMXcrosser(VAR a1,a2:array of integer);
   var
   u1,u2,u3,m: integer;
   temp1,temp2: array of integer;
   begin
      Randomize;
      repeat          {*随机选择交叉点 *}
        u1:=randomrange(low(a1)+1,high(a1));
        u2:=randomrange(low(a1)+1,high(a1));
      until (Abs(u1-u2)>0);
      if u1>u2 then
         begin
           m:=u1;
           u1:=u2;
           u2:=m;
         end;
       u3:=low(a1)+1;
     for m:=low(a1) to high(a1) do
        if a1[m]=0  then
           begin
             a1[m]:=u3+n;
             inc(u3);
           end;

     u3:=low(a2)+1;
     for m:=low(a2) to high(a2) do
        if a2[m]=0  then
           begin
             a2[m]:=u3+n;
             inc(u3);
           end;
      setlength(temp1,u2-u1+1);
      setlength(temp2,u2-u1+1);
      for u3:=u1 to u2 do         {*将选中部分放入数组存储*}
        begin
         temp1[u3-u1]:=a1[u3];
         temp2[u3-u1]:=a2[u3];

        end;
      for u3:=low(a1) to high(a1) do         {先替代重复部分 }
         for m:=low(temp2) to high(temp2) do
             if a1[u3]=temp2[m] then
                   begin
                      if judgein(a1[u3],temp1) then a1[u3]:=temp2[m];
                      break;
                   end;
      for u3:=low(a2) to high(a2) do         {先替代重复部分 }
         for m:=low(temp1) to high(temp1) do
             if a2[u3]=temp1[m] then
                   begin
                      if judgein(a2[u3],temp2) then a2[u3]:=temp1[m];
                      break;
                   end;
      for u3:=u1 to u2 do       {替换交叉部分}
        begin
          m:=a1[u3];
          a1[u3]:=a2[u3];
          a2[u3]:=m;
        end;
      for m:=low(a1) to high(a1) do            { 将原来替换为大于n的自然数的车场代号换回0}
        begin
           if a1[m]>n  then      a1[m]:=0;
           if a2[m]>n  then      a2[m]:=0;
        end;
   end;























  {*主程序开始*}

  BEGIN
                             {*read in data *}
    write('Which COVRP problem (1...5,11,12,15,16)? ');
    readln(problem);

  case problem of
    1 : AssignFile(data, 'C:\OVRP\vrpnc1.txt');
    2 : AssignFile(data, 'C:\OVRP\vrpnc2.txt');
    3 : AssignFile(data, 'C:\OVRP\vrpnc3.txt');
    4 : AssignFile(data, 'C:\OVRP\vrpnc4.txt');
    5 : AssignFile(data, 'C:\OVRP\vrpnc5.txt');
   11 : AssignFile(data, 'C:\OVRP\vrpnc11.txt');
   12 : AssignFile(data, 'C:\OVRP\vrpnc12.txt');
   15 : AssignFile(data, 'C:\OVRP\vrpnf11.txt');
   16 : AssignFile(data, 'C:\OVRP\vrpnf12.txt');
  end;

  reset(data);
  readln(data,n,C);
   setlength(x,n+1);
   setlength(y,n+1);
   setlength(d,n+1);
  for i:=0 to n do
    readln(data,x[i],y[i],d[i]);
  CloseFile(data);

  ds:=0;
  for i:=0 to n do
    ds:=ds+d[i];
  kmin:=trunc(ds/c)+1;

  setlength(costmatrix, n+1,n+1);
  for i:=0 to n do
      for j:=0 to n do
          if i=j then costmatrix[i,j]:=maxint
          else if i>j then costmatrix[i,j]:=costmatrix[j,i]
               else  costmatrix[i,j]:=sqrt(sqr(x[i]-x[j])+sqr(y[i]-y[j]));

   //output the distance matrix
  case problem of
    1 : AssignFile(Matrix,'c:\ovrp\MatrixC1.txt');
    2 : AssignFile(Matrix,'c:\ovrp\MatrixC2.txt');
    3 : AssignFile(Matrix,'c:\ovrp\MatrixC3.txt');
    4 : AssignFile(Matrix,'c:\ovrp\MatrixC4.txt');
    5 : AssignFile(Matrix,'c:\ovrp\MatrixC5.txt');
   11 : AssignFile(Matrix,'c:\ovrp\MatrixC11.txt');
   12 : AssignFile(Matrix,'c:\ovrp\MatrixC12.txt');
   15 : AssignFile(Matrix,'c:\ovrp\MatrixF11.txt');
   16 : AssignFile(Matrix,'c:\ovrp\MatrixF12.txt');
  end;
  rewrite(Matrix);
  write(Matrix,'   ');
  for i := 0 to n do write(Matrix,i:7);
  writeln(Matrix);
  for i := 0 to n do
  begin
    write(Matrix,i:3);
    for j := 0 to n do
    begin
      if i=j then write(Matrix,'0':7)          {*** t[i,i]=M **}
      else write(Matrix,costmatrix[i,j]:7:2);
    end;
    writeln(Matrix);
  end;
  CloseFile(Matrix);

{**下面开始随机生成初始可行解*}
    AssignFile(result,'c:\ovrp\result.txt');
    rewrite(result);
    writeln(result, 'Generate initial solution Randomly.');
    case problem of
     1 : writeln(result,'Problem C1');
     2 : writeln(result,'Problem C2');
     3 : writeln(result,'Problem C3');
     4 : writeln(result,'Problem C4');
     5 : writeln(result,'Problem C5');
    11 : writeln(result,'Problem C11');
    12 : writeln(result,'Problem C12');
    15 : writeln(result,'Problem F11');
    16 : writeln(result,'Problem F12');
    end;
   writeln(result,'Kmin=',Kmin);
   writeln(result);writeln(result,' ':12,'K':4,'cost':14,'CPU':8,'':8,'pop_int':7,'times':7);
   CloseFile(result);
FOR times:=1 to 20 do
BEGIN
     Initialize(s1);
     setlength(s1.route, 2*n);
     s1.sum_cost:=maxint;
     s1.sum_ec:=0;
     s1.sum_k:=n;
     for i:=0 to 30 do
      begin
        Initialize(s2);
        s2.sum_cost:=0;
        s2.sum_k:=0;
        setlength(s2.route, 2*n);
        RandomInitial_FS(n,s2);


        {**对返回的可行解s2进行判断，如优于s1中存放的解则取代之*}
        flag:=false;{**用flag标记是否对s1中的解进行修改*}
        if s2.sum_k<s1.sum_k then flag:=true
                             else if (s2.sum_k=s1.sum_k) and (s2.sum_cost<s1.sum_cost)
                                  then flag:=true;
        if flag then
        begin
             s1.sum_cost:=s2.sum_cost;
             s1.sum_k:=s2.sum_k;
             s1.route:=copy(s2.route,0,s2.sum_k+n);
        end {**if flag then*}

      end;{**for*}

     {**把最好的初始解写入最好解记录*}


       setlength(best.best_route,s1.sum_k+n);
       best.best_route:=copy(s1.route,0,s1.sum_k+n);
       best.best_cost:=s1.sum_cost;
       best.best_k:=s1.sum_k;
       best.cputime:=0;




     {**下面开始调用遗传算法 *}
     cpustart:=gettime;  {**开始计时*}
     change_k:= false;    {**表示最优解的车辆数是否改变*}
     POP_int:=0;         {**种群进化代数初始化*}
     nochange_s:=0;      {**最优解连续不改变的次数置0*}
     flag_nochange_s:=1 ; { **对解扰动控制变量置假  **}

     setlength(pop,pop_n);

     {**下面的repeat为整个遗传算法的核心部分，即运行至满足终止条件时算法终止*}
       repeat
         {**以当前最好解的车辆数重新生成初始解种群*}
         change_k:= false;    {**表示最优解的车辆数是否改变*}
         setlength(pop[0].route,length(best.best_route));
         pop[0].route:=copy(best.best_route);
         pop[0].sum_cost:=best.best_cost;
         pop[0].sum_k:=best.best_k;
         pop[0].adapt:=best.adapt;

         for i:=low(pop)+1 to high(pop) do
           begin
             setlength(pop[i].route,n+best.best_k);
             RandomInitial_half_FS(n,best.best_k,pop[i]);
             if 0.5>random then Sum_solution(pop[i].route,pop[i].sum_cost,pop[i].sum_ec,pop[i].adapt,pop[i].sum_k)
              else begin
                     Sum_solution2(pop[i].route,pop[i].sum_cost,pop[i].sum_ec,pop[i].adapt,pop[i].sum_k);
                   end;
           end;









         {**当最好解车辆数不发生改变时进行正常的遗传算法操作*}
         while (nochange_s<max_nochange_s)  and  (pop_int< max_popint) and (change_k=false) do
           begin



             punish_m:=20;{sqr(best.best_cost/best.best_k)**惩罚系数* }
             {**繁殖*}
             Reproduce2(pop,punish_m);



             {计算平均适应度和最大适应度}
             adaptmax:=0;
             adaptavg:=0;
             for i:=low(pop) to high(pop) do
             begin
               if pop[i].adapt>adaptmax then adaptmax:=pop[i].adapt;
               adaptavg:=adaptavg+pop[i].adapt;
             end;
             adaptavg:=adaptavg/length(pop);




             {**随机选出两个解作为父代即配对:方法是将解的序号0...（pop_n－1)组成
             的数组随机变换位置，然后每次取新数组的连续的2个元素值作为父
             代解的编号，如：
             假设 排序后的数组为 3 4 8 2 1 9 0 5 6 7
             则 原种群中的第3、4解为一对父代，第8、2解为一对父代.....第
             6、7为一对父代*}
             setlength(choice_array,pop_n);
             for i:=0 to (pop_n-1) do  choice_array[i]:=i;
             Randomline(choice_array,0,0);
             {**交叉*}



             randomize;
             crossAint:=Randomrange(3,8);
             i:=low(pop);

             while i<(high(pop)) do
                begin
                   pc:=0.8;{pcfunction(pop[choice_array[i]].adapt,pop[choice_array[i+1]].adapt);}
                   if pc>Random then
                     begin
                        tempadapt1:=pop[choice_array[i]].adapt;
                        tempadapt2:=pop[choice_array[i+1]].adapt;
                        if abs(tempadapt1-tempadapt2)=0  then
                                      begin
                                       if random<0.8 then
                                           begin
                                             singlecrosser(pop[choice_array[i]]);
                                             singlecrosser(pop[choice_array[i+1]]);
                                           end;
                                      end
                                      else  CrossA2(pop[choice_array[i]].route,pop[choice_array[i+1]].route,CrossAint)

                     end;
                   i:=i+2;


                end;

             {**突变*}

            if ((nochange_s+1) div 20)=0 then    { 打乱当前种群，用重新注入新解的方式获得基因}

               change_k:= true;

            for i:=low(pop) to high(pop) do
              begin
               pm:=pmfunction(pop[i].adapt);
               Sum_solution(pop[i].route,pop[i].sum_cost,pop[i].sum_ec,pop[i].adapt,pop[i].sum_k);
               if pm>Random then  downmutation2(pop[i]);
              end;

            Evaluation_s2(pop);

             for i:=low(pop) to high(pop) do
              begin
               pm:=pmfunction(pop[i].adapt);
               if pm>Random then  mutation(pop[i].route);
              end;

             

               


              
             for i:=low(pop) to high(pop) do
               if 1>random then Sum_solution(pop[i].route,pop[i].sum_cost,pop[i].sum_ec,pop[i].adapt,pop[i].sum_k)
                   else  Sum_solution2(pop[i].route,pop[i].sum_cost,pop[i].sum_ec,pop[i].adapt,pop[i].sum_k);




             Evaluation_s(pop);

             {**输出结果到屏幕便于监控*}
          {if ((pop_int mod 20)=1) then  }
            begin
              writeln(' ':12,best.best_k:4,best.best_cost:14:2,best.cputime:10:2,
                   '':3,pop_int:7,times:7);
              Append(result);

              writeln(result,' ':12,best.best_k:4,best.best_cost:14:2,best.cputime:10:2,
                   '':3,pop_int:7,times:7);
              CloseFile(result);
            end;

        cpuover:=gettime;
        best.cputime:=(cpuover-cpustart) * 86400;


        


         end;{**end while*}
       until (nochange_s>=max_nochange_s)  or  (pop_int>= max_popint);
       {**释放内存与最终管理*}


        {**下面两个分别用于统计未收敛和早熟收敛的次数*}
        if pop_int >= (max_popint-100) then noconvergence:=noconvergence+1;
        if pop_int < min_popint then prematurity:=prematurity+1;
        cpuover:=gettime;
        best.cputime:=(cpuover-cpustart) * 86400;
        Append(result);
        writeln(result);
        writeln(result,' ':12,'k=',best.best_k:2,'':10,'cost=',best.best_cost:7:2,'':10,'cpu=',best.cputime:5:2);
        j:=0;
        for i:=low(best.best_route) to high(best.best_route)  do
          begin
            if best.best_route[i]=0 then
                        begin
                           j:=j+1;
                           writeln(result);
                           write(result,'route':10,j,':':3,best.best_route[i]);
                        end
                                    else write(result,' ':8,best.best_route[i]:2);
          end;

        writeln(result);
     
        closefile(result);

END;{times}



readln;





end.