Main Content

Problema del viajante: Basado en problemas

Este ejemplo muestra cómo utilizar programación de enteros binaria para resolver el problema clásico del viajante. Este problema implica encontrar el trayecto (ruta) cerrado más corto pasando por un conjunto de paradas (ciudades). En este caso hay 200 paradas, pero puede cambiar fácilmente la variable nStops para obtener un tamaño de problema diferente. Resolverá el problema inicial y verá que la solución tiene subtrayectos. Esto significa que la solución óptima encontrada no proporciona una ruta continua por todos los puntos, sino que tiene varios bucles desconectados. Después utilizará un proceso iterativo para determinar los subtrayectos, añadiendo restricciones y volviendo a ejecutar la optimización hasta que los subtrayectos se eliminen.

Si desea ver el enfoque basado en solvers para este problema, consulte Problema del viajante: Basado en solvers.

Formulación del problema

Formule el problema del viajante para programación lineal de enteros de la siguiente forma:

  • Genere todos los viajes posibles, es decir, todas las parejas de paradas diferentes.

  • Calcule la distancia para cada viaje.

  • La función de coste que se desea minimizar es la suma de las distancias del viaje para cada viaje del trayecto.

  • Las variables de decisión son binarias y están asociadas a cada viaje, donde cada 1 representa un viaje que existe en el trayecto y cada 0 representa un viaje que no está en el trayecto.

  • Para garantizar que el trayecto incluya todas las paradas, añada la restricción lineal de que cada parada esté exactamente entre dos viajes. Esto implica una llegada y una salida en cada parada.

Generar paradas

Genere paradas aleatorias dentro de una representación poligonal en bruto de los EE. UU. continentales.

load('usborder.mat','x','y','xx','yy');
rng(3,'twister') % Makes stops in Maine & Florida, and is reproducible
nStops = 200; % You can use any number, but the problem size scales as N^2
stopsLon = zeros(nStops,1); % Allocate x-coordinates of nStops
stopsLat = stopsLon; % Allocate y-coordinates
n = 1;
while (n <= nStops)
    xp = rand*1.5;
    yp = rand;
    if inpolygon(xp,yp,x,y) % Test if inside the border
        stopsLon(n) = xp;
        stopsLat(n) = yp;
        n = n+1;
    end
end

Calcular distancias entre puntos

Dado que existen 200 paradas, hay 19.900 viajes, es decir, 19.900 variables binarias (# variables = 200 elegir 2).

Genere todos los viajes, es decir, todas las parejas de paradas.

idxs = nchoosek(1:nStops,2);

Calcule todas las distancias de los viajes, asumiendo que la tierra es plana para utilizar el teorema de Pitágoras.

dist = hypot(stopsLat(idxs(:,1)) - stopsLat(idxs(:,2)), ...
             stopsLon(idxs(:,1)) - stopsLon(idxs(:,2)));
lendist = length(dist);

Con esta definición del vector dist, la longitud de un trayecto es

dist'*trips

donde trips es el vector binario que representa los viajes que realiza la solución. Esta es la distancia de un trayecto que se intenta minimizar.

Crear una gráfica y diseñar un mapa

Represente el problema como una gráfica. Cree una gráfica donde todas las paradas sean nodos y los viajes sean bordes.

G = graph(idxs(:,1),idxs(:,2));

Muestre las paradas utilizando una gráfica. Represente los nodos sin los bordes de la gráfica.

figure
hGraph = plot(G,'XData',stopsLon,'YData',stopsLat,'LineStyle','none','NodeLabel',{});
hold on
% Draw the outside border
plot(x,y,'r-')
hold off

Crear variables y problemas

Cree un problema de optimización con variables de optimización binarias que representan los viajes potenciales.

tsp = optimproblem;
trips = optimvar('trips',lendist,1,'Type','integer','LowerBound',0,'UpperBound',1);

Incluya la función objetivo en el problema.

tsp.Objective = dist'*trips;

Restricciones

Cree las restricciones lineales de que cada parada tenga dos viajes asociados, ya que debe haber un viaje hasta cada parada y un viaje desde cada parada.

Utilice la representación de la gráfica para identificar todos los viajes que comienzan y terminan en una parada encontrando todos los bordes que están conectados a esa parada. Para cada parada, cree la restricción de que la suma de los viajes para esa parada es igual a dos.

constr2trips = optimconstr(nStops,1);
for stop = 1:nStops
    whichIdxs = outedges(G,stop); % Identify trips associated with the stop
    constr2trips(stop) = sum(trips(whichIdxs)) == 2;
end
tsp.Constraints.constr2trips = constr2trips;

Resolver el problema inicial

El problema está listo para resolverlo. Para suprimir la salida iterativa, desactive la visualización predeterminada.

opts = optimoptions('intlinprog','Display','off');
tspsol = solve(tsp,'options',opts)
tspsol = struct with fields:
    trips: [19900×1 double]

Visualizar la solución

Cree una nueva gráfica con los viajes de la solución como bordes. Para hacerlo, redondee la solución en caso de que algunos valores no sean exactamente enteros, y convierta los valores resultantes a logical.

tspsol.trips = logical(round(tspsol.trips));
Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2),[],numnodes(G));
% Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2)); % Also works in most cases

Superponga la nueva gráfica sobre la gráfica existente y destaque sus bordes.

hold on
highlight(hGraph,Gsol,'LineStyle','-')
title('Solution with Subtours')

Como puede verse en el mapa, la solución tiene varios subtrayectos. Las restricciones especificadas hasta ahora no evitan estos subtrayectos. Para evitar cualquier posible subtrayecto, necesitaría un número increíblemente amplio de restricciones de desigualdad.

Restricciones de subtrayecto

Dado que no se pueden añadir todas las restricciones de subtrayecto, adopte un enfoque iterativo. Detecte los subtrayectos en la solución actual y, después, añada restricciones de desigualdad para evitar esos subtrayectos en particular. Al hacerlo, se encuentra un trayecto adecuado en unas pocas iteraciones.

Elimine los subtrayectos con restricciones de desigualdad. Un ejemplo de cómo funciona es que, si tiene cinco puntos en un subtrayecto, tendrá cinco líneas conectando estos puntos para crear el subtrayecto. Elimine este subtrayecto implementando una restricción de desigualdad que diga que debe haber cuatro líneas o menos entre estos cinco puntos.

Aún mejor, encuentre todas las líneas entre estos cinco puntos y restrinja la solución para que no tenga más de cuatro de estas líneas. Se trata de una restricción correcta porque, si en una solución existen cinco o más líneas, la solución tendría un subtrayecto (una gráfica con n nodos y n bordes siempre contiene un ciclo).

Detecte el subtrayecto identificando los componentes conectados en Gsol, la gráfica creada con los bordes en la solución actual. conncomp devuelve un vector con el número del subtrayecto al que pertenece cada borde.

tourIdxs = conncomp(Gsol);
numtours = max(tourIdxs); % Number of subtours
fprintf('# of subtours: %d\n',numtours);
# of subtours: 27

Incluya las restricciones de desigualdad lineales para eliminar subtrayectos y llame repetidamente al solver hasta que solo quede un subtrayecto.

% Index of added constraints for subtours
k = 1;
while numtours > 1 % Repeat until there is just one subtour
    % Add the subtour constraints
    for ii = 1:numtours
        inSubTour = (tourIdxs == ii); % Edges in current subtour
        a = all(inSubTour(idxs),2); % Complete graph indices with both ends in subtour
        constrname = "subtourconstr" + num2str(k);
        tsp.Constraints.(constrname) = sum(trips(a)) <= (nnz(inSubTour) - 1);
        k = k + 1;        
    end
    
    % Try to optimize again
    [tspsol,fval,exitflag,output] = solve(tsp,'options',opts);
    tspsol.trips = logical(round(tspsol.trips));
    Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2),[],numnodes(G));
    % Gsol = graph(idxs(tspsol.trips,1),idxs(tspsol.trips,2)); % Also works in most cases
    
    % Plot new solution
    hGraph.LineStyle = 'none'; % Remove the previous highlighted path
    highlight(hGraph,Gsol,'LineStyle','-')
    drawnow

    % How many subtours this time?
    tourIdxs = conncomp(Gsol);
    numtours = max(tourIdxs); % Number of subtours
    fprintf('# of subtours: %d\n',numtours)    
end
# of subtours: 20
# of subtours: 7
# of subtours: 9
# of subtours: 9
# of subtours: 3
# of subtours: 2
# of subtours: 7
# of subtours: 2
# of subtours: 1
title('Solution with Subtours Eliminated');
hold off

Calidad de la solución

La solución representa un trayecto factible porque se trata de un bucle cerrado único. Pero, ¿se trata de un trayecto con coste mínimo? Una forma de descubrirlo es examinar la estructura de salida.

disp(output.absolutegap)
   1.8300e-04

El tamaño pequeño del intervalo absoluto implica que la solución es óptima o tiene una longitud total casi óptima.

Temas relacionados