!     TSP = Travelling Salesman Problem
!     The idea is to visit a set of N cities once and once only 
!     and return to the starting city.
!
!     This program generates all possible combinations of tours
!     with no duplicates.  Then it finds the cost of doing the tour
!     in either direction in case the cost matrix is not symmetric.
!
!     The algorithm uses a N x N matrix (actually half of it) to generate
!     and store cities not visited so far.  This is equivalent to generating
!     the complete tree of tours and then traversing the tree in order.
!
!	Programmed in FORTRAN 90 (ELF90 from Lahey - check www.lahey.com)
!     Copyright Gerald Ratzer, February 1999
!
!      Approximate timing on 366MHz Pentium II
!      11 city tour - 20 seconds - 3,628,800 tours
!      12 city tour - 250 seconds - 39,916,800 tours
!      13 city tour - 53 minutes - 479,001,600 tours
!      14 city tour - 12 hours - 6.22 billion tours
!
!      20 city tour - over 25,000 years (estimated of course!)

      program tspn
      implicit none
   
      integer :: i,j,k,m,n,ntours
      integer :: cost(20,20),city(20,20)
      integer :: visit(21),tour(21),best(21)
      integer :: cost1,cost2,mincost
      real:: notours,nobest
!
      write(*,*) 'This is Program >> TSPN  - Travelling Salesman'
      write(*,*) 'General Program for TSP tours of size N'
!     READ IN cost matrix size
      read (*,*)n
      write(*,*)'TSP tours with n=',n
      if(n < 3 .or. n > 20) then
        write(*,*)'Invalid value of n - stopping - Bye now'
        stop
      end if
!
!     READ IN cost matrix and zero city matrix
L1:   DO i=1,n
            READ (*,*) (cost(i,j),j=1,n)
            do j=1,n
            city(i,j)=0
            end do
      END DO L1
!
!     PRINT cost MATRIX
      write(*, 100)
!
L2:   DO i=1,n
            write(*, 101) (cost(I,j),j=1,n)
      END DO L2
101   FORMAT (20i4)
100   FORMAT(/' Input Cost MATRIX '/)

      ntours=0          ! Counter from 1 to 50,000
      notours=0         ! Counter for greater than 2 billion
      mincost=1E8       ! Some big number
!
!     Initialization of city data structure
!     Generate first row of city matrix
!
L3:   do i=1,n-1
         city(1,i)=i+1  ! As in 2 3 4 5 . . n  in Row 1
      end do L3
!
!      Main Calculation loop starts here
!
L4:   do while (city(1,1) /= 0)
!     Generate - expands the matrix with cities not visited so far
!
L5:      do k=2,n-1
           if(city(k,1) == 0) then
L6:          do i=1,n-1             ! Generate all cities to visit
                visit(i)=i+1        ! visit is 2 3 4 5 .. n
             end do L6
!
L61:         do i=1,n-1             !  Run down Column 1
                if(city(i,1) /= 0 ) visit(city(i,1)-1)=0
             end do L61      !  Build up list of places visited                
                             !  What's left - cities to be visited
             m=0
L7:          do i=1,n-1
                if(visit(i) /= 0) then
                   m=m+1
                   city(k,m)=visit(i)
                end if
             end do L7
             m=m+1
             city(k,m)=0   ! Mark end of this row with a 0
           end if          ! We have generated 1 new row
         end do L5         ! Repeat for rows below
!
!      Extract the tour from column 1 of city matrix
!
       tour(1)=1        ! Start from home
L8:    do i=1,n-1
          tour(i+1)=city(i,1)   ! Straight down column 1
       end do L8
       tour(n+1)=1      ! Return home 
!
       cost1=0      ! cost going
       cost2=0      ! cost coming back - in case cost matrix is not
                    ! symmetric
L9:    do i=1,n
          j=tour(i)          ! i and i+1 point to adjacent pairs 
          k=tour(i+1)        ! in the tour array
          cost1=cost1+cost(j,k)    ! j is a start city of 1 leg
          cost2=cost2+cost(k,j)    ! k is the next city
       end do L9
       ntours=ntours+1      ! Total number of tours so far
       if (cost1 < mincost) then
           mincost=cost1    ! Find minimum cost so far
           nobest=notours+ntours
           do i=1,n+1       ! Remember the best tour 
              best(i)=tour(i)
           end do
       end if  
       if (cost2 < mincost) then
           mincost=cost2
           nobest=notours+ntours
           do i=1,n+1       ! Remember the best tour 
              best(i)=tour(i)
           end do
        end if

       if(mod(ntours,50000) == 0) then
          notours=notours+50000.
          ntours=0
       write(*,99)notours,(tour(j),j=1,n+1),cost1,cost2,mincost
99     format(' Tour# ',e11.3, 20i3)
       end if 
!
!      Compress city matrix starting from the foot and move up
!
        city(n-1,1)=0    ! Zero out the last city visited
        m=n-1
L27:    do while (m > 1)
          if(city(m,1) == 0) then
!         Compress the (m-1)th row - the one above.
            i=1
            m=m-1        ! Shift previous row 1 column left
L28:        do while(city(m,i) /= 0)
               city(m,i) = city(m,i+1)
               i=i+1
            end do L28
          end if
        if(city(m,1) /= 0) exit
        end do L27     

!       Repeat the main loop over again        
        end do L4

        write(*,*)'Total number of tours = ',notours+ntours, &
         & ' and best cost is =',mincost
        Write(*,*)'The (first) best tour is number',nobest,' and is -'
        write(*,127)(best(i),i=1,n+1)
127     Format(20i4)
        stop
        end program tspn

