Solves the Concurrent Flow variant of the Multi-Commodity Flow Problem with GLPK via LEMON interface.

mcfp_cf_lp_lemon_glpk.cc

/*
 * mcfp_cf_lp_lemon_glpk.cc
 * LEMON provides a simple interface to known solvers! In
 * this example the Concurrent Flow variant of the
 * Multi-Commodity Flow Problem is solved with the glpk
 * solver.
 *
 * The following libraries are needed to run this application:
 *
 * lemon
 * glpk
 *
 * Also note that the -std=c++11 flag is needed, as well.
 *
 * Please see the following links for more information on
 * how the problem is formulated:
 *
 *
 *  Created on: Oct 9, 2016
 *      Author: andrewandrepowell2
 */
#include <iostream>
#include <stdexcept>
#include <lemon/list_graph.h>
#include <lemon/lp.h>
#include <lemon/lp_base.h>
#include <utility>
#include <vector>
#include <map>
#include <string>

/* Just so this application looks cleaner
all the declarations are stored in the following
header file. */
#include "mcfp_cf_lp_lemon_glpk.h"

int main()
try
{
	/* It's important lp is declared locally.
	An error occurs when declared globally. */
	LinProb lp;
	lpp = &lp;

	/* Define the graph. */
	Node A = addNode("A");
	Node B = addNode("B");
	Node C = addNode("C");
	addArc(A,B,1.0,"AB");
	addArc(B,C,1.0,"BC");
	addArc(C,A,1.0,"CA");

	/* Define the Commodities. */
	addCommodity(A,C,0.5,"k0");
	addCommodity(C,B,0.5,"k1");

	/* Define the objective. */
	addObjective();

	/* Define the constraints. */
	addCapacityConstraint();
	addConserveConstraint();

	/* Solve the LP problem. Not sure what solver is applied,
	but I'm assuming it's Simplex. */
	lp.solve();

	/* Display the results. */
	displayResults();

	return 0;
}
catch ( std::exception& e )
{
	std::cerr << e.what() << std::endl;
}
catch (...)
{
	std::cerr << "A serious error occurred" << std::endl;
}

mcfp_cf_lp_lemon_glpk.h

template <typename T_> using Set = std::vector<T_>;
typedef lemon::Lp LinProb;
typedef LinProb::Col Variable;
typedef double Weight;
typedef Weight Demand;
typedef Weight Capacity;
typedef Variable Flow;
typedef lemon::ListDigraph Graph;
typedef Graph::Node Node;
typedef Graph::Arc Arc;
typedef Graph::ArcMap<Capacity> Capacities;
typedef int Commodity;
typedef Set<Commodity> Commodities;
typedef std::map<Commodity,Node> Sources;
typedef std::map<Commodity,Node> Sinks;
typedef std::map<Commodity,Demand> Demands;
typedef std::map<std::pair<Commodity,Arc>,Flow> Flows;
typedef std::string Label;
typedef Graph::NodeMap<Label> NodeLabels;
typedef Graph::ArcMap<Label> ArcLabels;
typedef std::map<Commodity,Label> ComLabels;

LinProb* lpp;
Graph g;
Capacities cs(g);
Commodities ks;
Sources kss;
Sinks kts;
Demands kds;
Flows kfs;

NodeLabels nls(g);
ArcLabels als(g);
ComLabels kls;

/* Adds Node with corresponding label. */
Node addNode( Label l )
{
	Node n = g.addNode();
	nls[n] = l;
	return n;
}

/* Adds Arc with corresponding capacity and label. */
Arc addArc( Node s, Node t, Capacity c, Label l )
{
	Arc a = g.addArc(s,t);
	cs[a] = c;
	als[a] = l;
	return a;
}

/* Adds commodity with corresponding source node,
sink node, demand, and label. This function also
ensures the non-negative constraint is added to the
LP problem. */
Commodity addCommodity( Node s, Node t, Demand d, Label l )
{
	Commodity k = ( ks.size() ) ? ks.back()+1 : 0;
	ks.push_back(k);
	kss[k] = s;
	kts[k] = t;
	kds[k] = d;
	kls[k] = l;
	for ( Graph::ArcIt a(g); a!=lemon::INVALID; ++a )
	{
		Variable flow = lpp->addCol();
		kfs[std::make_pair(k,a)] = flow;
		lpp->addRow( flow >= 0 );
	}
	return k;
}

/* Adds the objective to the LP problem. */
void addObjective()
{
	Variable z = lpp->addCol();
	for ( Commodities::iterator k=ks.begin(); k!=ks.end(); ++k )
	{
		LinProb::Expr expr;
		Demand d = kds[*k];
		Node s = kss[*k];
		for ( Graph::NodeIt v(g); v!=lemon::INVALID; ++v )
		{
			Arc a = lemon::findArc(g,s,v);
			if ( a!=lemon::INVALID )
				expr += kfs[std::make_pair(*k,a)]/d;
		}
		lpp->addRow( expr >= z );
	}
	lpp->max();
	lpp->obj( z );
}

/* Adds the Capacity Constraint to the LP problem. */
void addCapacityConstraint()
{
	for ( Graph::ArcIt a(g); a!=lemon::INVALID; ++a )
	{
		LinProb::Expr expr;
		for ( Commodities::iterator k=ks.begin(); k!=ks.end(); ++k )
			expr += kfs[std::make_pair(*k,a)];
		lpp->addRow( expr <= cs[a] );
	}
}

/* Adds the Conservation Constraint to the LP problem. */
void addConserveConstraint()
{
	for ( Commodities::iterator k=ks.begin(); k!=ks.end(); ++k )
	{
		Demand d = kds[*k];
		Node s = kss[*k];
		Node t = kts[*k];
		for ( Graph::NodeIt v(g); v!=lemon::INVALID; ++v )
		{
			LinProb::Expr expr;
			for ( Graph::NodeIt w(g); w!=lemon::INVALID; ++w )
			{
				Arc pos = lemon::findArc(g,v,w);
				Arc neg = lemon::findArc(g,w,v);
				if ( pos!=lemon::INVALID )
					expr += kfs[std::make_pair(*k,pos)];
				if ( neg!=lemon::INVALID )
					expr -= kfs[std::make_pair(*k,neg)];
			}
			if ( v!=s && v!=t )
			{
				lpp->addRow( expr==0 );
			}
			else if ( v==s )
			{
				lpp->addRow( expr==d );
			}
			else if ( v==t )
			{
				lpp->addRow( expr== -d );
			}
		}
	}
}

/* Displays the results after applying the solver. */
void displayResults()
{
	if ( lpp->primalType() == LinProb::OPTIMAL )
	{
		std::cout << "Displaying the flows..." << std::endl;
		std::cout << "Objective Value: " << lpp->primal() << std::endl;
		for ( Commodities::iterator k=ks.begin(); k!=ks.end(); ++k )
		{
			std::cout << "For Commodity " << kls[*k] << "..." << std::endl;
			for ( Graph::ArcIt a(g); a!=lemon::INVALID; ++a )
				std::cout << "Arc " << als[a] << ": " <<
				lpp->primal(kfs[std::make_pair(*k,a)]) << std::endl;
		}
	}
	else
	{
		std::cout << "Optimal value could not be found in feasible region!" << std::endl;
	}
}