Lindo Systems

MODEL:
 !  Traveling Salesman Problem. Simple version.(TSPCUTx)
   Find the shortest tour that visits each city exactly once. 
   Subtour elimination method. We add a row/cut/constraint
   to cut off a subtour until a full tour is found;
 ! Keywords: TSP, Traveling Salesman, Cut Generation, Tour, Subtour Elimination;
 SETS:
  CITY; 
  LINK( CITY, CITY):
       DIST,  ! The distance matrix;
          Y;  ! Y( I, J) = 1 if link I, J is in tour;
  SUBTOUR: TOURSIZE;
  SXC(SUBTOUR,CITY): FLAG;
 ENDSETS
DATA:
  SUBTOUR = 1..30;  ! Max subtour cuts allowed;
  CITY = 
   NYK  ATL  CHI  CIN  HOU  LAX  MON  PHL  PIT  STL  SND  SNF;
   ! Distance matrix, it need not be symmetric;
  DIST =
     0  864  845  664 1706 2844  396   92  386 1002 2892 3032 !NYK;
   864    0  702  454  842 2396 1196  772  714  554 2363 2679 !ATL;
   845  702    0  324 1093 2136  764  764  459  294 2184 2187 !CHI;
   664  454  324    0 1137 2180  798  572  284  338 2228 2463 !CIN;
  1706  842 1093 1137    0 1616 1857 1614 1421  799 1521 2021 !HOU;
  2844 2396 2136 2180 1616    0 2900 2752 2464 1842   95  405 !LAX;
   396 1196  764  798 1857 2900    0  424  514 1058 2948 2951 !MON;
    92  772  764  572 1614 2752  424    0  305  910 2800 2951 !PHL;
   386  714  459  284 1421 2464  514  305    0  622 2512 2646 !PIT;
  1002  554  294  338  799 1842 1058  910  622    0 1890 2125 !STL;
  2892 2363 2184 2228 1521   95 2948 2800 2512 1890    0  500 !SND;
  3032 2679 2187 2463 2021  405 2951 2951 2646 2125  500    0;!SNF;
ENDDATA
!-----------------------------------------------------------;
SUBMODEL TSP_CUT:
! Minimize total distance traveled;
 MIN = TOURLEN;
    TOURLEN = @SUM( LINK: DIST * Y);

! The Assignment constraints;
  @FOR( CITY( K):
  !  It must be entered;
    @SUM( CITY( I)| I #NE# K: Y( I, K)) = 1;
  !  It must be departed;
    @SUM( CITY( J)| J #NE# K: Y( K, J)) = 1;
  !  Cannot visit self;        Y( K, K) = 0;
   );

! Subtour cuts added so far;
  @FOR( SUBTOUR(t) | t #LE# ICUT:
     @SUM( CITY(I) | FLAG(t,i) #EQ# 1:
      @SUM( CITY(J) | FLAG(t,j) #EQ# 1: Y(i,j))) <= TOURSIZE(t) - 1;
      );

 ! The Y(i,j) must be 0 or 1;
  @FOR( LINK(i,j): @BIN( Y(i,j)));
ENDSUBMODEL

CALC:
  @SET("TERSEO",2);  ! Turn off default output;
  N = @SIZE( CITY);
  MXCUTS = @SIZE(SUBTOUR); ! Max cuts we have space for;
  ICUT = 0;
  MORECUTS = 1; ! =0 if no more cuts;
  @WHILE( MORECUTS :
! Loop over subtour cuts, ICUT;
!  Solve current version;
   @SOLVE( TSP_CUT);
   ICUT = ICUT + 1; 
   @WRITE( @NEWLINE(1),'#',ICUT,' Obj value = ', TOURLEN, @NEWLINE(1));

   ! Turn this on to see intermediate solutions;
   !@FOR( LINK(i,j) | Y(i,j) #GT# .01: @WRITE(' ',i,' ',j,' ',Y(i,j),@NEWLINE(1)) );

!  Find subtour if any; 
   StartCity = 1; ! Start at city 1;
   KURSTOP = StartCity;
   @FOR( CITY(k): FLAG(ICUT,k) = 0); ! Start with no cities in cut;
   @WRITE('  Begin tour  at     ', CITY(KURSTOP), @NEWLINE(1));
   TOURSIZE(ICUT) = 1;
!  Loop over cities KURSTOP to find subtour including 1;		
   FLAG(ICUT, KURSTOP) = 1; ! City KURSTOP is in the cut;
   NotHome = 1;
   @WHILE( NotHome:	
     @FOR(CITY(j) | Y(KURSTOP, j) #GT# 0.5:  ! Find next stop;
      @WRITE('  Next stop on tour= ', CITY(j),@NEWLINE(1));
	  @IFC( j #EQ# startCity:
             NotHome = 0;   ! Back home to startCity;
          @ELSE         
             TOURSIZE(ICUT) = TOURSIZE(ICUT) + 1;
             KURSTOP = j;                         
             FLAG(ICUT,KURSTOP) = 1; ! KURSTOP is in the cut;
             ); ! end of IFC;
	 ); ! end of for loop;
      ); ! end of while loop;

    @IFC( TOURSIZE(ICUT) #EQ# N:
      @WRITE(' Above is complete tour, so min length tour found.');
      MORECUTS = 0;
     @ELSE
      @WRITE(' Add constraint to cut off above subtour.',@NEWLINE(1));
        );
    @IFC( ICUT #EQ# @SIZE(SUBTOUR):
      @WRITE(' Exhausted memory.');
      MORECUTS = 0;
        );
    ); ! End loop over add cuts;
!   (There are various refinements/complications of this basic idea 
     to make order of magnitude performance improvements);
       
ENDCALC