%% Solve the Time Division Multiplexing variant of the Multi Commodity Flow Problem in MATLAB. % % The goal of this script is to solve the Multi Commodity Flow Problem % modified for Time Division Multiplexing ( TDM ) as seen in Contention % Free Routing ( CFR ). The Multi Commodity Flow Problem is solved for the % Concurrent Flow variant ( MCFP-CF ) with mixed-integer linear programming % ( MILP ). % % The documentation generated from this script is a tex file. The PDF can % be published 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. % % % \begin{align*} % G & = \text{Directed Graph} \\ % & = (V,E) \\ % V & = \text{Node Set} \\ % E & = \text{Arc Set} \\ % e & = (v_e\in V,w_e\in V) \\ % & \in E \\ % v_e & = \text{Source Node of Arc $e\in E$} \\ % w_e & = \text{Sink Node of Arc $e\in E$} \\ % \end{align*} % % Define the Nodes. V = {'A','B','C','D'}; % Define the Arcs. tuples = { { 'eAB' 'A' 'B' }, ... { 'eBC' 'B' 'C' }, ... { 'eCA' 'C' 'A' } }; E = {}; for tuple=tuples e = tuple{1}{1}; E{end+1} = e; E_s.(e) = tuple{1}{2}; E_t.(e) = tuple{1}{3}; end % Plot the graph. figure; G = digraph(struct2cell(E_s),struct2cell(E_t)); plot(G); title( 'Digraph G' ); %% Define the Commodities. % % The commodities are defined as follows. % % % \begin{align*} % K & = \text{Commodity Set} \\ % k & = (s_k\in V,t_k\in V,d_k\in \mathbf{Z}_{\ge 0} : s_k \neq t_k) \\ % & \in K \\ % s_k & = \text{Source Node of Commodity $k\in K$} \\ % t_k & = \text{Sink Node of Commodity $k\in K$} \\ % d_k & = \text{Demand of Commodity $k\in K$} \\ % \end{align*} % % Define the commodities. tuples = { { 'k1' 'A' 'C' 1 },... { 'k2' 'B' 'A' 2 } }; 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 %% Define the Slots. % % An important distinction between TDM MCFP-CF over regular MCFP-CF is that % the arcs are divided into several slots. % % % \begin{align*} % L & = \text{Slot Sequence} \\ % & = (1,2,...,c) \\ % c & = \text{Capacity of each Arc $e\in E$ } \\ % & = \sum_{k\in K} d_k % \end{align*} % % Define the Slots. c = sum(struct2array(K_d)); L = 1:c; %% Define variables needed for script. % Variable related declaration. col = 1; intcont = []; % The following anonymous function makes string comparison easier. cs = @(s1,s2)strcmp(s1,s2); %% Define the Flows and Objective Function. % % The objective function has the following form. The goal is to maximize % the smallest supply-demand ratio of all the commodities. % % % \begin{align*} % \text{Maximize} & : \min\left\{ \frac{\sum_{v\in V,l \in L}f_{k,(s_k,v),l}}{d_k}: k \in K \right\} \\ % f_{k,e,l} & = \text{Flow of Commodity $k\in K$ on Edge $e\in E$ for Slot $l\in L$} \\ % & \in \{0,1\} \\ % \end{align*} % % Define the flows with their constraints. for k=K for e=E for l=L f.(k{1}).(e{1})(l) = col; bu(col) = 1; bl(col) = 0; intcont(end+1) = col; col = col+1; end end end % Define the variable related to the objective. z = col; bu(col) = Inf; bl(col) = 0; col = col+1; % Constraints related declaration. rowin = 1; roweq = 1; Ain = zeros(1,col-1); bin = zeros(1,1); Aeq = zeros(1,col-1); beq = zeros(1,1); % Define the objective function. for k=K s = K_s.(k{1}); d = K_d.(k{1}); for v=V for e=E if ( cs(s,E_s.(e{1})) && cs(v{1},E_t.(e{1})) ) for l=L Ain(rowin,f.(k{1}).(e{1})(l)) = -1/d; end end end end Ain(rowin,z) = 1; bin(rowin) = 0; rowin = rowin+1; end o = zeros(1,col-1); o(z) = -1; %% Define the Slot Capacity Constraint. % % A slot can only be filled by a single commodity % % % \begin{align*} % 1 & \ge \sum_{k\in K} f_{k,e,l} : \forall e\in E,\forall l \in L \\ % \end{align*} % % Define the Slot Capacity Constraint. for e=E for l=L for k=K Ain(rowin,f.(k{1}).(e{1})(l)) = 1; end bin(rowin) = 1; rowin = rowin+1; end end %% Define the TDM Constraint. % % As defined in CFR, each node can buffer a single flow unit, thus the flow % of a commodity should shift to an adjacent slot between edges who share the % same node. However, this constraint should be broken for the commodity's % source and sink nodes. % % % \begin{align*} % 0 & = \sum_{w \in V} f_{k,(w,v),l} - \sum_{w \in V} f_{k,(v,w),(l\text{ mod}|L|)+1} \\ % & : \forall k\in K,\forall v\in \{w\in V: w\neq s_k,t_k \},\forall l \in L \\ % \end{align*} % % Define the TDM Constraint. for k=K k_s = K_s.(k{1}); k_t = K_t.(k{1}); for v=V if ( cs(v{1},k_s) || cs(v{1},k_t) ) continue; end for l=L for w=V for e=E if ( cs(w{1},E_s.(e{1})) && cs(v{1},E_t.(e{1})) ) Aeq(roweq,f.(k{1}).(e{1})(l)) = 1; end if ( cs(v{1},E_s.(e{1})) && cs(w{1},E_t.(e{1})) ) Aeq(roweq,f.(k{1}).(e{1})(mod(l,c)+1)) = -1; end end end beq(roweq) = 0; roweq = roweq+1; end end end %% Define the Supply Equals Demand Constraint. % % The supply of a commodity's source node should equal the demand of the % commodity's sink node. % % % \begin{align*} % 0 & = \sum_{v\in V,l\in L} f_{k,(s_k,v),l} - \sum_{v\in V,l\in L} f_{k,(v,t_k),l} \\ % & : \forall k\in K \\ % \end{align*} % % Define the Supply Equals Demand Constraint. for k=K k_s = K_s.(k{1}); k_t = K_t.(k{1}); for e=E e_s = E_s.(e{1}); e_t = E_t.(e{1}); if ( cs(e_s,k_s) && cs(e_t,k_t) ) Aeq(roweq,f.(k{1}).(e{1})(L(:))) = zeros(1,c); elseif ( cs(e_s,k_s) ) Aeq(roweq,f.(k{1}).(e{1})(L(:))) = ones(1,c); elseif ( cs(e_t,k_t) ) Aeq(roweq,f.(k{1}).(e{1})(L(:))) = -ones(1,c); end end beq(roweq) = 0; roweq = roweq+1; end %% Define the Demand Always Met Constraint. % % Clearly, the supply of each commodity must equal or exceed the % commodity's demand. % % % \begin{align*} % d_k & \le \sum_{w\in V,l \in L} f_{k,(s_k,w),l} \\ % & : \forall k \in K % \end{align*} % % Define the Demand Always Met Constraint. for k=K k_s = K_s.(k{1}); k_d = K_d.(k{1}); for e=E e_s = E_s.(e{1}); if ( cs(e_s,k_s) ) Ain(rowin,f.(k{1}).(e{1})(L(:))) = -ones(1,c); end end bin(rowin) = -k_d; rowin = rowin+1; end %% Solve with MATLAB's intlinprog. % % Optimize for the flows and objective variable with MATLAB's function for % performing MILP. % Solve for the variables with MILP. [x,fval,exitflag,output] = intlinprog(o,intcont,Ain,bin,Aeq,beq,bl,bu); %% Display the Table Results. % Build table with results. Variables = {}; Variables{z} = 'z'; for k=K for e=E for l=L Variables{f.(k{1}).(e{1})(l)} = ['f_',k{1},e{1},num2str(l)]; end end end T = table; T.Variables = Variables'; T.Values = x; % Display the results. output T %% Generate the Visualized Graph results. % Declarations. fn = 'temp.gif'; % Add an index value to the commodity tuple. k_q = 1; for k=K K_q.(k{1}) = k_q; k_q = k_q+1; end % Generate figure for each slot. fsi = []; for l=L % Create figure and grab necessary values. fh = figure; ph = plot(G); pxlim = fh.CurrentAxes.XAxis.Limits; pylim = fh.CurrentAxes.YAxis.Limits; title(['Digraph G for Time Slot ',num2str(l)]); % Highlight the flow of each commodity during the current time slot. for k=K k_q = K_q.(k{1}); hls = {}; for e=E if (x(f.(k{1}).(e{1})(l))) hls{end+1} = E_s.(e{1}); hls{end+1} = E_t.(e{1}); end end col = [0.5,k_q/numel(K),0]; highlight(ph,hls,'EdgeColor',col,'LineWidth',2); % Add legend. text( ... pxlim(1), ... pylim(2)*.9-((pylim(2)-pylim(1))*(k_q-1)/numel(K)), ... ['---' k{1}],'Color',col); end axis square; % Update figure; drawnow if isempty(fsi) fsi = fh.Position; else fh.Position = fsi; refresh; drawnow; end % Build gid for animation. im = frame2im(getframe(fh)); [imind,cm] = rgb2ind(im,256); if ( l==1 ) imwrite(imind,cm,fn,'gif','Loopcount',inf); else imwrite(imind,cm,fn,'gif','WriteMode','append'); end end