Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Solves the Concurrent Flow variant of the Multi-Commodity Flow Problem in MATLAB with linprog.
%% Solve the Multi-Commodity Flow Problem in MATLAB.
%
% The goal of this script is to solve the Concurrent-Flow variant of the
% Multi-Commodity Flow Problem ( MCFP ) with MATLAB's Linear Programming
% ( LP ) tools.
%
% The documentation generated from this script is a tex file. The PDF can
% be generated from the tex, but \usagepackage{amsmath} needs to be added
% to the tex file's preamble.
% It's important the workspace is cleared
% before executing this script!
clear all;
close all;
%% Define the Flow Network.
%
% The flow network is defined as the following.
%
% <latex>
% \begin{align*}
% V & = \text{Node Set}\\
% E & = \text{Directed Edge ( Arc ) Set}\\
% e & = (u,v: u,v\in V) \\
% & \in E \\
% G & = \text{Directed Graph} \\
% & = (V,E) \\
% c & = \text{Capacity, } \forall e\in E \\
% \end{align*}
% </latex>
% Define the Nodes.
V = {'C','F','B'};
% Define the Arcs.
tuples = { ...
{ 'eAB' 'F' 'B' 1 }, ...
{ 'eBC' 'B' 'C' 1 }, ...
{ 'eCA' 'C' 'F' 1 } };
E = {};
for tuple=tuples
e = tuple{1}{1};
E{end+1} = e;
E_s.(e) = tuple{1}{2};
E_t.(e) = tuple{1}{3};
E_c.(e) = tuple{1}{4};
end
% Plot the graph.
figure;
plot(digraph(struct2cell(E_s),struct2cell(E_t)));
title( 'Digraph G' );
%% Define the Commodities.
%
% The commodities are each represented as a 3-tuple, containing a Source Node,
% a Target Node, and a Demand.
%
% <latex>
% \begin{align*}
% K & = \text{Commodity Set} \\
% & \subseteq V \times V \times \mathbf{R}_{\ge 0}\\
% s_k & = \text{Source Node of $k\in K$}\\
% t_k & = \text{Target Node of $k\in K$}\\
% d_k & = \text{Demand of $k\in K$}\\
% k & = (s_k,t_k,d_k) \\
% & \in K
% \end{align*}
% </latex>
% Define the commodities.
tuples = {
{ 'k0' 'F' 'C' 0.5 },...
{ 'k1' 'B' 'F' 0.5 } };
K = {};
for tuple=tuples
k = tuple{1}{1};
K{end+1} = k;
K_s.(k) = tuple{1}{2};
K_t.(k) = tuple{1}{3};
K_d.(k) = tuple{1}{4};
end
%% Structure MCFP.
%
% In order to run MATLAB's linprog, the problem needs to be adapted to fit
% under the following form.
%
% <latex>
% \begin{align*}
% \text{Objective: } & \min_x o^T x \\
% \text{Constraints: } & A_{in}x \le b_{in},\\
% & A_{eq} x = b_{eq}, \\
% & b_l \le x \le b_h
% \end{align*}
% </latex>
%
% $o$, $b_{in}$, and $b_{eq}$ are all constant vectors. $A_{in}$ and $A_{eq}$ are
% constant matrices. The bounds $b_l$ and $b_u$ can include constants or
% infinity. Finally, $x$ is the vector of variables, most of which in MCFP
% are the flow of commodity over each edge.
% Declaration starting indices.
col = 1;
rowin = 1;
roweq = 1;
% The following anonymous function makes string comparison easier.
cs = @(s1,s2)strcmp(s1,s2);
%% MCFP Variables.
%
% The MCFP problem has the following MCFP variables.
%
% <latex>
% \begin{align*}
% f_{k,e} & = \text{Flow of Commodity $k$ at Edge $e$}, \forall k\in K, \forall e\in E \\
% & \ge 0 \\
% z & = \text{Objective Variable} \\
% & \ge 0
% \end{align*}
% </latex>
% Define the necessary objective variable with its constraint.
z = col;
bu(col) = Inf;
bl(col) = 0;
col = col+1;
% Define the flows with their constraints.
for k=K
for e=E
f.(k{1}).(e{1}) = col;
bu(col) = Inf;
bl(col) = 0;
col = col+1;
end
end
%% Define the Objective.
%
% Say Percentage of Injected Flow ( PIF ) describes the ratio of a
% commodity's flow injected into the network over the commodity's demand.
% The objective of MCFP can then be described as maximizing the smallest
% PIF of all the commodities.
%
% <latex>
% \begin{align*}
% z_k & = \text{PIF of Commodity $k\in K$} \\
% & = \frac{\sum_{v\in V} f_{k,(s_k,v)}}{d_k} \\
% \text{Objective } & : \text{ Maximize } \min\{ z_k : k \in K\}
% \end{align*}
% </latex>
% Define the objective constraints.
for k=K
s = K_s.(k{1});
d = K_d.(k{1});
for v=V
for e=E
if ( cs(E_s.(e{1}),s) && cs(E_t.(e{1}),v{1}) )
Ain(rowin,f.(k{1}).(e{1})) = -1/d;
end
end
end
Ain(rowin,z) = 1;
rowin = rowin+1;
end
% Define the objective function.
o = zeros(1,col-1);
o(z) = -1;
%% Define the Capacity Constraint
%
% The Capacity Constraint can be summarized as the sum of flow at any edge
% must be less than or equal to the capacity of the edge.
%
% <latex>
% \begin{align*}
% c \ge \sum_{k \in K} f_{k,e}, \forall e \in E
% \end{align*}
% </latex>
% Define the capacity constraint.
for e=E
for k=K
Ain(rowin,f.(k{1}).(e{1})) = 1;
end
bin(rowin) = E_c.(e{1});
rowin = rowin+1;
end
%% Define the Conservation Constraints.
%
% Say net flow describes for a commodity the sum of flow going into and out
% of a node. The conservation constraint can then be defined as net flow
% should equal zero if the node isn't any commodity's source or target.
% But, if the node is a source, the netflow should equal the demand. If a
% target, netflow should equal negative demand.
%
% <latex>
% \begin{align*}
% \Delta_{k,v} &= \text{Net Flow for Commodity $k$ of Edge $e$}\\
% & = \sum_{w\in V} f_{k,(v,w)} - \sum_{w\in V} f_{k,(w,v)}\\
% & = 0, \forall k \in K, \forall v \in \{w \in V: w\neq s_k, t_k\} \\
% & = d_k, \forall k \in K, \forall v \in \{w \in V: w = s_k\} \\
% & = -d_k, \forall k \in K, \forall v \in \{w \in V: w = t_k\} \\
% \end{align*}
% </latex>
% Define the conservation constraint.
for k=K
s = K_s.(k{1});
t = K_t.(k{1});
d = K_d.(k{1});
for v=V
for w=V
for e=E
if ( cs(E_s.(e{1}),v{1}) && cs(E_t.(e{1}),w{1}) )
Aeq(roweq,f.(k{1}).(e{1})) = 1;
end
if ( cs(E_s.(e{1}),w{1}) && cs(E_t.(e{1}),v{1}) )
Aeq(roweq,f.(k{1}).(e{1})) = -1;
end
end
end
if ( ~cs(v{1},s) && ~cs(v{1},t) )
beq(roweq) = 0;
elseif ( v{1}==s )
beq(roweq) = d;
elseif ( v{1}==t )
beq(roweq) = -d;
end
roweq = roweq+1;
end
end
%% Optimized with MATLAB's linprog.
%
% Finally, the LP problem is solved with MATLAB's linprog.
% Apply MATLAB's linprog function to solve the LP problem.
[ x, fval, exitflag, output ] = ...
linprog( o, Ain, bin, Aeq, beq, bl, bu );
%% Display the results.
% Build table with results.
Variables = {};
Variables{z} = 'z';
for k=K
for e=E
Variables{f.(k{1}).(e{1})} = ['f_',k{1},e{1}];
end
end
T = table;
T.Variables = Variables';
T.Values = x;
% Display the results.
output
T
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment