Este ejemplo muestra cómo utilizar la programación de enteros binarios para resolver el problema clásico del vendedor ambulante. Este problema implica encontrar el recorrido cerrado más corto (ruta) a través de un conjunto de paradas (ciudades). En este caso hay 200 paradas, pero puede cambiar fácilmente la variable para obtener un tamaño de problema diferente.nStops
Resolverás el problema inicial y verás que la solución tiene subtours. Esto significa que la solución óptima encontrada no da una ruta continua a través de todos los puntos, pero en su lugar tiene varios bucles desconectados. A continuación, utilizará un proceso iterativo para determinar los subrecorridos, agregar restricciones y volver a ejecutar la optimización hasta que se eliminen los subrecorridos.
Para el enfoque basado en problemas, vea.Problema del vendedor viajante: basado en problemas
Genere paradas aleatorias dentro de una representación poligonal bruta de los Estados Unidos continentales.
figure; load('usborder.mat','x','y','xx','yy'); rng(3,'twister') % makes a plot with 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 plot(x,y,'Color','red'); % draw the outside border hold on % Add the stops to the map plot(stopsLon,stopsLat,'*b')
Formular el problema del vendedor ambulante para la programación lineal de enteros de la siguiente manera:
Genere todos los viajes posibles, lo que significa que todos los pares de paradas diferentes.
Calcule la distancia para cada viaje.
La función de coste para minimizar es la suma de las distancias de viaje para cada viaje en el recorrido.
Las variables de decisión son binarias, y asociadas con cada viaje, donde cada 1 representa un viaje que existe en el Tour, y cada 0 representa un viaje que no está en el Tour.
Para asegurarse de que el tour incluye cada parada, incluya la restricción lineal que cada parada está en exactamente dos viajes. Esto significa una llegada y una salida de la parada.
Debido a que hay 200 paradas, hay 19.900 viajes, que significa 19.900 variables binarias (# variables = 200 elegir 2).
Generar todos los viajes, lo que significa todos los pares de paradas.
idxs = nchoosek(1:nStops,2);
Calcular todas las distancias de viaje, suponiendo que la tierra es plana con el fin de utilizar la regla 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, la duración de un recorrido esdist
dist'*x_tsp
donde está el vector de solución binaria.x_tsp
Esta es la distancia de un recorrido que intenta minimizar.
El problema tiene dos tipos de restricciones de igualdad. La primera impone que debe haber 200 viajes en total. El segundo impone que cada parada debe tener dos viajes adjuntos a él (debe haber un viaje a cada parada y un viaje saliendo de cada parada).
Especifique el primer tipo de restricción de igualdad, que debe tener viajes, en el formulario.nStops
Aeq*x_tsp = beq
Aeq = spones(1:length(idxs)); % Adds up the number of trips beq = nStops;
Para especificar el segundo tipo de restricción de igualdad, que debe haber dos viajes asociados a cada parada, extender la matriz como disperso.Aeq
Aeq = [Aeq;spalloc(nStops,length(idxs),nStops*(nStops-1))]; % allocate a sparse matrix for ii = 1:nStops whichIdxs = (idxs == ii); % find the trips that include stop ii whichIdxs = sparse(sum(whichIdxs,2)); % include trips where ii is at either end Aeq(ii+1,:) = whichIdxs'; % include in the constraint matrix end beq = [beq; 2*ones(nStops,1)];
Todas las variables de decisión son binarias. Ahora, establezca el argumento en el número de variables de decisión, coloque un límite inferior de 0 en cada uno y un límite superior de 1.intcon
intcon = 1:lendist; lb = zeros(lendist,1); ub = ones(lendist,1);
El problema está listo para la solución. Para suprimir la salida iterativa, desactive la visualización predeterminada.
opts = optimoptions('intlinprog','Display','off'); [x_tsp,costopt,exitflag,output] = intlinprog(dist,intcon,[],[],Aeq,beq,lb,ub,opts);
segments = find(x_tsp); % Get indices of lines on optimal path lh = zeros(nStops,1); % Use to store handles to lines on plot lh = updateSalesmanPlot(lh,x_tsp,idxs,stopsLon,stopsLat); title('Solution with Subtours');
Como se puede ver en el mapa, la solución tiene varios subtours. Las restricciones especificadas hasta ahora no impiden que estos subcamino ocurran. Con el fin de evitar cualquier posible subtour de suceder, se necesitaría un número increíblemente grande de restricciones de desigualdad.
Dado que no puede agregar todas las restricciones de subrecorrido, tome un enfoque iterativo. Detecte los subrecorridos en la solución actual y, a continuación, agregue restricciones de desigualdad para evitar que esos subcamino particulares ocurran. Al hacer esto, encontrará un recorrido adecuado en unas cuantas iteraciones.
Elimine subrecorridos con restricciones de desigualdad. Un ejemplo de cómo funciona esto es si tienes cinco puntos en un subtour, entonces tienes cinco líneas conectando esos puntos para crear el subtour. Eliminar este subtour implementando una restricción de desigualdad para decir que debe haber menos o igual a cuatro líneas entre estos cinco puntos.
Aún más, encuentre todas las líneas entre estos cinco puntos, y restrinja la solución para no tener más de cuatro de estas líneas presentes. Esta es una restricción correcta porque si existían cinco o más líneas en una solución, la solución tendría un subtour (un gráfico con
La función analiza la solución y devuelve una matriz de vectores de celdas.detectSubtours
Cada vector en la matriz de celdas contiene las paradas involucradas en ese subrecorrido en particular.
tours = detectSubtours(x_tsp,idxs); numtours = length(tours); % number of subtours fprintf('# of subtours: %d\n',numtours);
# of subtours: 27
Incluya las restricciones de desigualdad lineales para eliminar subrecorridos, y llame repetidamente al solucionador, hasta que solo quede un subrecorrido.
A = spalloc(0,lendist,0); % Allocate a sparse linear inequality constraint matrix b = []; while numtours > 1 % repeat until there is just one subtour % Add the subtour constraints b = [b;zeros(numtours,1)]; % allocate b A = [A;spalloc(numtours,lendist,nStops)]; % a guess at how many nonzeros to allocate for ii = 1:numtours rowIdx = size(A,1)+1; % Counter for indexing subTourIdx = tours{ii}; % Extract the current subtour % The next lines find all of the variables associated with the % particular subtour, then add an inequality constraint to prohibit % that subtour and all subtours that use those stops. variations = nchoosek(1:length(subTourIdx),2); for jj = 1:length(variations) whichVar = (sum(idxs==subTourIdx(variations(jj,1)),2)) & ... (sum(idxs==subTourIdx(variations(jj,2)),2)); A(rowIdx,whichVar) = 1; end b(rowIdx) = length(subTourIdx)-1; % One less trip than subtour stops end % Try to optimize again [x_tsp,costopt,exitflag,output] = intlinprog(dist,intcon,A,b,Aeq,beq,lb,ub,opts); % Visualize result lh = updateSalesmanPlot(lh,x_tsp,idxs,stopsLon,stopsLat); % How many subtours this time? tours = detectSubtours(x_tsp,idxs); numtours = length(tours); % 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
La solución representa un recorrido factible, ya que es un único bucle cerrado. Pero, ¿es un tour de costo mínimo? Una forma de averiguarlo es examinando la estructura de salida.
disp(output.absolutegap)
0
La pequeñez de la brecha absoluta implica que la solución es o bien óptima o tiene una longitud total que está cerca de óptima.