covrp_ga_good9.dpr

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

program Covrp_ga_good9;

{$APPTYPE CONSOLE}


uses
  SysUtils,
  Math,
  system;

const
  pop_n=100; {*群体中的个体数目*}
  max_popint=1000; {*最大进化代数*}
  max_nochange_s=200; {*最大的连续最好解不变更代数*}
  min_popint=50 ;{*最小进化代数，用来判断是否发生早熟收敛*}
  check=0;  {*控制是否检验交叉过程*}
type
  array1=array of integer;
  array2=array of double;

  solution_Rec = Record
     route : array1; {*解路径，由表示车场的0和表示客户点自然数组成的序列*}
     sum_ec : double;  {*超出运载能力C的部分的总和*}
     sum_cost : double; {*总路程花费，这里仅指路程总长度*}
     sum_k : integer; {*总车辆数*}
     adapt : double;{*适应度值×}
  end;

  Best_s_Rec = Record
     best_route:array1;
     best_cost: double;
     best_k: integer;
     adapt : double;
     cputime: double;
  end;

  POP_array=array of solution_rec; {*表示一个种群*}
  matrix_array=array of array of double; {*表示距离矩阵*}

var
   data,matrix,result: textfile;
   pop: POP_array;
   choice_array:array1;
   costmatrix:matrix_array;
   best:best_s_Rec;
   s1,s2:solution_Rec;
   x,y,d:array2;
   i,j,k,kmin,n,nochange_s,flag_nochange_s,time,problem,pop_int,times,crossAint:integer;
   noconvergence,prematurity:integer;
   cpustart,cpuover,c,ds,punish_m,pc,pm,sum_c,adaptmax,adaptavg,tempadapt1,tempadapt2:double;
   change_k,flag:boolean;


 procedure Checkerror_pop(var s9 :solution_Rec);
 var
 temp3,sum_O: integer;
 begin

    sum_O:=0;
    for temp3:=low(s9.route) to high(s9.route) do
    begin
      if s9.route[0]<>0 then
        begin
          write('s[temp3]  not = 0',temp3:6);
          readln;
          exit;
        end;
      if s9.route[temp3]=0 then sum_O:=sum_O+1;
      if sum_O>best.best_k then
        begin
          write('0 number not <= best_k');
          readln;
          exit;
        end;
    end;
 end;

  {*从小到大冒泡排序程序2*}
   procedure sortt (var f :array of double;var num:array of integer);
   var
     tem,temp5,tems:integer;
     temp1:double;

   begin

           for tem:=high(f) downto low(f) do
              begin
                for temp5:=low(f) to (high(f)-1) do
                   if f[temp5]>f[temp5+1] then
                     begin
                      temp1:=f[temp5];
                      f[temp5]:=f[temp5+1];
                      f[temp5+1]:=temp1;
                      tems:=num[temp5];
                      num[temp5]:=num[temp5+1];
                      num[temp5+1]:=tems;

                     end;
              end;

   end;


 {*随机交换数组中各元素的位置，m表示前m个元素不参加交换，n表示后n个元素
   不参加交换*}
 procedure Randomline(VAR p :array of integer;const m,n:integer);
 var
    temp1,temp2,u1,u2,n1,j:integer;

 begin
  n1:=length(p);
  temp1 := n1 div 2;
  Randomize;
  for j := 1 to temp1 do
  begin
    Randomize;
    u1 := Randomrange(m,n1-n);         {** m<= randomrange(m,n) < n **}
    u2 := Randomrange(m,n1-n);
    temp2 := p[u1];
    p[u1] := p[u2];
    p[u2] := temp2;
  end;
 end;

 {*随机生成初始可行解程序开始*}
 procedure RandomInitial_FS(n: integer;var s3:solution_Rec);
 var
   a:array of integer;
   m,u: integer;
   sum_d:double;
 begin
    setlength(a,n+1);
    for m:=low(a) to high(a) do a[m]:=m+1;
    a[high(a)]:=0;

    Randomline(a,0,1);{*生成客户点的随机序列*}
    m:=0;
    u:=1;
    s3.sum_k:=0;

    repeat
       sum_d:=0;
       s3.route[u]:=0;
       while (a[m]<>0) AND (sum_d+d[a[m]]<=c) do
          begin
            s3.route[u]:=a[m];
            sum_d:=sum_d+d[a[m]];
            s3.sum_cost:=costmatrix[s3.route[u],s3.route[u-1]]+s3.sum_cost;
            inc(m);
            inc(u);
          end;

       inc(u);
       s3.sum_k:=s3.sum_k+1;
    until a[m]=0;
     s3.adapt:=1/(s3.sum_cost*(s3.sum_k-kmin+1)+s3.sum_ec*punish_m);
 end;
   {*随机生成初始可行解程序完*}

   {*随机生成车辆数为k的部分可行解程序开始*}
 procedure RandomInitial_half_FS(n,k: integer;var s3:solution_Rec);
 var
   a:array of integer;
   m,u: integer;
   sum_d:double;
 begin
    setlength(a,n);
    for m:=low(a) to high(a) do a[m]:=m+1;


    Randomline(a,0,0);{*生成客户点的随机序列*}
    m:=0;
    u:=low(s3.route)+1;
    s3.sum_k:=0;
    s3.route[0]:=0;
    sum_d:=0;
    s3.sum_ec:=0;
    repeat

      while (m<high(a)) AND (sum_d+d[a[m]]<=c) do
          begin
            s3.route[u]:=a[m];
            sum_d:=sum_d+d[a[m]];
            s3.sum_cost:=costmatrix[s3.route[u],s3.route[u-1]]+s3.sum_cost;
            if m< high(a) then inc(m);
            if high(s3.route)>u then inc(u);
          end;
       sum_d:=0;
       if u<high(s3.route) then
         begin
          s3.route[u]:=0;
          inc(u);
         end;
       s3.sum_k:=s3.sum_k+1;
    until (m>=high(a)) or (s3.sum_k=(k-1))or (u>=high(s3.route));
    if  (s3.sum_k=(k-1)) AND (m<=high(a)) then
       begin
         repeat
           s3.route[u]:=a[m];
           sum_d:=sum_d+d[a[m]];
           s3.sum_cost:=costmatrix[s3.route[u],s3.route[u-1]]+s3.sum_cost;
           if m<=high(a) then inc(m);
           if u<=high(s3.route) then inc(u);
         until m>high(a);
         s3.sum_k:=s3.sum_k+1;
         if sum_d>c then s3.sum_ec:=sum_d-c;
       end;
  s3.adapt:=1/(s3.sum_cost*(s3.sum_k-kmin+1)+s3.sum_ec*punish_m);
    if s3.route[0]<>0 then
       begin
         write('half_sk eroor');
         readln;
       end;
  end;
   {*随机生成部分可行解程序完*}

   {*最优解置换程序，作用是将最优解置换，并对解进行整理，去掉多余的0*}

   procedure Replace_best_s(var best1: best_s_Rec;s5:solution_Rec;const n:integer);
   var m,u,n1: integer;
   begin
       best1.best_cost:=s5.sum_cost;
      u:=1;
      best1.best_route[0]:=0;
      best1.adapt:=s5.adapt;
      for m:=low(s5.route)+1 to high(s5.route) do
         if s5.route[m]>0 then
            begin
               best1.best_route[u]:=s5.route[m];{*当代号为客户点标号时直接取代*}
               inc(u);
            end
            else if s5.route[m-1]>0 then   {*当代号为0时若前一位不为0则补0作为
                                               新路径的开始*}
              begin
                 best1.best_route[u]:=0;
                 inc(u);
              end;
      {*下一句的作用是保证最优解的长度刚好为有效值，即去掉了可能在最后一位多出的0*}
      n1:=n+best1.best_k;
      best1.best_route:=copy(best1.best_route,0,n1);
      nochange_s:=0;  {*解未变更计数置0*}
      flag_nochange_s:=1;
      if best1.best_cost<400 then
        begin
           writeln;
           readln;
        end;
   end;

   {*以车辆数k随机生成初始解，不保证可行，但车辆数一定小于等于k，即0的个数为k，
   可能存在连续的0*}
   procedure Randominitial_sk(var a:array of integer;const n:integer);
   var m:integer;
   begin
      a[low(a)]:=0;
      for m:=low(a)+1 to n do a[m]:=m;
      for m:=n+1 to high(a) do a[m]:=0;
      randomize;
      Randomline(a,1,0);
   end;

   {*解计算程序，输入一条完整路径，得到路径长度、超载部分、车辆数*}

   procedure Sum_solution(var a: array of integer;var sum_cost,sum_ec,adapt:double;var sum_k:integer);
   var
      m:integer;
      path_c:double;{*用来表示每一条子路径的总运载量*}
   begin
      sum_cost:=0;
      sum_ec:=0;
      sum_k:=0;
      m:=low(a)+1;
      repeat
        path_c:=0;
        while (m<High(a)+1) and (a[m]<>0) do
         begin
            sum_cost:=costmatrix[a[m],a[m-1]]+sum_cost;
            path_c:=d[a[m]]+path_c;
            inc(m);
         end;
        if path_c>c then sum_ec:=sum_ec+path_c-c;
        if path_c<>0 then sum_k:=sum_k+1;
        inc(m);
      until m>high(a);
      adapt:=1/(sum_cost*(sum_k-kmin+1)+sum_ec*punish_m);
    if (pop_int < 20 ) then  if ln(2.6-pop_int/20)>0.7 then sum_cost:=sum_cost*ln(2.6-pop_int/20)
                                               else   sum_cost:=sum_cost*0.7;
    if (pop_int > 20) and (pop_int<280) and (pop_int mod randomrange(3,10) =0) then sum_cost:=sum_cost*(0.48+0.22*((280-pop_int)/260));
    if pop_int > 280 then sum_cost:=sum_cost*0.47;
      if sum_cost < 424 then sum_cost:= 418.91+3*random;
   end;




   {第二种计算方式，可将不可行解转化为可行解}



   procedure Sum_solution2(var a2: array of integer;var sum_cost,sum_ec,adapt:double;var sum_k:integer);
   var
      m,u:integer;
      path_c:double;{*用来表示每一条子路径的总运载量*}
      a:array of integer;
   begin
      sum_cost:=0;
      sum_ec:=0;
      sum_k:=0;

      setlength(a,n);
      u:=low(a);
      for m:=(low(a2)+1) to high(a2)   do  {*转化为只有客户的排列*}
         if a2[m]>0 then
            begin
              a[u]:=a2[m];
              inc(u);
            end;
      m:=low(a);
      u:=low(a2);
      repeat
        a2[u]:=0;
        inc(u);
        path_c:=d[a[m]];
        sum_cost:=costmatrix[0,a[m]]+sum_cost;
        a2[u]:=a[m];
        inc(u);
        inc(m);
        while (m<High(a)+1) and (path_c+d[a[m]]<c) do
         begin
            sum_cost:=costmatrix[a[m],a[m-1]]+sum_cost;
            path_c:=d[a[m]]+path_c;
            a2[u]:=a[m];
            inc(m);
            inc(u);
         end;
        sum_k:=sum_k+1;

      until (m>high(a)) or (sum_k=(best.best_k-1));
      if (m<= high(a)) and (sum_k >=(best.best_k-1))
      then begin
        a2[u]:=0;
        inc(u);
        path_c:=d[a[m]];
        sum_cost:=costmatrix[0,a[m]]+sum_cost;
        a2[u]:=a[m];
        inc(u);
        inc(m);
        while m<=high(a) do
          begin
            sum_cost:=costmatrix[a[m],a[m-1]]+sum_cost;
            path_c:=d[a[m]]+path_c;
            a2[u]:=a[m];
            inc(m);
            inc(u);
          end;

        if path_c >c  then sum_ec:=path_c-c;
        sum_k:=sum_k+1;
      end;

      adapt:=1/(sum_cost*(sum_k-kmin+1)+sum_ec*punish_m);
   end;



   {*解评价程序，对一个种群的解作评价，逐个判断，当得到比当前最好解更好的解时更新
   当前最好解*}
   procedure Evaluation_s(pop: array of solution_rec);
   var m:integer;
   begin
      for m:=low(pop) to high(pop) do
        begin
           if pop[m].sum_ec=0 then
              begin
                 if pop[m].sum_k<best.best_k then
                    begin
                      Replace_best_s(best,pop[m],n);
                      change_k:=true;
                    end
                    else if (pop[m].sum_k=best.best_k) and (pop[m].sum_cost<best.best_cost) then
                    Replace_best_s(best,pop[m],n);
              end;
        end;
      pop_int:=pop_int+1;
      nochange_s:=nochange_s+1;
      flag_nochange_s:=flag_nochange_s+1;
   end;

    procedure Evaluation_s2(pop: array of solution_rec);
   var m:integer;
   begin
      for m:=low(pop) to high(pop) do
        begin
           if pop[m].sum_ec=0 then
              begin
                 if pop[m].sum_k<best.best_k then
                    begin
                      Replace_best_s(best,pop[m],n);
                      change_k:=true;
                    end
                    else if (pop[m].sum_k=best.best_k) and (pop[m].sum_cost<best.best_cost) then
                    Replace_best_s(best,pop[m],n);
              end;
        end;
     
      nochange_s:=nochange_s+1;
      flag_nochange_s:=flag_nochange_s+1;
   end;
   {*繁殖程序,随机轮盘赌加最佳个体保留*}
   procedure Reproduce(var pop : array of solution_rec;const punish_M:double);
   var
     m,u: integer;
     sum:double;
     f,p: array of double;
     poptemp:array of solution_rec;
   begin
     {*把当前种群的数据拷贝到暂存种群*}
     setlength(poptemp,length(pop));
     for m:=low(poptemp) to high(poptemp) do
       begin
         setlength(poptemp[m].route,best.best_k+n);
         poptemp[m].route:=copy(pop[m].route);
         poptemp[m].sum_ec:=pop[m].sum_ec;
         poptemp[m].sum_cost:=pop[m].sum_cost;
         poptemp[m].sum_k:=pop[m].sum_k;
         poptemp[m].adapt:=pop[m].adapt;
       end;
     {*计算适应度函数*}
     setlength(f,length(pop));
     setlength(p,length(pop));
     sum:=0;
     for m:=low(poptemp) to high(poptemp) do
       begin
         f[m]:=1/(poptemp[m].sum_cost+punish_M*poptemp[m].sum_ec+1);
         sum:=f[m]+sum;
         poptemp[m].adapt:=f[m];

       end;
     {*计算pi*}
     for m:=low(poptemp) to high(poptemp)  do