Capacitated Facility Location (CFLP)

Principles learned

  • Create a set decision variable

  • Use a lambda expression to compute a sum on a set

Problem

../_images/capacitatedfacilitylocation.svg

The capacitated facility location problem is the problem of locating a number of facilities which have to serve sites, at minimum cost, where each site has a given demand and each facility has a given capacity. The cost for a facility is defined as the sum of a fixed opening price and the allocation prices due to transportation between the facility and every site it provides. This problem is also known as capacitated warehouse location problem, or capacitated p-median problem.

Download the example


Data

The data files cap61, cap62, cap63 and cap64 provided come from the OR-LIB. They are the test problems from Tables 1, 2 and 3 of J.E.Beasley “An algorithm for solving large capacitated warehouse location problems” European Journal of Operational Research 33 (1988) 314-325.

Each instance file contains the number of potential facilities, the number of sites, the table of capacities and opening prices for each potential facility, the table of demand for each site, and the matrix of allocation prices between each facility and each site.

Program

In this problem, you need to decide the set of sites each facility will provide. If this set is empty, the facility will remain close, otherwise it will open and the associated cost will be the opening price plus the allocation costs. The capacity constraint is that for every facility, the sum of demand from the sites it provides can not exceed its capacity.

Thus, the only decision variables are the sets facilityAssignment[f] that contains all the sites provided by facility f. The expressions costs[f] are created to store cost due to facility f. This cost is :

  • 0 if facilityAssignment[f] is an empty set.

  • opening price + the sum of allocation prices between f and the sites it provides otherwise

We want to minimize the total cost.

Execution:
localsolver capacitated_facility_location.lsp inFileName=instances/cap61 [lsTimeLimit=] [solFileName=]
use io;

/* Reads instance data */
function input() {
    local inFile = io.openRead(inFileName);
    nbMaxFacilities = inFile.readInt();
    nbSites = inFile.readInt();

    for [f in 0..nbMaxFacilities - 1] {
        //List of facilities capacities
        capacity[f] = inFile.readDouble();
        //List of fixed costs induced by the facilities opening
        openningPrice[f] = inFile.readDouble();
    }

    //Demand of each site
    demand[0..nbSites - 1] = inFile.readDouble();

    //Allocation Price btw sites and facilities
    allocationPrice[0..nbMaxFacilities - 1][0..nbSites - 1] = inFile.readDouble();
}

/* Declares the optimization model. */
function model(){
    //Facilities are represented by the set of sites they provide
    facilityAssignments[0..nbMaxFacilities - 1] <- set(nbSites);

    //Each site is covered by exactly one facility -> partition :
    constraint partition[f in 0..nbMaxFacilities - 1](facilityAssignments[f]);

    for [f in 0..nbMaxFacilities-1] {
        local facility <- facilityAssignments[f];
        local size <- count(facility);

        //capacity constraint
        constraint sum(facility, i => demand[i]) <= capacity[f];

        //Cost (Allocation price + Openning price)
        cost[f] <- sum(facility, i => allocationPrice[f][i]) + openningPrice[f] * (size > 0);
    }
    totalCost <- sum[f in 0..nbMaxFacilities - 1](cost[f]);

    //Objective : minimize total cost
    minimize totalCost;
}

/* Parameterizes the solver. */
function param() {
    if (lsTimeLimit == nil) lsTimeLimit = 20;
}

/* Writes the solution in a file with the following format:
   - value of the objective
   - indices of the open facilities followed by all the sites they provide */
function output() {
    if (solFileName == nil) return;
    local outfile = io.openWrite(solFileName);
    outfile.println(totalCost.value);
    for [f in 0..nbMaxFacilities - 1 : cost[f].value > 0]{
        outfile.println(f);
        for [s in facilityAssignments[f].value]
            outfile.print(s, " ");
        outfile.println();
    }
}
Execution (Windows)
set PYTHONPATH=%LS_HOME%\bin\python
python capacitated_facility_location.py instances\cap61
Execution (Linux)
export PYTHONPATH=/opt/localsolver_11_0/bin/python
python capacitated_facility_location.py instances/cap61
########## capacitated_facility_location.py ##########

import localsolver
import sys


def main(instanceFile, strTimeLimit, solFile):

    #
    # Reads instance data
    #
    nbMaxFacilities, nbSites, capacity, openingPrice, demand, allocationPrice = read_data(
        instanceFile)

    with localsolver.LocalSolver() as ls:

        #
        # Declares the optimization model
        #
        model = ls.model

        # Facilities are represented by the set of sites they provide
        facilityAssignments = [model.set(nbSites)
                               for f in range(nbMaxFacilities)]

        # Each site is covered by exactly one facility -> partition :
        model.constraint(model.partition(facilityAssignments))

        # Converting demand and allocationPrice into model array
        demandArray = model.array(demand)
        allocationPriceArray = model.array()
        for f in range(nbMaxFacilities):
            allocationPriceArray.add_operand(model.array(allocationPrice[f]))

        cost = [None] * nbMaxFacilities
        for f in range(nbMaxFacilities):
            facility = facilityAssignments[f]
            size = model.count(facility)

            # Capacity constraint
            demandSelector = model.lambda_function(lambda i: demandArray[i])
            model.constraint(
                model.sum(facility, demandSelector) <= capacity[f])

            # Cost (Allocation price + openning price)
            costSelector = model.lambda_function(
                lambda i: model.at(allocationPriceArray, f, i))
            cost[f] = model.sum(facility, costSelector) + \
                openingPrice[f] * (size > 0)

        # Objective : minimizing total cost
        totalCost = model.sum(cost)
        model.minimize(totalCost)

        model.close()

        #
        # Parameterizes the solver
        #
    
        ls.param.time_limit = int(strTimeLimit)

        ls.solve()

        #
        # Writes the solution in a file following the following format:
        # - value of the objective
        # - indices of the open facilities followed by all the sites they provide
        #
        if solFile:
            with open(solFile, 'w') as outputFile:
                outputFile.write("%d" % totalCost.value)
                for f in range(nbMaxFacilities):
                    if cost[f].value > 0:
                        outputFile.write("%d\n"% f)
                        for site in facilityAssignments[f].value:
                            outputFile.write("%d " % site)
                        outputFile.write("\n")


def read_elem(filename):
    with open(filename) as f:
        return [str(elem) for elem in f.read().split()]


def read_data(filename):
    file_it = iter(read_elem(filename))

    nbMaxFacilities = int(next(file_it))
    nbSites = int(next(file_it))

    capacity = []
    openingPrice = []
    demand = []
    allocationPrice = []

    for f in range(nbMaxFacilities):
        capacity.append(float(next(file_it)))
        openingPrice.append(float(next(file_it)))
        allocationPrice.append([])

    for s in range(nbSites):
        demand.append(float(next(file_it)))

    for f in range(nbMaxFacilities):
        for s in range(nbSites):
            allocationPrice[f].append(float(next(file_it)))

    return nbMaxFacilities, nbSites, capacity, openingPrice, demand, allocationPrice


if __name__ == '__main__':
    instanceFile = sys.argv[1]
    solFile = sys.argv[2] if len(sys.argv) > 2 else None
    strTimeLimit = sys.argv[3] if len(sys.argv) > 3 else "20"

    main(instanceFile, strTimeLimit, solFile)
Compilation / Execution (Windows)
cl /EHsc capacitated_facility_location.cpp -I%LS_HOME%\include /link %LS_HOME%\bin\localsolver110.lib
capacitated_facility_location instances\cap61
Compilation / Execution (Linux)
g++ capacitated_facility_location.cpp -I/opt/localsolver_11_0/include -llocalsolver110 -lpthread -o capacitated_facility_location
./capacitated_facility_location instances/pmed1.in
#include <iostream>
#include <fstream>
#include <vector>
#include <string.h>
#include "localsolver.h"

using namespace localsolver;
using namespace std;

class FacilityLocationCapacity {
public :

    //Data 
    int nbMaxFacilities;
    int nbSites;
    vector<lsdouble> capacity;
    vector<lsdouble> openingPrice;
    vector<lsdouble> demand;
    vector<vector<lsdouble>> allocationPrice;

    //Solver
    LocalSolver localsolver;

    vector<LSExpression> facilityAssignments;
    vector<LSExpression> cost;
    LSExpression totalCost;

    void solve(int limit) {
        // Declare the optimization model
        LSModel model = localsolver.getModel();

        // Table of sets representing the sites provided by each facility
        facilityAssignments.resize(nbMaxFacilities);
        for (int f = 0; f < nbMaxFacilities; f++) {
            facilityAssignments[f] = model.setVar(nbSites);
        }

        // Converting demand and allocationPrice into localsolver array
        LSExpression demandArray = model.array(demand.begin(), demand.end());
        LSExpression allocationPriceArray = model.array();
        for (int f = 0; f < nbMaxFacilities; f++) {
            allocationPriceArray.addOperand(model.array(allocationPrice[f].begin(), allocationPrice[f].end()));
        }

        // Partition constraint, as each site has to be provided by exactly one facility
        model.constraint(model.partition(facilityAssignments.begin(), facilityAssignments.end()));
        
        // Cost induced by facilities
        cost.resize(nbMaxFacilities);

        for (int f = 0; f < nbMaxFacilities; f++) {
            LSExpression facility = facilityAssignments[f];
            LSExpression size = model.count(facility);

            LSExpression demandSelector = model.createLambdaFunction([&](LSExpression i) { return demandArray[i]; });
            model.constraint(model.sum(facility, demandSelector) <= capacity[f]);

            LSExpression costSelector = model.createLambdaFunction([&](LSExpression i) { return model.at(allocationPriceArray, f, i); });
            cost[f] = model.sum(facility, costSelector) + (size > 0)*openingPrice[f];
        }
        totalCost = model.sum(cost.begin(), cost.end());

        // Objective : minimize total cost 
        model.minimize(totalCost);
        model.close();

        // Parametrizes the solver
        localsolver.getParam().setTimeLimit(limit);

        localsolver.solve();
    }

    void readInstance(const string& fileName) {
        ifstream infile;
        infile.exceptions(ifstream::failbit | ifstream::badbit);
        infile.open(fileName.c_str());

        infile >> nbMaxFacilities;
        infile >> nbSites;
        capacity.resize(nbMaxFacilities);
        openingPrice.resize(nbMaxFacilities);
        demand.resize(nbSites);
        allocationPrice.resize(nbMaxFacilities,vector<lsdouble>(nbSites));

        for (int f = 0; f < nbMaxFacilities; ++f) {
            infile >> capacity[f];
            infile >> openingPrice[f];
        }
        for (int s = 0; s < nbSites; ++s){
            infile >> demand[s];
        }
        for (int f = 0; f < nbMaxFacilities; ++f) {
            for (int s = 0; s < nbSites; ++s){
                infile >> allocationPrice[f][s];
            }
        }
        infile.close();
    }

    // Writes the solution in a file with the following format:
    //  - value of the objective
    //  - indices of the open facilities followed by all the sites they provide
    void writeSolution(const string& fileName) {
        ofstream outfile;
        outfile.exceptions(ofstream::failbit | ofstream::badbit);
        outfile.open(fileName.c_str());
        outfile << totalCost.getDoubleValue() << endl;
        for (int f = 0; f < nbMaxFacilities; f++) {
            if (cost[f].getDoubleValue() > 0) {
                outfile << f << endl;
                LSCollection sites_assigned = facilityAssignments[f].getCollectionValue();
                for (lsint s = 0; s < sites_assigned.count(); s++) {
                    outfile << sites_assigned[s] << " ";
                }
                outfile << endl;
            }
        }
    }
};

int main(int argc, char** argv){
    if (argc < 2) {
        cerr << "Usage: facilitylocation_capacity 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 {
        FacilityLocationCapacity model;
        model.readInstance(instanceFile);
        model.solve(atoi(strTimeLimit));
        if (solFile != NULL) model.writeSolution(solFile);
    } catch (const exception& e) {
        cerr << "An error occured " << e.what() << endl;
        return -1;
    }
}
Compilation / Execution (Windows)
copy %LS_HOME%\bin\localsolvernet.dll .
csc CapacitatedFacilityLocation.cs /reference:localsolvernet.dll
CapacitatedFacilityLocation instances\cap61
using System;
using System.IO;
using System.Collections.Generic;
using localsolver;

public class CapacitatedFacilityLocation : IDisposable
{
    //Solver
    LocalSolver localsolver = new LocalSolver();

    //Data
    int nbMaxFacilities;
    int nbSites;
    double[] capacity;
    double[] openingPrice;
    double[] demand;
    double[,] allocationPrice;

    //Variables
    LSExpression[] facilityAssignments;
    LSExpression[] cost;
    LSExpression totalCost;

    private string[] SplitInput(StreamReader input)
    {
        string line = input.ReadLine();
        if (line == null) return new string[0];
        return line.Split(new[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
    }

    void ReadInstance(string fileName)
    {
        using (StreamReader input = new StreamReader(fileName))
        {
            string[] splitted;
            splitted = SplitInput(input);
            nbMaxFacilities = int.Parse(splitted[0]);
            nbSites = int.Parse(splitted[1]);

            // Table initialization
            capacity = new double[nbMaxFacilities];
            openingPrice = new double[nbMaxFacilities];
            demand = new double[nbSites];
            allocationPrice = new double[nbMaxFacilities, nbSites];

            for (int f = 0; f < nbMaxFacilities; f++)
            {
                splitted = SplitInput(input);
                capacity[f] = double.Parse(splitted[0]);
                openingPrice[f] = double.Parse(splitted[1]);
            }
            splitted = SplitInput(input);
            for (int s = 0; s < nbSites; s++)
            {
                demand[s] = double.Parse(splitted[s]);
            }
            for (int f = 0; f < nbMaxFacilities; f++)
            {
                splitted = SplitInput(input);
                for (int s = 0; s < nbSites; s++)
                {
                    allocationPrice[f, s] = double.Parse(splitted[s]);
                }
            }
        }
    }

    public void Dispose()
    {
        if (localsolver != null)
            localsolver.Dispose();
    }


    void Solve(int limit)
    {

        LSModel model = localsolver.GetModel();

        //Facilities are represented by the sets of sites they provide
        facilityAssignments = new LSExpression[nbMaxFacilities];
        for (int f = 0; f < nbMaxFacilities; f++)
        {
            facilityAssignments[f] = model.Set(nbSites);
        }
        //Each site is provided by exactly one facility, so partition constraint
        model.Constraint(model.Partition(facilityAssignments));

        // Cost define the cost induced by each facility
        cost = new LSExpression[nbMaxFacilities];

        LSExpression demandArray = model.Array(demand);
        LSExpression allocationPriceArray = model.Array(allocationPrice);

        for (int f = 0; f < nbMaxFacilities; f++)
        {
            LSExpression facility = facilityAssignments[f];
            LSExpression size = model.Count(facility);

            // Capacity constraint
            LSExpression demandSelector = model.LambdaFunction(i => demandArray[i]);
            model.Constraint(model.Sum(facility, demandSelector) <= capacity[f]);

            // Cost calcul
            LSExpression cost_selector = model.LambdaFunction(i => allocationPriceArray[f, i]);
            cost[f] = model.Sum(facility, cost_selector) + model.If(size > 0, openingPrice[f], 0);
        }
        //Objective : minimizing total cost
        totalCost = model.Sum(cost);
        model.Minimize(totalCost);

        model.Close();

        // Parameterize the solver.
        localsolver.GetParam().SetTimeLimit(limit);

        localsolver.Solve();
    }

    // Writes the solution in a file with the following format:
    //  - value of the objective
    //  - indices of the open facilities followed by all the sites they provide
    void WriteSolution(string fileName)
    {
        using (StreamWriter output = new StreamWriter(fileName))
        {
            output.WriteLine(totalCost.GetDoubleValue());
            for (int f = 0; f < nbMaxFacilities; f++)
            {
                if (cost[f].GetDoubleValue() > 0)
                {
                    output.WriteLine(f);
                    LSCollection assigned_sites = facilityAssignments[f].GetCollectionValue();
                    for (int s = 0; s < assigned_sites.Count(); s++)
                    {
                        output.Write(assigned_sites[s] + " ");
                    }
                    output.WriteLine();
                }
            }
        }
    }

    public static void Main(string[] args)
    {
        string instanceFile = args[0];
        string outputFile = args.Length > 1 ? args[1] : null;
        string strTimeLimit = args.Length > 2 ? args[2] : "20";

        using (CapacitatedFacilityLocation model = new CapacitatedFacilityLocation())
        {
            model.ReadInstance(instanceFile);
            model.Solve(int.Parse(strTimeLimit));
            if (outputFile != null)
                model.WriteSolution(outputFile);
        }
    }
}
Compilation / Execution (Windows)
javac CapacitatedFacilityLocation.java -cp %LS_HOME%\bin\localsolver.jar
java -cp %LS_HOME%\bin\localsolver.jar;. CapacitatedFacilityLocation instances\cap61
Compilation / Execution (Linux)
javac CapacitatedFacilitylocation.java -cp /opt/localsolver_11_0/bin/localsolver.jar
java -cp /opt/localsolver_11_0/bin/localsolver.jar:. CapacitatedFacilityLocation instances/cap61
import java.util.*;
import java.io.*;
import localsolver.*;

public class CapacitatedFacilityLocation {

    // Data from the problem
    private int nbMaxFacilities;
    private int nbSites;
    private double[] capacity;
    private double[] openingPrice;
    private double[] demand;
    private double[][] allocationPrice;

    // Solver
    private final LocalSolver localsolver;

    // Variable
    private LSExpression[] facilityAssignments;

    // Objective
    private LSExpression[] cost;
    private LSExpression totalCost;

    private CapacitatedFacilityLocation(LocalSolver localsolver) {
        this.localsolver = localsolver;
    }

    private void readInstance(String fileName) throws IOException {
        try (Scanner input = new Scanner(new File(fileName))) {
            nbMaxFacilities = input.nextInt();
            nbSites = input.nextInt();
            capacity = new double[nbMaxFacilities];
            openingPrice = new double[nbMaxFacilities];
            demand = new double[nbSites];
            allocationPrice = new double[nbMaxFacilities][nbSites];

            for (int f = 0; f < nbMaxFacilities; f++) {
                capacity[f] = input.nextDouble();
                openingPrice[f] = input.nextDouble();
            }

            for (int s = 0; s < nbSites; s++) {
                demand[s] = input.nextDouble();
            }

            for (int f = 0; f < nbMaxFacilities; f++) {
                for (int s = 0; s < nbSites; s++) {
                    allocationPrice[f][s] = input.nextDouble();
                }
            }
        }
    }

    private void solve(int limit) {
        // Declare the optimization model
        LSModel model = localsolver.getModel();

        // Variable : table of sets
        facilityAssignments = new LSExpression[nbMaxFacilities];
        for (int f = 0; f < nbMaxFacilities; f++) {
            facilityAssignments[f] = model.setVar(nbSites);
        }

        // Table of cost induced by facililties
        cost = new LSExpression[nbMaxFacilities];

        // Constraint : sets of assignment must form a partition of sites
        model.constraint(model.partition(facilityAssignments));

        for (int f = 0; f < nbMaxFacilities; f++) {
            LSExpression facility = facilityAssignments[f];
            LSExpression size = model.count(facility);

            // Capacity constraint
            LSExpression demandSelector = model.lambdaFunction(i -> model.at(model.array(demand), i));
            model.constraint(model.leq(model.sum(facility, demandSelector), capacity[f]));

            // Calculating cost induced by facility f
            double[] facilityAllocationPrice = allocationPrice[f];
            LSExpression costSelector = model.lambdaFunction(i -> model.at(model.array(facilityAllocationPrice), i));
            cost[f] = model.sum(model.sum(facility, costSelector), model.iif(model.gt(size, 0), openingPrice[f], 0));
        }

        // Objective : minimize total cost
        totalCost = model.sum(cost);
        model.minimize(totalCost);

        model.close();
        localsolver.getParam().setTimeLimit(limit);
        localsolver.solve();
    }

    // Writes the solution in a file with the following format:
    // - value of the objective
    // - indices of the open facilities followed by all the sites they provide
    private void writeSolution(String file_name) throws IOException {
        try (PrintWriter output = new PrintWriter(file_name)) {
            output.println(totalCost.getDoubleValue());
            for (int f = 0; f < nbMaxFacilities; f++) {
                if (cost[f].getDoubleValue() > 0) {
                    output.println(f);
                    LSCollection sites_assigned = facilityAssignments[f].getCollectionValue();
                    for (int s = 0; s < sites_assigned.count(); s++) {
                        output.print(sites_assigned.get(s) + " ");
                    }
                    output.println();
                }
            }
        }
    }

    public static void main(String[] args) {
        try (LocalSolver localsolver = new LocalSolver()) {
            String instanceFile = args[0];
            String outputFile = args.length > 1 ? args[1] : null;
            String strTimeLimit = args.length > 2 ? args[2] : "20";

            CapacitatedFacilityLocation model = new CapacitatedFacilityLocation(localsolver);
            model.readInstance(instanceFile);
            model.solve(Integer.parseInt(strTimeLimit));
            if (outputFile != null) {
                model.writeSolution(outputFile);
            }
        } catch (Exception ex) {
            System.err.println(ex);
            ex.printStackTrace();
            System.exit(1);
        }
    }
}