# WgraphTSPfile.pl 
# copyright 2018 Dan Gusfield

# call on a command line in a terminal window with:
#
# perl WgraphTSPfile.pl edge-weight-file max-weight-for-a-real-edge

# Program WgraphTSPfile.pl creates the MTZ ILP formulation to compute a TS *path* from an edge-weight file
# for a graph G,  with no specific start node.
# It adds a 0 node with edges between 0 and all other nodes 1 through n, with cost in the real data. 
# 0. That transformation creates graph G' from G. So, the ILP computes a TS tour in G'.

# Also, this version of the program handles asymmetric data - as long as for every pair of nodes (i,j)
# both (i,j) and (j,i) are allowed, i.e., with cost below threshold, or both are disallowed.
#
#

#
open (IN, "$ARGV[0]");
$file = $ARGV[0];
open (TSP, ">W$file.lp");
$bound = $ARGV[1];

@D = ();
%bin = ();
%generals = ();
$X0 = "";
$inflow = "";
$outflow = "";
$conn_constraints = "";
$objective = "";
$Xconstraints = "";

$i = 1;
while ($line = <IN>) {
    @linedata = split (/ /, $line);
    $ncols = @linedata;
       foreach $j (1 .. $ncols) {
             $D[$i][$j] = $linedata[$j-1];  # works for both symmetric and assymetric data
          }
    $i++;
}
$n = $i -1; # the number of rows in the input matrix for the graph.


foreach $p (0..$n) {
    $generals{"U$p"} = 1;
    foreach $q (0..$n) {
        
       if ($p != $q) {    
         $Xconstraints .= "X$p,$q + X$q,$p <= 1\n";
       	 $outflow .= "+ X$p,$q ";
	 $inflow .= "+ X$q,$p ";
	 $bin{"X$p,$q"} = 1;
	 $bin{"X$q,$p"} = 1;
         if (($D[$p][$q] < $bound ) and ($p > 0) and ($q > 0)) {   # check if distance is for a real edge.
#             print "Distance for $p, $q, $D[$p][$q] \n";

            $mult100 = int ($D[$p][$q] * 100); # the effect of this is to round to two decimal places after the .
            $D[$p][$q] = $mult100/100;
            $mult100 = int ($D[$q][$p] * 100); # the effect of this is to round to two decimal places after the .
            $D[$q][$p] = $mult100/100;
#
            $objective .= "+ $D[$p][$q] X$p,$q ";
          }
          if (($D[$p][$q] >= $bound) and ($p > 0) and ($q > 0)) {   
#              print "Excluded edge $p, $q, $D[$p][$q]\n";
              $X0 .= "X$p,$q = 0\n";
          }
          $nm1 = $n - 1;
          if ($q > 0) {
              $conn_constraints .= "U$p - U$q + $n X$p,$q <= $nm1\n";
          }
#          if (($p > 0) and ($q > 0)) {
#               $conn_constraints .= "U$p - U$q + $n X$p,$q <= $nm1\n";
#          }
      }
    }
    $outflow .= "= 1\n";
    $inflow .= "= 1\n";
}



print TSP "minimize\n";
print TSP "$objective\n";
print TSP "such that\n";
print TSP "$inflow";
print TSP "$outflow";
print TSP "$conn_constraints";
print TSP "$Xconstraints";
print TSP "$X0";
print TSP "binaries\n";
foreach $var (sort keys %bin) {
  print TSP "$var\n";
}
print TSP "generals\n";
foreach $var (sort keys %generals) {
  print TSP "$var\n";
}

print TSP "end";
close (TSP);
print "\n The ILP formulation is in file W$file.lp \n";