function [ order_best, cost_best ] = tsp_anneal ( x, temp, rate )

%*****************************************************************************80
%
%% tsp_anneal() solves the traveling salesperson problem (TSP) using simulated annealing.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    14 June 2026
%
%  Author:
%
%    Original R code by James Howard;
%    This version by John Burkardt.
%
%  Reference:
%
%    James Howard,
%    Computational Methods for Numerical Analysis with R,
%    CRC Press, 2017,
%    ISBN13: 978-1-4987-2363-3.
%
%  Input:
%
%    real x(n,*): the coordinates of the cities.
%    A second dimension of 2 is typical, but a 3D problem would also work.
%
%    real temp: the annealing temperature.
%
%    real rate: the annealing rate.
%    A default value is 0.0001
%
%  Output:
%
%    integer order_best(n): the optimal ordering.
%
%    real cost_best: the length of the optimal round trip.
%
  trial_num = 5;

  step = 1.0 - rate;
%
%  N is the number of cities.
%
  n = size ( x, 1 );
%
%  Create n-vectors order_best, order_next, order_curr.
%
  order_best = randperm ( n );
  order_next = order_best;
  order_curr = order_best;
%
%  Compute the lengths cost_best, cost_next, cost_curr.
%
  cost_best = tsp_tour_cost ( n, x, order_best );
  cost_next = cost_best;
  cost_curr = cost_best;
%
%  As the temperature decreases, repeatedly adjust the current path.
%
  while ( 1.0 < temp )
%
%  Try up to NTRY variations.
%
    for trial = 1 : trial_num
%
%  Pick a path component i to swap with component i-1.
%
      i = randi ( n );
      im1 = i - 1;
      if ( im1 == 0 )
        im1 = n;
      end
%
%  Create the next path candidate, by swapping positions i and i-1.
%
      order_next = order_curr;
%
%  Swap positions of cities i and i-1.
%
      temporary = order_next(i);
      order_next(i) = order_next(im1);
      order_next(im1) = temporary;
%
%  Compute the path length.
%
      cost_next = tsp_tour_cost ( n, x, order_next );
%
%  The candidate solution is accepted if
%  A) the cost is lower than the best cost seen so for, OR
%  B) the cost is at most slightly greater than the best cost,
%     where "slightly" depends on the temperature.
%
      accept = exp ( - ( cost_next - cost_curr ) / temp );

      if ( cost_next < cost_curr || rand ( ) < accept )
        order_curr = order_next;
        cost_curr = cost_next;
      end

      if ( cost_next < cost_best )
        order_best = order_curr;
        cost_best = cost_curr;
      end

    end
%
%  Decrease the temperature for the next round.
%
    temp = temp * step;

  end

  return
end
function cost = tsp_tour_cost ( n, x, order )

%*****************************************************************************80
%
%% tsp_tour_cost() evaluates the cost of a TSP tour.
%
%  Licensing:
%
%    This code is distributed under the MIT license.
%
%  Modified:
%
%    14 June 2026
%
%  Author:
%
%    John Burkardt
%
%  Input:
%
%    integer n: the number of cities.
%
%    real x(n,:): the city coordinates.
%
%    integer order(n): the order in which the cities are visited.
%
%  Output:
%
%    real cost: the cost of the route.
%
  cost = 0.0;
  i1 = n;
  for i2 = 1 : n
    cost = cost + norm ( x(order(i1),:) - x(order(i2),:) );
    i1 = i2;
  end

  return
end
