Traveling Salesman Prob with Local Optimization(ga1)

/****************************************************************/
/*                                                              */
/*             S A S   S A M P L E   L I B R A R Y              */
/*                                                              */
/*    NAME: ga1                                                 */
/*   TITLE: Traveling Salesman Prob with Local Optimization(ga1)*/
/* PRODUCT: OR                                                  */
/*  SYSTEM: ALL                                                 */
/*    KEYS: OR                                                  */
/*   PROCS: GA                                                  */
/*    DATA:                                                     */
/*                                                              */
/* SUPPORT:                             UPDATE:                 */
/*     REF:                                                     */
/*    MISC: Example 1 from the GA chapter of Local Search       */
/*          Optimization.                                       */
/*                                                              */
/****************************************************************/

/* 20 random locations for a  Traveling Salesman Problem */
data locations;
   input x y;
   datalines;
0.0333692 0.9925079
0.6020896 0.0168807
0.1532083 0.7020444
0.3181124 0.1469288
0.1878440 0.8679120
0.9786112 0.4925364
0.7918010 0.7943144
0.5145329 0.0363478
0.5500754 0.8324617
0.3893757 0.6635483
0.9641841 0.6400201
0.7718126 0.5463923
0.7549037 0.4584584
0.2837881 0.7733415
0.3308411 0.1974851
0.7977221 0.1193149
0.3221207 0.7930478
0.9201035 0.1186234
0.2397964 0.1448552
0.3967470 0.6716172
;


proc ga data1 = locations seed = 5554;
   call SetEncoding('S20');
   ncities = 20;
   array distances[20,20] /nosym;
   do i = 1 to 20;
      do j = 1 to i;
         distances[i,j] = sqrt((x[i] - x[j])**2 + (y[i] - y[j])**2);
         distances[j,i] = distances[i,j];
      end;
   end;
   call SetObj('TSP',0,'distances', distances);
   call SetCross('Order');
   call SetMut('Invert');
   call SetMutProb(0.05);
   call SetCrossProb(0.8);
   call SetElite(1);
   call Initialize('DEFAULT',200);
   call ContinueFor(140);
run;

proc ga data1 = locations seed = 5554;
   call SetEncoding('S20');
   ncities = 20;
   array distances[20,20] /nosym;
   do i = 1 to 20;
      do j = 1 to i;
         distances[i,j] = sqrt((x[i] - x[j])**2 + (y[i] - y[j])**2);
         distances[j,i] = distances[i,j];
      end;
   end;
   /* Objective function with local optimization */
   function TSPSwap(selected[*],ncities,distances[*,*]);
   array s[1] /nosym;
   call dynamic_array(s,ncities);
   call ReadMember(selected,1,s);
   /* First try to improve solution by swapping adjacent cities */
   do i = 1 to ncities;
      city1 = s[i];
      inext = 1 + mod(i,ncities);
      city2 = s[inext];
      if i=1 then
         before = s[ncities];
      else
         before = s[i-1];
         after = s[1 + mod(inext,ncities)];
      if (distances[before,city1]+distances[city2,after]) >
         (distances[before,city2]+distances[city1,after]) then do;
         s[i] = city2;
         s[inext] = city1;
      end;
   end;
   call WriteMember(selected,1,s);
   /* Now compute distance of tour */
   distance = distances[s[ncities],s[1]];
   do i = 1 to (ncities - 1);
      distance + distances[s[i],s[i+1]];
   end;
   return(distance);
   endsub;
   call SetObjFunc('TSPSwap',0);
   call SetCross('Order');
   call SetMut('Invert');
   call SetMutProb(0.05);
   call SetCrossProb(0.8);
   call SetElite(1);
   call Initialize('DEFAULT',200);
   call ContinueFor(35);
run;