Disjunction does not work with priority stops

[TeodoraBogdan]

Hello!
I am trying to add disjunction for allowing the algorithm to drop stops when a solution with all of them does not exists.
I have this example with 9 stops and 1 vehicle.
If I don't set a disjunction, I receive the solution: 8, 1, 2, 7, 0, 3, 6, 5, 4, 8

Input:

Vehicles = 1
Stops = 9
Start Node = [8]
End Node = [8]
Distance Matrix:
[[0, 5846123, 4823294, 5625907, 6288410, 5926199, 5256152, 6114973, 25076],
[6189881, 0, 8327115, 9576363, 16337688, 16247708, 5352708, 3399940, 6214958],
[5383585, 8225098, 0, 15991840, 2347220, 9526869, 3370806, 10396249, 5408662],
[5371109, 8867346, 16252502, 0, 15681982, 9219884, 8014188, 3173896, 5396187],
[8166318, 17530262, 3080488, 15533248, 0, 9068277, 7454251, 19420644, 8191395],
[5762927, 16439689, 9746438, 9274403, 9175918, 0, 14122843, 13161801, 5788004],
[5748356, 5798339, 3116301, 9302443, 5529019, 15132525, 0, 7885111, 5773433],
[6408878, 3058719, 18710926, 3629279, 18593018, 12130919, 6912235, 0, 6433955],
[1420784, 5821045, 4798217, 5600830, 6263333, 5901122, 5231074, 6089896, 0]]
Time Matrix:
[[0, 971040, 757077, 912002, 1053181, 1018343, 782161, 847962, 6018],
[949704, 0, 1231799, 952045, 1280706, 1336106, 798873, 499006, 955722],
[865086, 1179297, 0, 1231678, 452079, 920324, 532147, 1341104, 871104],
[898240, 734710, 1272133, 0, 1110179, 815725, 1206672, 527548, 904258],
[1164197, 1310897, 569918, 1230873, 0, 919519, 998967, 1432325, 1170215],
[1015304, 1222643, 977476, 828801, 815522, 0, 1342327, 1030253, 1021323],
[855223, 808536, 510295, 1245850, 921039, 1358738, 0, 976943, 861241],
[987570, 459100, 1369467, 535718, 1217977, 923524, 1035151, 0, 993588],
[245912, 965021, 751059, 905984, 1047162, 1012324, 776143, 841944, 0]]
Number of Pieces = [30, 30, 20, 20, -10, -20, -20, -30, 0]      (positive numbers are pieces that need to be collected, negative numbers are pieces that need to be delivered)
Priorities = [1, 1, 1, 1, 0, 0, 0, 0, 0]             (1=HighPriority 0=LowPriority)
Max Number of Pieces per Vehicle = [150]
Start Time = [1730937600000]
End Time = [1733615940000]

The algorithm can run with multiple distance and time matrices for different vehicle types.

For the disjunction I take the biggest value in the distance matrix, if the vrp is optimized by distance, or in time matrix, if the vrp is optimized by time and multiply it by 10, to increase it.

int64 max = 0;
typedef gp::vector<gp::vector<uint32>> VRPMatrix;
gp::map<EVehicleType, VRPMatrix> matrix = m_distanceMatrices;
if (m_eOptimCriterion == EOptimizationCriterion::eTime)
	matrix = m_timeMatrices;

for (const auto& pair : matrix)	{
	VRPMatrix m = pair.second;
	for (const auto& row : m)	{
		int64 row_max = *std::max_element(row.begin(), row.end());
		if (row_max > max)
			max = row_max;
	}
}
max = max * 10;
int64 disjunctionVariable = max < INT64_MAX ? max : INT64_MAX;

for (int32 i = 0; i < manager->num_nodes(); i++) {
	int32 nodeIndex = manager->IndexToNode(i).value();
	if (std::find(startNodes.begin(), startNodes.end(), nodeIndex) == startNodes.end() &&
		std::find(endNodes.begin(), endNodes.end(), nodeIndex) == endNodes.end())
		m_routing->AddDisjunction({ nodeIndex }, disjunctionVariable);
}

This is how I set the matrices:

for (int32 vehicle = 0; vehicle < manager->num_vehicles(); vehicle++) {
	if (!m_vServiceTime.empty() || !m_timeMatrices.empty())
		timeEvaluatorIndices.push_back(m_routing->RegisterTransitCallback(
			[this, vehicle](int64 from_index, int64 to_index) -> int64 {
				auto from_node = manager->IndexToNode(from_index);
				auto to_node = manager->IndexToNode(to_index);
					int64 time = 0;
				if (!m_vServiceTime.empty())
					time = m_vServiceTime[from_node.value()];
				if (!m_timeMatrices.empty())
					time += m_timeMatrices[m_vehicleTypes[vehicle]][from_node.value()][to_node.value()];
				return time;
			}));
	if (!m_distanceMatrices.empty()) {
		distanceEvaluatorIndices.push_back(m_routing->RegisterTransitCallback(
			[this, vehicle](int64 from_index, int64 to_index) -> int64 {
				auto from_node = manager->IndexToNode(from_index);
				auto to_node = manager->IndexToNode(to_index);
					return m_distanceMatrices[m_vehicleTypes[vehicle]][from_node.value()][to_node.value()];
			}));
	}
}
m_routing->AddDimensionWithVehicleTransitAndCapacity(timeEvaluatorIndices, m_stops->GetMaxWaitTime(), m_vehicles->GetEndTimes(), false, "Time");
m_routing->AddDimensionWithVehicleTransitAndCapacity(distanceEvaluatorIndices, 0, m_vehicles->GetMaxDistForVehicles(), true, "Distance");
//set the cost function of the solution and try to minimize their values
const operations_research::RoutingDimension* dimension;
std::vector<int> evaluator;
if (m_eOptimCriterion == EOptimizationCriterion::eTime) {
	dimension = m_routing->GetMutableDimension("Time");
	evaluator = timeEvaluatorIndices;
} else {
	dimension = m_routing->GetMutableDimension("Distance");
	evaluator = distanceEvaluatorIndices;
}

for (int32 vehicle = 0; vehicle < manager->num_vehicles(); vehicle++)
	m_routing->AddVariableMinimizedByFinalizer(dimension->CumulVar(manager->IndexToNode(manager->GetEndIndex(vehicle)).value()));

This is how I set the priority to the stops:

std::vector<int64> highPriorityStops, otherPriorityStops;
gp::vector<EVRPOrderPriority> priorities = m_stops->GetPriorities();

if (std::count(priorities.begin(), priorities.end(), (EVRPOrderPriority)0) == priorities.size())
	return;
std::vector<VRPTimeLimit> timeWindows = m_stops->GetTimeWindows();
int32 stopsNumber = m_stops->GetSize() - m_vehicles->GetStartIndexes().size() - m_vehicles->GetEndIndexes().size();
operations_research::RoutingDimension* time_dimension = m_routing->GetMutableDimension("Time");

int64 maxTWHighPriority = 0;

for (int32 i = 0; i < m_stops->GetSize(); i++) {
	if (priorities[i] == EVRPOrderPriority::eHigh)	{
		highPriorityStops.push_back(i);
		if (!timeWindows.empty() && timeWindows[i].second != GP_INT64_MAX && timeWindows[i].second > maxTWHighPriority)
				maxTWHighPriority = timeWindows[i].second;
	} else
		otherPriorityStops.push_back(i);
}
if (timeWindows.empty()) {
	std::vector<int64> endTimes = m_vehicles->GetEndTimes();
	maxTWHighPriority = endTimes[0];
	if (endTimes.size() > 1)
		maxTWHighPriority = *std::min_element(endTimes.begin(), endTimes.end());
}

//if all have the same priority
int32 startAndEndsNumber = m_vehicles->GetStartIndexes().size() + m_vehicles->GetEndIndexes().size();
if (highPriorityStops.empty() || otherPriorityStops.size() == startAndEndsNumber)
	return;

std::vector<uint64> startTimes = m_vehicles->GetStartTimes();
int64 minStart = startTimes[0];
if (startTimes.size() > 1)
	minStart = *std::min_element(startTimes.begin(), startTimes.end());

for (int32 i = highPriorityStops.size()-1; i >=0; i--) {
	auto node = manager->IndexToNode(highPriorityStops[i]);
	auto index = manager->NodeToIndex(node);
	time_dimension->SetCumulVarSoftLowerBound(index, minStart, 1000);
}
for (int32 i = 0; i < otherPriorityStops.size(); i++) {
	auto node = manager->IndexToNode(otherPriorityStops[i]);
	auto index = manager->NodeToIndex(node);
	time_dimension->SetCumulVarSoftLowerBound(index, maxTWHighPriority, 1000);
}

And this is how I set the pieces:
*All the deliveries are made from the start node, that's why I created the 'deliveryPiecesDimension'

//add pieces constraint
std::vector<int32> numberOfPieces = m_stops->GetNumberOfPieces();

if (numberOfPieces.empty() || std::count(numberOfPieces.begin(), numberOfPieces.end(), 0) == numberOfPieces.size())
	return;

const int32 numberOfPiecesCallback = m_routing->RegisterTransitCallback(
	[numberOfPieces, this](int64 from_index, int64 to_index) -> int64
	{
		auto from_node = manager->IndexToNode(from_index);
		auto to_node = manager->IndexToNode(to_index);

		return numberOfPieces[from_node.value()];
	});

const int32 numberOfPiecesDeliveriesCallback = m_routing->RegisterTransitCallback(
	[numberOfPieces, this](int64 from_index, int64 to_index) -> int64
	{
		auto from_node = manager->IndexToNode(from_index);
		auto to_node = manager->IndexToNode(to_index);

		if (numberOfPieces[from_node.value()] < 0)
			return numberOfPieces[from_node.value()];
		return 0;
	});
m_routing->AddDimensionWithVehicleCapacity(numberOfPiecesCallback, 0, m_vehicles->GetMaxNumberOfPieces(), false, "Pieces");
m_routing->AddDimensionWithVehicleCapacity(numberOfPiecesDeliveriesCallback, 0, m_vehicles->GetMaxNumberOfPieces(), false, "DeliveryPieces");

const operations_research::RoutingDimension* piecesDimension = m_routing->GetMutableDimension("Pieces");
const operations_research::RoutingDimension* deliveryPiecesDimension = m_routing->GetMutableDimension("DeliveryPieces");

for (int32 i = 0; i < manager->num_nodes(); i++)
{
	int32 nodeIndex = manager->IndexToNode(i).value();
	piecesDimension->CumulVar(i)->SetMin(0);
}

operations_research::Solver* solver = m_routing->solver();

for (int32 i = 0; i < m_vehiclesSize; i++)
{
	piecesDimension->CumulVar(manager->GetStartIndex(i))->SetMin(0);
	piecesDimension->CumulVar(manager->GetEndIndex(i))->SetMin(0);

	auto startIndex = m_routing->Start(i);
	solver->AddConstraint(solver->MakeEquality(deliveryPiecesDimension->CumulVar(startIndex), piecesDimension->CumulVar(startIndex)));
}

For firstSolutionStrategies I use: LOCAL_CHEAPEST_ARC, GLOBAL_CHEAPEST_ARC, LOCAL_CHEAPEST_INSERTION, SAVINGS, PATH_CHEAPEST_ARC, PATH_MOST_CONSTRAINED_ARC, PARALLEL_CHEAPEST_INSERTION, SEQUENTIAL_CHEAPEST_INSERTION.
And for localSearchMetaheuristics I use: AUTOMATIC and GREEDY_DESCENT.
I try to find a solution by running all of the combinations between firstSolutionStrategies and localSearchMetaheuristics and I take the best one.

I tried lowering the disjunctionValue, but then it will drop stops from other examples where solution with all stops can be found.

What I want is that the algorithm to drop stops only when there is absolutely no way to visit them all.
Do you have any idea on how to update this?

[Mizux]

Did you try to enable the GLS ?

search_parameters.set_local_search_metaheuristic(
      LocalSearchMetaheuristic::GUIDED_LOCAL_SEARCH);
// GLS will run forever unless a time limit is set (here 20s)
 search_parameters.mutable_time_limit()->set_seconds(20);

Why using some near-overflow data value ?

Start Time = [1730937600000]
End Time = [1733615940000]

Can't you offset your problem to zero for start time ?

[TeodoraBogdan]

I don't want to use the GLS, because it will keep the project blocked for X seconds and to determinate the needed amount of seconds depending on the size of the input, might be tricky (for example if I have an input with 100 stops, 1000 stops or 20 stops, the amount of time needed might be different). I want to let the algorithm decide how much time it needs to give me a solution.

I tried setting the start time to 0, but I still don't receive a solution.

Do you know if there is another way to allow the algorithm to drop stops?

[TeodoraBogdan]

Or maybe do you know another way to implement the priority stops?

[TeodoraBogdan]

I changed the implementation of the SetCumulVarSoftLowerBound/SetCumulVarSoftUpperBound for time dimenstion for priority stops, and now seems to be ok.

Another issue I have now is that if I have different vehicle types with different matrices, after I set the disjunction, the stops order in the solution is different than how it was before. The new solution is worse, having a total duration greater than the one before. Why is that and how can I fix it?

[jmarca]

First, in reference to your earlier comment, GLS should be considered a requirement for all but the most trivial of problems. Without it the solver will most likely stop at a local minimum. Perhaps that relates to this problem you raise now, with the solution apparently being worse. On the other hand, you haven't really shared how you implemented the soft lower/upper bounds, and the different per-vehicle matrices, so it is impossible to say what might be going wrong here. I for one am not that good of a debugger! James James E. Marca Activimetrics LLC

[TeodoraBogdan]

"First, in reference to your earlier comment, GLS should be considered a requirement for all but the most trivial of problems. Without it the solver will most likely stop at a local minimum." GLS should be always used or is absolutely needed in this type of problem (with disjunction)?
If so, what is a good approach to set a time limit corelated to the size of the problem? For example an input with 100 stops might need more time than an input with 20 stops or maybe the input has more constraints, such as time window, priorities, capacities, etc. How would I know how much time to set?

"On the other hand, you haven't really shared how you implemented the soft lower/upper bounds, and the different per-vehicle matrices, so it is impossible to say what might be going wrong here." I showed how I set the matrices. Do you want to see how i set different constraints? If so, which ones?