Location Routing Problem (LRP)¶
Principles learned¶
Combine list and set decision variables
Use the find operator
Use a lambda expression to define a recursive array
Problem¶
In the Location Routing Problem (LRP), which is an extension of the Capacitated Vehicle Routing Problem (CVRP), a fleet of delivery vehicles with uniform capacity must service customers with known demand for a single commodity. Contrary to the CVRP, which treats the LRP problem with only one depot location, different depot locations are available, each of them having its own opening cost and capacity. The vehicles end their routes at their starting depot. Each customer can only be served by one vehicle. The objectives are to determine which depots to open and routes to take while minimizing the total cost, such that all customers are served, and the total demand served by each depot and truck does not exceed its capacity.
Download the exampleData¶
The instances used in this example come from the S. Barreto instances.
The format of the data files is as follows:
The number of customers
The number of depots
The x and y coordinates of the depots and the customers
The capacity of the delivery vehicles
The capacity of each depot
The demand for every customer
The opening cost of each depot
The opening cost of a route
The areCostDouble boolean value indicating if costs should be rounded
Program¶
This LocalSolver model defines a sequence for nbTrucks
trucks as a list
variable (customersSequences[r]
). It corresponds to the sequence of
customers visited. To ensure that all customers must be served, all the list
variables are constrained to form a partition by using to the “partition”
operator. The model also defines a group of sequences
associated to every depot as a set variable (depots[d]
). The depots are
also constrained to form a partition, to ensure that every sequence is
associated to a depot.
The number of sequences nbTrucks
is determined as the minimum number of
trucks necessary to serve every customer without overloading the vehicles,
multiplied by a factor 1.5 to ensure feasibility. The quantity delivered to a
customer has to meet the corresponding demand. The quantity served in a
sequence quantityServed[r]
is defined with a sum and the access to the
demands
array with the “at” operator
and must not exceed its capacity.
To determine the depot assigned to a route, we use the “find” operator.
A depot is considered opened if at least one (non-empty) sequence uses it.
This condition is computed with the “count” operator. The depotsCost
are
defined as the sum of the opening costs of the opened depots. Similarly, a
sequence is considered opened if at least one customer is served for this
sequence. We calculate the total distance for every sequence by adding the
distance from the depot to the first customer, the distances between every
customer and the distance back to the depot from the last customer. The sum of
sequenceLength
and openingCostRoute
corresponds to the routingCost
.
The objective is to minimize the sum of depotsCost
and routingCost
.
- Execution:
- localsolver location_routing_problem.lsp inFileName=instances/coordChrist100.dat [lsTimeLimit=] [solFileName=]
use io;
/* Read instance data */
function input() {
local usage = "Usage : localsolver location_routing_problem.lsp "
+ "inFileName=inputFile [solFileName=outputFile] [lsTimeLimit=timeLimit]";
if (inFileName == nil) throw usage;
readInputLrp();
computeDistances();
minNbTrucks = ceil(sum[c in 0...nbCustomers](demands[c]) / truckCapacity);
nbTrucks = ceil(1.5 * minNbTrucks);
}
/* Declare the optimization model */
function model() {
// A route is represented as a list containing the customers in the order they are
// visited
customersSequences[r in 0...nbTrucks] <- list(nbCustomers);
// A depot is represented as a set containing the associated customersSequences
depots[d in 0...nbDepots] <- set(nbTrucks);
// All customers should be assigned to a route
constraint partition(customersSequences);
// All the customersSequences should be assigned to a depot
constraint partition(depots);
for [r in 0...nbTrucks] {
local sequence <- customersSequences[r];
local c <- count(sequence);
quantityServed[r] <- sum(sequence, j => demands[j]);
// The quantity needed in each route must not exceed the vehicle capacity
constraint quantityServed[r] <= truckCapacity;
// A route is used if it serves at least one customer
sequenceUsed[r] <- c > 0;
// The "find" function gets the depot that is assigned to the route
associatedDepot[r] <- find(depots, r);
sequenceLength[r] <- sum(0...c-1, i => distCustomers[sequence[i]][sequence[i + 1]])
+ (sequenceUsed[r] ? distCustomersDepots[sequence[0]][associatedDepot[r]]
+ distCustomersDepots[sequence[c - 1]][associatedDepot[r]] : 0);
// The route cost is the sum of the opening cost and the route length
routeCost[r] <- sequenceLength[r] + sequenceUsed[r] * openingRouteCost;
}
for [d in 0...nbDepots] {
// The total demand served by a depot must not exceed its capacity
constraint sum(depots[d], i => quantityServed[i]) <= depotsCapacity[d];
// A depot is open if at least a route starts from there
isDepotOpen[d] <- count(depots[d]) > 0;
}
depotsCost <- sum[d in 0...nbDepots](isDepotOpen[d] * openingCostDepots[d]);
routingCost <- sum[r in 0...nbTrucks](routeCost[r]);
totalCost <- routingCost + depotsCost;
minimize totalCost;
}
function param() {
if (lsTimeLimit == nil) lsTimeLimit = 20;
}
/* Write the solution in a file */
function output() {
if (solFileName == nil) return;
local resFile = io.openWrite(solFileName);
resFile.println("File name: ", inFileName + "; totalCost = " + totalCost.value);
for [r in 0...nbTrucks] {
if (sequenceUsed[r].value) {
resFile.print("Route ", r, ", assigned to depot ", associatedDepot[r].value, ": ");
for [customer in customersSequences[r].value] {
resFile.print(customer, " ");
}
resFile.println();
}
}
resFile.close();
}
function readInputLrp() {
if (inFileName.endsWith(".dat")) readInputDat();
else throw "Unknown file format";
}
function readInputDat() {
local inFile = io.openRead(inFileName);
nbCustomers = inFile.readInt();
nbDepots = inFile.readInt();
coordinatesDepots[d in 0...nbDepots][l in 0...2] = inFile.readDouble();
coordinatesCustomers[c in 0...nbCustomers][l in 0...2] = inFile.readDouble();
truckCapacity = inFile.readInt();
depotsCapacity[d in 0...nbDepots] = inFile.readInt();
demands[c in 0...nbCustomers] = inFile.readInt();
local tempOpeningCostDepots[d in 0...nbDepots] = inFile.readDouble();
local tempOpeningCostRoute = inFile.readDouble();
areCostsDouble = inFile.readInt();
if (areCostsDouble == 1) {
openingCostDepots[d in 0...nbDepots] = tempOpeningCostDepots[d];
openingRouteCost = tempOpeningCostRoute;
} else {
openingCostDepots[d in 0...nbDepots] = round(tempOpeningCostDepots[d]);
openingRouteCost = round(tempOpeningCostRoute);
}
}
function computeDistances() {
local tempDistanceCustomers[c0 in 0...nbCustomers][c1 in 0...nbCustomers] =
sqrt(pow((coordinatesCustomers[c0][0] - coordinatesCustomers[c1][0]), 2)
+ pow((coordinatesCustomers[c0][1] - coordinatesCustomers[c1][1]), 2));
local tempDistanceCustomersDepots[c in 0...nbCustomers][d in 0...nbDepots] =
sqrt(pow((coordinatesCustomers[c][0] - coordinatesDepots[d][0]), 2)
+ pow((coordinatesCustomers[c][1] - coordinatesDepots[d][1]), 2));
if (areCostsDouble == 1) {
distCustomers[c0 in 0...nbCustomers][c1 in 0...nbCustomers] =
tempDistanceCustomers[c0][c1];
distCustomersDepots[c in 0...nbCustomers][d in 0...nbDepots] =
tempDistanceCustomersDepots[c][d];
} else {
distCustomers[c0 in 0...nbCustomers][c1 in 0...nbCustomers] =
ceil(100 * tempDistanceCustomers[c0][c1]);
distCustomersDepots[c in 0...nbCustomers][d in 0...nbDepots] =
ceil(100 * tempDistanceCustomersDepots[c][d]);
}
}
- Execution (Windows)
- set PYTHONPATH=%LS_HOME%\bin\pythonpython location_routing_problem.py instances/coordChrist100.dat
- Execution (Linux)
- export PYTHONPATH=/opt/localsolver_12_5/bin/pythonpython location_routing_problem.py instances/coordChrist100.dat
import localsolver
import sys
import math
def read_elem(filename):
with open(filename) as f:
return [str(elem) for elem in f.read().split()]
def main(instance_file, str_time_limit, sol_file):
#
# Read instance data
#
nb_customers, nb_depots, vehicle_capacity, opening_route_cost, demands_data, \
capacity_depots, opening_depots_cost, dist_matrix_data, dist_depots_data = \
read_input_lrp(instance_file)
min_nb_trucks = int(math.ceil(sum(demands_data) / vehicle_capacity))
nb_trucks = int(math.ceil(1.5 * min_nb_trucks))
with localsolver.LocalSolver() as ls:
#
# Declare the optimization model
#
m = ls.model
# A route is represented as a list containing the customers in the order they are
# visited
customers_sequences = [m.list(nb_customers) for _ in range(nb_trucks)]
# All customers should be assigned to a route
m.constraint(m.partition(customers_sequences))
# A depot is represented as a set containing the associated sequences
depots = [m.set(nb_trucks) for _ in range(nb_depots)]
# All the sequences should be assigned to a depot
m.constraint(m.partition(depots))
route_costs = [None] * nb_trucks
sequence_used = [None] * nb_trucks
dist_routes = [None] * nb_trucks
associated_depot = [None] * nb_trucks
# Create LocalSolver arrays to be able to access them with "at" operators
demands = m.array(demands_data)
dist_matrix = m.array()
dist_depots = m.array()
quantity_served = m.array()
for i in range(nb_customers):
dist_matrix.add_operand(m.array(dist_matrix_data[i]))
dist_depots.add_operand(m.array(dist_depots_data[i]))
for r in range(nb_trucks):
sequence = customers_sequences[r]
c = m.count(sequence)
# A sequence is used if it serves at least one customer
sequence_used[r] = c > 0
# The "find" function gets the depot that is assigned to the sequence
associated_depot[r] = m.find(m.array(depots), r)
# The quantity needed in each sequence must not exceed the vehicle capacity
demand_lambda = m.lambda_function(lambda j: demands[j])
quantity_served.add_operand(m.sum(sequence, demand_lambda))
m.constraint(quantity_served[r] <= vehicle_capacity)
# Distance traveled by each truck
dist_lambda = m.lambda_function(
lambda i: m.at(dist_matrix, sequence[i], sequence[i + 1]))
depot = associated_depot[r]
dist_routes[r] = m.sum(m.range(0, c - 1), dist_lambda) + m.iif(
sequence_used[r],
m.at(dist_depots, sequence[0], depot)
+ m.at(dist_depots, sequence[c - 1], depot),
0)
# The sequence cost is the sum of the opening cost and the sequence length
route_costs[r] = sequence_used[r] * opening_route_cost + dist_routes[r]
depot_cost = [None] * nb_depots
for d in range(nb_depots):
# A depot is open if at least one sequence starts from there
depot_cost[d] = (m.count(depots[d]) > 0) * opening_depots_cost[d]
# The total demand served by a depot must not exceed its capacity
depot_lambda = m.lambda_function(lambda r: quantity_served[r])
depot_quantity = m.sum(depots[d], depot_lambda)
m.constraint(depot_quantity <= capacity_depots[d])
depots_cost = m.sum(depot_cost)
routing_cost = m.sum(route_costs)
totalCost = routing_cost + depots_cost
m.minimize(totalCost)
m.close()
# Parameterize the solver
ls.param.time_limit = int(str_time_limit)
ls.solve()
if sol_file != None:
with open(sol_file, 'w') as file:
file.write("File name: %s; totalCost = %d \n" % (instance_file, totalCost.value))
for r in range(nb_trucks):
if sequence_used[r].value:
file.write("Route %d, assigned to depot %d: " % (r, associated_depot[r].value))
for customer in customers_sequences[r].value:
file.write("%d " % customer)
file.write("\n")
def read_input_lrp_dat(filename):
file_it = iter(read_elem(filename))
nb_customers = int(next(file_it))
nb_depots = int(next(file_it))
x_depot = [None] * nb_depots
y_depot = [None] * nb_depots
for i in range(nb_depots):
x_depot[i] = int(next(file_it))
y_depot[i] = int(next(file_it))
x_customer = [None] * nb_customers
y_customer = [None] * nb_customers
for i in range(nb_customers):
x_customer[i] = int(next(file_it))
y_customer[i] = int(next(file_it))
vehicle_capacity = int(next(file_it))
capacity_depots = [None] * nb_depots
for i in range(nb_depots):
capacity_depots[i] = int(next(file_it))
demands = [None] * nb_customers
for i in range(nb_customers):
demands[i] = int(next(file_it))
temp_opening_cost_depot = [None] * nb_depots
for i in range(nb_depots):
temp_opening_cost_depot[i] = float(next(file_it))
temp_opening_route_cost = int(next(file_it))
are_cost_double = int(next(file_it))
opening_depots_cost = [None] * nb_depots
if are_cost_double == 1:
opening_depots_cost = temp_opening_cost_depot
opening_route_cost = temp_opening_route_cost
else:
opening_route_cost = round(temp_opening_route_cost)
for i in range(nb_depots):
opening_depots_cost[i] = round(temp_opening_cost_depot[i])
distance_customers = compute_distance_matrix(x_customer, y_customer, are_cost_double)
distance_customers_depots = compute_distance_depot(x_customer, y_customer,
x_depot, y_depot, are_cost_double)
return nb_customers, nb_depots, vehicle_capacity, opening_route_cost, demands, \
capacity_depots, opening_depots_cost, distance_customers, distance_customers_depots
# Compute the distance matrix
def compute_distance_matrix(customers_x, customers_y, are_cost_double):
nb_customers = len(customers_x)
dist_customers = [[None for _ in range(nb_customers)] for _ in range(nb_customers)]
for i in range(nb_customers):
dist_customers[i][i] = 0
for j in range(nb_customers):
dist = compute_dist(customers_x[i], customers_x[j],
customers_y[i], customers_y[j], are_cost_double)
dist_customers[i][j] = dist
dist_customers[j][i] = dist
return dist_customers
# Compute the distance depot matrix
def compute_distance_depot(customers_x, customers_y, depot_x, depot_y, are_cost_double):
nb_customers = len(customers_x)
nb_depots = len(depot_x)
distance_customers_depots = [[None for _ in range(nb_depots)] for _ in range(nb_customers)]
for i in range(nb_customers):
for d in range(nb_depots):
dist = compute_dist(customers_x[i], depot_x[d],
customers_y[i], depot_y[d], are_cost_double)
distance_customers_depots[i][d] = dist
return distance_customers_depots
def compute_dist(xi, xj, yi, yj, are_cost_double):
dist = math.sqrt(math.pow(xi - xj, 2) + math.pow(yi - yj, 2))
if are_cost_double == 0:
dist = math.ceil(100 * dist)
return dist
def read_input_lrp(filename):
if filename.endswith(".dat"):
return read_input_lrp_dat(filename)
else:
raise Exception("Unknown file format")
if __name__ == '__main__':
if len(sys.argv) < 2:
print("Usage: python location_routing_problem.py input_file \
[output_file] [time_limit]")
sys.exit(1)
instance_file = sys.argv[1]
sol_file = sys.argv[2] if len(sys.argv) > 2 else None
str_time_limit = sys.argv[3] if len(sys.argv) > 3 else "20"
main(instance_file, str_time_limit, sol_file)
- Compilation / Execution (Windows)
- cl /EHsc location_routing_problem.cpp -I%LS_HOME%\include /link %LS_HOME%\bin\localsolver125.liblocation_routing_problem instances/coordChrist100.dat
- Compilation / Execution (Linux)
- g++ location_routing_problem.cpp -I/opt/localsolver_12_5/include -llocalsolver125 -lpthread -o location_routing_problem./location_routing_problem instances/coordChrist100.dat
#include "localsolver.h"
#include <cmath>
#include <cstring>
#include <fstream>
#include <iostream>
#include <numeric>
#include <vector>
using namespace localsolver;
using namespace std;
class LocationRoutingProblem {
public:
// LocalSolver
LocalSolver localsolver;
// Number of customers
int nbCustomers;
// Customers coordinates
vector<int> xCustomers;
vector<int> yCustomers;
// Customers demands
vector<double> demandsData;
// Number of depots
int nbDepots;
// Depots coordinates
vector<double> xDepots;
vector<double> yDepots;
// Capacity of depots
vector<double> depotsCapacity;
// Cost of opening a depot
vector<double> openingDepotsCost;
// Number of trucks
int nbTrucks;
// Capacity of trucks
int truckCapacity;
// Cost of opening a route
int openingRouteCost;
// Is the route used ?
vector<LSExpression> sequenceUsed;
// What is the depot of the route ?
vector<LSExpression> associatedDepot;
// Distance matrixes
vector<vector<double>> distMatrixData;
vector<vector<double>> distDepotsData;
int areCostDouble;
// Decision variables
vector<LSExpression> customersSequences;
vector<LSExpression> depots;
// Sum of all the costs
LSExpression totalCost;
void readInstance(const char* fileName) { readInputLrp(fileName); }
void solve(const char* limit) {
// Declare the optimization model
LSModel m = localsolver.getModel();
int minNbTrucks = ceil(accumulate(demandsData.begin(), demandsData.end(), 0) / truckCapacity);
nbTrucks = ceil(1.5 * minNbTrucks);
// A sequence is represented as a list containing the customers in the order they are visited
customersSequences.resize(nbTrucks);
for (int i = 0; i < nbTrucks; ++i) {
customersSequences[i] = m.listVar(nbCustomers);
}
// All customers should be assigned to a sequence
m.constraint(m.partition(customersSequences.begin(), customersSequences.end()));
// A depot is represented as a set containing the associated customersSequences
depots.resize(nbDepots);
for (int d = 0; d < nbDepots; ++d) {
depots[d] = m.setVar(nbTrucks);
}
// All the customersSequences should be assigned to a depot
m.constraint(m.partition(depots.begin(), depots.end()));
vector<LSExpression> distRoutes;
vector<LSExpression> routeCosts;
distRoutes.resize(nbTrucks);
sequenceUsed.resize(nbTrucks);
routeCosts.resize(nbTrucks);
associatedDepot.resize(nbTrucks);
// Create LocalSolver arrays to be able to access them with "at" operators
LSExpression quantityServed = m.array();
LSExpression demands = m.array(demandsData.begin(), demandsData.end());
LSExpression distMatrix = m.array();
LSExpression distDepots = m.array();
for (int i = 0; i < nbCustomers; ++i) {
distMatrix.addOperand(m.array(distMatrixData[i].begin(), distMatrixData[i].end()));
distDepots.addOperand(m.array(distDepotsData[i].begin(), distDepotsData[i].end()));
}
for (int r = 0; r < nbTrucks; ++r) {
LSExpression sequence = customersSequences[r];
LSExpression c = m.count(sequence);
// A sequence is used if it serves at least one customer
sequenceUsed[r] = c > 0;
// The "find" function gets the depot assigned to the sequence
associatedDepot[r] = m.find(m.array(depots.begin(), depots.end()), r);
LSExpression demandLambda = m.lambdaFunction([&](LSExpression j) { return demands[j]; });
quantityServed.addOperand(m.sum(sequence, demandLambda));
// The quantity needed in each sequence must not exceed the vehicle capacity
m.constraint(quantityServed[r] <= truckCapacity);
LSExpression distLambda =
m.lambdaFunction([&](LSExpression i) { return m.at(distMatrix, sequence[i], sequence[i + 1]); });
distRoutes[r] = m.iif(sequenceUsed[r],
m.at(distDepots, sequence[0], associatedDepot[r]) +
m.at(distDepots, sequence[c - 1], associatedDepot[r]),
0) +
m.sum(m.range(0, c - 1), distLambda);
// The sequence cost is the sum of the opening cost and the sequence length
routeCosts[r] = sequenceUsed[r] * openingRouteCost + distRoutes[r];
}
vector<LSExpression> depotCost;
depotCost.resize(nbDepots);
for (int d = 0; d < nbDepots; ++d) {
// A depot is open if at least a sequence starts from there
depotCost[d] = openingDepotsCost[d] * (m.count(depots[d]) > 0);
LSExpression depotLambda = m.lambdaFunction([&](LSExpression r) { return quantityServed[r]; });
LSExpression depotQuantity = m.sum(depots[d], depotLambda);
// The total demand served by a depot must not exceed its capacity
m.constraint(depotQuantity <= depotsCapacity[d]);
}
LSExpression depotsCost = m.sum(depotCost.begin(), depotCost.end());
LSExpression routingCost = m.sum(routeCosts.begin(), routeCosts.end());
totalCost = routingCost + depotsCost;
m.minimize(totalCost);
m.close();
localsolver.getParam().setTimeLimit(atoi(limit));
localsolver.solve();
}
/* Write the solution in a file */
void writeSolution(const char* inFile, const string& solFile) {
ofstream file;
file.exceptions(ofstream::failbit | ofstream::badbit);
file.open(solFile.c_str());
file << "File name: " << inFile << "; total cost = " << totalCost.getDoubleValue() << endl;
for (int r = 0; r < nbTrucks; ++r) {
if (sequenceUsed[r].getValue()) {
file << "Sequence " << r << ", assigned to depot " << associatedDepot[r].getValue() << " : ";
LSCollection customersCollection = customersSequences[r].getCollectionValue();
for (lsint i = 0; i < customersCollection.count(); ++i) {
file << customersCollection[i] << " ";
}
file << endl;
}
}
}
private:
void readInputLrp(const char* fileName) {
string file = fileName;
ifstream infile(file.c_str());
if (!infile.is_open()) {
throw std::runtime_error("File cannot be opened.");
}
infile >> nbCustomers;
xCustomers.resize(nbCustomers);
yCustomers.resize(nbCustomers);
demandsData.resize(nbCustomers);
distMatrixData.resize(nbCustomers);
distDepotsData.resize(nbCustomers);
infile >> nbDepots;
xDepots.resize(nbDepots);
yDepots.resize(nbDepots);
depotsCapacity.resize(nbDepots);
openingDepotsCost.resize(nbDepots);
for (int i = 0; i < nbDepots; ++i) {
infile >> xDepots[i];
infile >> yDepots[i];
}
for (int i = 0; i < nbCustomers; ++i) {
infile >> xCustomers[i];
infile >> yCustomers[i];
}
infile >> truckCapacity;
for (int i = 0; i < nbDepots; ++i) {
infile >> depotsCapacity[i];
}
for (int i = 0; i < nbCustomers; ++i) {
infile >> demandsData[i];
}
vector<double> tempOpeningCostDepots;
tempOpeningCostDepots.resize(nbDepots);
for (int i = 0; i < nbDepots; ++i) {
infile >> tempOpeningCostDepots[i];
}
int tempOpeningCostRoute;
infile >> tempOpeningCostRoute;
infile >> areCostDouble;
infile.close();
if (areCostDouble == 1) {
openingRouteCost = tempOpeningCostRoute;
for (int i = 0; i < nbDepots; ++i) {
openingDepotsCost[i] = tempOpeningCostDepots[i];
}
} else {
openingRouteCost = round(tempOpeningCostRoute);
for (int i = 0; i < nbDepots; ++i) {
openingDepotsCost[i] = round(tempOpeningCostDepots[i]);
}
}
computeDistanceMatrix();
}
void computeDistanceMatrix() {
for (int i = 0; i < nbCustomers; ++i) {
distMatrixData[i].resize(nbCustomers);
distDepotsData[i].resize(nbDepots);
for (int j = 0; j < nbCustomers; ++j) {
distMatrixData[i][j] =
computeDist(xCustomers[i], yCustomers[i], xCustomers[j], yCustomers[j], areCostDouble);
}
for (int d = 0; d < nbDepots; ++d) {
distDepotsData[i][d] = computeDist(xCustomers[i], yCustomers[i], xDepots[d], yDepots[d], areCostDouble);
}
}
}
double computeDist(int xi, int yi, int xj, int yj, int areCostDouble) {
double dist = sqrt(pow(xi - xj, 2) + pow(yi - yj, 2));
if (areCostDouble == 0) {
dist = ceil(dist * 100);
}
return dist;
}
};
int main(int argc, char** argv) {
if (argc < 2) {
cerr << "Usage: ./lrp inputFile [outputFile] [timeLimit]" << endl;
return 1;
}
const char* instanceFile = argv[1];
const char* solFile = argc > 2 ? argv[2] : NULL;
const char* strTimeLimit = argc > 3 ? argv[3] : "20";
try {
LocationRoutingProblem model;
model.readInstance(instanceFile);
model.solve(strTimeLimit);
if (solFile != NULL)
model.writeSolution(instanceFile, solFile);
} catch (const std::exception& e) {
std::cerr << "An error occured: " << e.what() << endl;
}
return 0;
}
- Compilation / Execution (Windows)
- copy %LS_HOME%\bin\localsolvernet.dll .csc LocationRoutingProblem.cs /reference:localsolvernet.dllLocationRoutingProblem instances\coordChrist100.dat
using System;
using System.IO;
using System.Collections.Generic;
using localsolver;
public class LocationRoutingProblem : IDisposable
{
// LocalSolver
LocalSolver localsolver;
// Number of customers
int nbCustomers;
// Customers coordinates
int[] xCustomers;
int[] yCustomers;
// Customers demands
int[] demandsData;
// Number of depots
int nbDepots;
// Depots coordinates
int[] xDepots;
int[] yDepots;
// Capacity of depots
int[] depotsCapacity;
double[] openingDepotsCost;
// Number of trucks
int nbTrucks;
// Capacity of a truck
int truckCapacity;
// Cost of opening a route
int openingRouteCost;
// Is the route used ?
LSExpression[] sequenceUsed;
// What is the depot of the route ?
LSExpression[] associatedDepot;
// Distance matrixes
double[][] distMatrixData;
double[][] distDepotsData;
int areCostDouble;
// Decision variables
LSExpression[] customersSequences;
LSExpression[] depots;
// Sum of all the costs
LSExpression totalCost;
public LocationRoutingProblem()
{
localsolver = new LocalSolver();
}
/* Read instance data */
void ReadInstance(string fileName)
{
using (StreamReader input = new StreamReader(fileName))
{
string[] line;
line = ReadNextLine(input);
nbCustomers = int.Parse(line[0]);
xCustomers = new int[nbCustomers];
yCustomers = new int[nbCustomers];
demandsData = new int[nbCustomers];
distMatrixData = new double[nbCustomers][];
line = ReadNextLine(input);
nbDepots = int.Parse(line[0]);
xDepots = new int[nbDepots];
yDepots = new int[nbDepots];
depotsCapacity = new int[nbDepots];
openingDepotsCost = new double[nbDepots];
distDepotsData = new double[nbCustomers][];
line = ReadNextLine(input);
for (int i = 0; i < nbDepots; ++i)
{
line = ReadNextLine(input);
xDepots[i] = int.Parse(line[0]);
yDepots[i] = int.Parse(line[1]);
}
line = ReadNextLine(input);
for (int i = 0; i < nbCustomers; ++i)
{
line = ReadNextLine(input);
xCustomers[i] = int.Parse(line[0]);
yCustomers[i] = int.Parse(line[1]);
}
line = ReadNextLine(input);
line = ReadNextLine(input);
truckCapacity = int.Parse(line[0]);
line = ReadNextLine(input);
for (int i = 0; i < nbDepots; ++i)
{
line = ReadNextLine(input);
depotsCapacity[i] = int.Parse(line[0]);
}
line = ReadNextLine(input);
for (int i = 0; i < nbCustomers; ++i)
{
line = ReadNextLine(input);
demandsData[i] = int.Parse(line[0]);
}
line = ReadNextLine(input);
double[] tempOpeningCostDepots = new double[nbDepots];
for (int i = 0; i < nbDepots; ++i)
{
line = ReadNextLine(input);
tempOpeningCostDepots[i] = double.Parse(
line[0],
System.Globalization.CultureInfo.InvariantCulture
);
}
line = ReadNextLine(input);
line = ReadNextLine(input);
double tempOpeningCostRoute = double.Parse(line[0]);
line = ReadNextLine(input);
line = ReadNextLine(input);
areCostDouble = int.Parse(line[0]);
if (areCostDouble == 1)
{
openingRouteCost = (int)tempOpeningCostRoute;
for (int i = 0; i < nbDepots; ++i)
openingDepotsCost[i] = tempOpeningCostDepots[i];
}
else
{
openingRouteCost = (int)Math.Round(tempOpeningCostRoute);
for (int i = 0; i < nbDepots; ++i)
openingDepotsCost[i] = Math.Round(tempOpeningCostDepots[i]);
}
DistanceMatrixes();
}
}
void DistanceMatrixes()
{
for (int i = 0; i < nbCustomers; ++i)
{
distMatrixData[i] = new double[nbCustomers];
for (int j = 0; j < nbCustomers; ++j)
{
distMatrixData[i][j] = ComputeDistance(
xCustomers[i],
yCustomers[i],
xCustomers[j],
yCustomers[j]
);
}
distDepotsData[i] = new double[nbDepots];
for (int d = 0; d < nbDepots; ++d)
{
distDepotsData[i][d] = ComputeDistance(
xCustomers[i],
yCustomers[i],
xDepots[d],
yDepots[d]
);
}
}
}
double ComputeDistance(int xi, int yi, int xj, int yj)
{
double dist = Math.Sqrt(Math.Pow(xi - xj, 2) + Math.Pow(yi - yj, 2));
if (areCostDouble == 0)
dist = Math.Ceiling(dist * 100);
return dist;
}
public void Dispose()
{
if (localsolver != null)
localsolver.Dispose();
}
void Solve(int limit)
{
// Declare the optimization model
LSModel m = localsolver.GetModel();
int totalDemand = 0;
foreach (int dem in demandsData)
totalDemand += dem;
int minNbTrucks = (int)Math.Ceiling((double)totalDemand / truckCapacity);
nbTrucks = (int)Math.Ceiling(1.5 * minNbTrucks);
customersSequences = new LSExpression[nbTrucks];
depots = new LSExpression[nbDepots];
// A route is represented as a list containing the customers in the order they are visited
for (int r = 0; r < nbTrucks; ++r)
customersSequences[r] = m.List(nbCustomers);
// All customers should be assigned to a route
m.Constraint(m.Partition(customersSequences));
// A depot is represented as a set containing the associated sequences
for (int d = 0; d < nbDepots; ++d)
depots[d] = m.Set(nbTrucks);
// All the sequences should be assigned to a depot
m.Constraint(m.Partition(depots));
// Create LocalSolver arrays to be able to access them with "at" operators
LSExpression demands = m.Array(demandsData);
LSExpression distMatrix = m.Array(distMatrixData);
LSExpression distDepots = m.Array(distDepotsData);
sequenceUsed = new LSExpression[nbTrucks];
LSExpression[] routeCosts = new LSExpression[nbTrucks];
LSExpression[] distRoutes = new LSExpression[nbTrucks];
associatedDepot = new LSExpression[nbTrucks];
LSExpression quantityServed = m.Array();
for (int r = 0; r < nbTrucks; ++r)
{
LSExpression sequence = customersSequences[r];
LSExpression c = m.Count(sequence);
// A sequence is used if it serves at least one customer
sequenceUsed[r] = m.Gt(c, 0);
// The "find" function gets the depot that is assigned to the sequence
associatedDepot[r] = m.Find(m.Array(depots), r);
LSExpression demandLambda = m.LambdaFunction(j => demands[j]);
quantityServed.AddOperand(m.Sum(sequence, demandLambda));
// The quantity needed in each sequence must not exceed the vehicle capacity
m.Constraint(quantityServed[r] <= truckCapacity);
LSExpression distLambda = m.LambdaFunction(
i => m.At(distMatrix, m.At(sequence, i), m.At(sequence, i + 1))
);
distRoutes[r] =
m.Sum(m.Range(0, c - 1), distLambda)
+ m.If(
sequenceUsed[r],
m.At(distDepots, m.At(sequence, 0), associatedDepot[r])
+ m.At(distDepots, m.At(sequence, c - 1), associatedDepot[r]),
0
);
// The sequence cost is the sum of the opening cost and the sequence length
routeCosts[r] = distRoutes[r] + openingRouteCost * sequenceUsed[r];
}
LSExpression[] depotCost = new LSExpression[nbDepots];
for (int d = 0; d < nbDepots; ++d)
{
// A depot is open if at least a sequence starts from there
depotCost[d] = openingDepotsCost[d] * (m.Count(depots[d]) > 0);
LSExpression depotLambda = m.LambdaFunction(r => m.At(quantityServed, r));
LSExpression depotQuantity = m.Sum(depots[d], depotLambda);
// The total demand served by a depot must not exceed its capacity
m.Constraint(depotQuantity <= depotsCapacity[d]);
}
LSExpression depotsCost = m.Sum(depotCost);
LSExpression routingCost = m.Sum(routeCosts);
totalCost = routingCost + depotsCost;
m.Minimize(totalCost);
m.Close();
localsolver.GetParam().SetTimeLimit(limit);
localsolver.Solve();
}
/* Write the solution in a file */
void WriteSolution(string infile, string outfile)
{
using (StreamWriter output = new StreamWriter(outfile))
{
output.WriteLine(
"File name: " + infile + "; total cost = " + totalCost.GetDoubleValue()
);
for (int r = 0; r < nbTrucks; ++r)
{
if (sequenceUsed[r].GetValue() != 0)
{
output.Write(
"Route " + r + ", assigned to depot " + associatedDepot[r].GetValue() + ": "
);
LSCollection customersCollection = customersSequences[r].GetCollectionValue();
for (int i = 0; i < customersCollection.Count(); ++i)
output.Write(customersCollection[i] + " ");
output.WriteLine();
}
}
}
}
String[] ReadNextLine(StreamReader input)
{
return input.ReadLine().Split(new[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
}
public static void Main(string[] args)
{
if (args.Length < 1)
{
Console.WriteLine("Usage: LocationRoutingProblem inputFile [outputFile] [timeLimit]");
Environment.Exit(1);
}
string instanceFile = args[0];
string strOutput = args.Length > 1 ? args[1] : null;
string strTimeLimit = args.Length > 2 ? args[2] : "20";
using (LocationRoutingProblem model = new LocationRoutingProblem())
{
model.ReadInstance(instanceFile);
model.Solve(int.Parse(strTimeLimit));
if (strOutput != null)
model.WriteSolution(instanceFile, strOutput);
}
}
}
- Compilation / Execution (Windows)
- javac LocationRoutingProblem.java -cp %LS_HOME%\bin\localsolver.jarjava -cp %LS_HOME%\bin\localsolver.jar;. LocationRoutingProblem instances\coordChrist100.dat
- Compilation / Execution (Linux)
- javac LocationRoutingProblem.java -cp /opt/localsolver_12_5/bin/localsolver.jarjava -cp /opt/localsolver_12_5/bin/localsolver.jar:. LocationRoutingProblem instances/coordChrist100.dat
import localsolver.*;
import java.io.*;
import java.io.File;
import java.io.FileNotFoundException;
import java.util.Scanner;
public class LocationRoutingProblem {
// LocalSolver
private final LocalSolver localsolver;
// Number of customers
private int nbCustomers;
// Customers coordinates
private int[] xCustomers;
private int[] yCustomers;
// Customers demands
private int[] demandsData;
// Number of depots
private int nbDepots;
// Depots coordinates
private int[] xDepots;
private int[] yDepots;
// Capacity of depots
private int[] depotsCapacity;
// Cost of opening a depot
private double[] openingDepotsCost;
// Number of trucks
private int nbTrucks;
// Capacity of trucks
private int truckCapacity;
// Cost of opening a route
private int openingRouteCost;
// Is the sequence used ?
private LSExpression[] sequenceUsed;
// What is the depot of the sequence ?
private LSExpression[] associatedDepot;
// Distance matrixes
private double[][] distMatrixData;
private double[][] distDepotsData;
private int areCostDouble;
// Decision variables
private LSExpression[] customersSequences;
private LSExpression[] depots;
// Objective value
private LSExpression totalCost;
private LocationRoutingProblem(LocalSolver localsolver) {
this.localsolver = localsolver;
}
private void solve(int limit) {
// Declare the optimization model
LSModel m = localsolver.getModel();
int totalDemand = 0;
for (int d : demandsData) {
totalDemand += d;
}
int minNbTrucks = (int) Math.ceil(totalDemand / truckCapacity);
int nbTrucks = (int) Math.ceil(1.5 * minNbTrucks);
customersSequences = new LSExpression[nbTrucks];
depots = new LSExpression[nbDepots];
// A route is represented as a list containing the customers in the order they are visited
for (int i = 0; i < nbTrucks; ++i) {
customersSequences[i] = m.listVar(nbCustomers);
}
// All customers should be assigned to a route
m.constraint(m.partition(customersSequences));
// A depot is represented as a set containing the associated customersSequences
for (int d = 0; d < nbDepots; ++d) {
depots[d] = m.setVar(nbTrucks);
}
// All the customersSequences should be assigned to a depot
m.constraint(m.partition(depots));
LSExpression demands = m.array(demandsData);
LSExpression distMatrix = m.array(distMatrixData);
LSExpression distDepots = m.array(distDepotsData);
sequenceUsed = new LSExpression[nbTrucks];
LSExpression[] routeCosts = new LSExpression[nbTrucks];
LSExpression[] distRoutes = new LSExpression[nbTrucks];
associatedDepot = new LSExpression[nbTrucks];
LSExpression quantityServed = m.array();
for (int r = 0; r < nbTrucks; ++r) {
LSExpression sequence = customersSequences[r];
LSExpression c = m.count(sequence);
// A sequence is used if it serves at least one customer
sequenceUsed[r] = m.gt(c, 0);
// The "find" function gets the depot that is assigned to the sequence
associatedDepot[r] = m.find(m.array(depots), r);
LSExpression demandLambda = m.lambdaFunction(j -> m.at(demands, j));
quantityServed.addOperand(m.sum(sequence, demandLambda));
// The quantity needed in each sequence must not exceed the vehicle capacity
m.constraint(m.leq(m.at(quantityServed, r), truckCapacity));
LSExpression distLambda = m
.lambdaFunction(i -> m.at(distMatrix, m.at(sequence, i), m.at(sequence, m.sum(i, 1))));
distRoutes[r] = m.sum(m.sum(m.range(0, m.sub(c, 1)), distLambda),
m.iif(sequenceUsed[r], m.sum(m.at(distDepots, m.at(sequence, 0), associatedDepot[r]),
m.at(distDepots, m.at(sequence, m.sub(c, 1)), associatedDepot[r])), 0));
// The sequence cost is the sum of the opening cost and the sequence length
routeCosts[r] = m.sum(distRoutes[r], m.prod(openingRouteCost, sequenceUsed[r]));
}
LSExpression[] depotCost = new LSExpression[nbDepots];
for (int d = 0; d < nbDepots; ++d) {
// A depot is open if at least a sequence starts from there
depotCost[d] = m.prod(openingDepotsCost[d], m.gt(m.count(depots[d]), 0));
LSExpression depotLambda = m.lambdaFunction(r -> m.at(quantityServed, r));
LSExpression depotQuantity = m.sum(depots[d], depotLambda);
// The total demand served by a depot must not exceed its capacity
m.constraint(m.leq(depotQuantity, depotsCapacity[d]));
}
LSExpression depotsCost = m.sum(depotCost);
LSExpression routingCost = m.sum(routeCosts);
totalCost = m.sum(routingCost, depotsCost);
m.minimize(totalCost);
m.close();
localsolver.getParam().setTimeLimit(limit);
localsolver.solve();
}
/* Write the solution in a file */
private void writeSolution(String infile, String outfile) throws IOException {
try (PrintWriter output = new PrintWriter(outfile)) {
output.println("File name: " + infile + "; total cost = " + totalCost.getDoubleValue());
for (int r = 0; r < nbTrucks; ++r) {
if (sequenceUsed[r].getIntValue() != 0) {
output.print("Route " + r + ", assigned to depot " + associatedDepot[r].getIntValue() + " : ");
LSCollection customersCollection = customersSequences[r].getCollectionValue();
for (int i = 0; i < customersCollection.count(); ++i) {
output.print(customersCollection.get(i) + " ");
}
output.println();
}
}
} catch (FileNotFoundException e) {
System.out.println("An error occurred.");
e.printStackTrace();
}
}
private void readInstanceLrp(String fileName) throws IOException {
try (Scanner infile = new Scanner(new File(fileName))) {
nbCustomers = infile.nextInt();
xCustomers = new int[nbCustomers];
yCustomers = new int[nbCustomers];
demandsData = new int[nbCustomers];
distMatrixData = new double[nbCustomers][];
distDepotsData = new double[nbCustomers][];
nbDepots = infile.nextInt();
xDepots = new int[nbDepots];
yDepots = new int[nbDepots];
depotsCapacity = new int[nbDepots];
openingDepotsCost = new double[nbDepots];
for (int i = 0; i < nbDepots; ++i) {
xDepots[i] = infile.nextInt();
yDepots[i] = infile.nextInt();
}
for (int i = 0; i < nbCustomers; ++i) {
xCustomers[i] = infile.nextInt();
yCustomers[i] = infile.nextInt();
}
truckCapacity = infile.nextInt();
for (int i = 0; i < nbDepots; ++i) {
depotsCapacity[i] = infile.nextInt();
}
for (int i = 0; i < nbCustomers; ++i) {
demandsData[i] = infile.nextInt();
}
double[] tempOpeningCostDepots = new double[nbDepots];
for (int i = 0; i < nbDepots; ++i) {
if (infile.hasNext()) {
tempOpeningCostDepots[i] = Double.parseDouble(infile.next());
} else if (infile.hasNextInt()) {
tempOpeningCostDepots[i] = infile.nextInt();
}
}
int tempOpeningCostRoute = infile.nextInt();
areCostDouble = infile.nextInt();
if (areCostDouble == 1) {
openingRouteCost = tempOpeningCostRoute;
for (int i = 0; i < tempOpeningCostDepots.length; ++i) {
openingDepotsCost[i] = tempOpeningCostDepots[i];
}
} else {
openingRouteCost = Math.round(tempOpeningCostRoute);
for (int i = 0; i < tempOpeningCostDepots.length; ++i) {
openingDepotsCost[i] = Math.round(tempOpeningCostDepots[i]);
}
}
computeDistanceMatrix();
} catch (FileNotFoundException e) {
System.out.println("An error occurred.");
e.printStackTrace();
}
}
void computeDistanceMatrix() {
for (int i = 0; i < nbCustomers; ++i) {
distMatrixData[i] = new double[nbCustomers];
for (int j = 0; j < nbCustomers; ++j) {
distMatrixData[i][j] = computeDist(xCustomers[i], xCustomers[j], yCustomers[i], yCustomers[j],
areCostDouble);
}
distDepotsData[i] = new double[nbDepots];
for (int d = 0; d < nbDepots; ++d) {
distDepotsData[i][d] = computeDist(xCustomers[i], xDepots[d], yCustomers[i], yDepots[d], areCostDouble);
}
}
}
private double computeDist(int xi, int xj, int yi, int yj, int areCostDouble) {
double dist = Math.sqrt(Math.pow(xi - xj, 2) + Math.pow(yi - yj, 2));
if (areCostDouble == 0) {
dist = Math.ceil(100 * dist);
}
return dist;
}
public static void main(String[] args) {
if (args.length < 1) {
System.err.println("Usage: java LocationRoutingProblem inputFile [outputFile] [timeLimit]");
System.exit(1);
}
try (LocalSolver localsolver = new LocalSolver()) {
String instanceFile = args[0];
String strOutfile = args.length > 1 ? args[1] : null;
String strTimeLimit = args.length > 2 ? args[2] : "20";
LocationRoutingProblem model = new LocationRoutingProblem(localsolver);
model.readInstanceLrp(instanceFile);
model.solve(Integer.parseInt(strTimeLimit));
if (strOutfile != null)
model.writeSolution(instanceFile, strOutfile);
} catch (Exception ex) {
System.err.println(ex);
ex.printStackTrace();
System.exit(1);
}
}
}