2 %> Algorithm 21.1: Generation of a simple cycle flow from a circulation. Implementation of algorithm 21.1 of \cite Bier15-book
4 %> @note Tested with \ref runSimpleCycle.m
8 %> @author Michel Bierlaire
9 %> @date Fri Mar 27 17:09:10 2015
10 %> @ingroup Algorithms
13 %> Extract a simple cycle flow from a circulation
14 %> @param adj the adjacency matrix of the network. It is a \f$m \times m\f$ matrix, such that the element at row i and column j corresponds to the
id of the arc (i,j). The numbering should be from 1 to n.
15 %> @param circ the circulation flow. It is a vector of \f$\mathbb{R}^n\f$.
16 %> @param printlevel If different from 0, the cycle is printed. (Default: 0)
17 %> @
return cycleflow a simple cycle flow. It is a vector of \f$\mathbb{R}^n\f$.
18 function cycleflow =
simpleCycle(adj,circ,printlevel=0)
19 cycleflow = zeros(size(circ)) ;
21 if (columns(adj) != nnodes)
22 error("Adjacency matrix must be square") ;
25 nonzeros = find(circ) ;
26 if(length(nonzeros)==0)
27 printf("The circ flow vector is zero")
30 # Find an arc carrying positive flow 31 theArc = nonzeros(1) ;
32 theLoopArcs = [theArc] ;
36 layer{layerNumber} = arcs(theArc,2) ;
37 otherNode = arcs(theArc,1) ;
38 theLoopOrientation = [1] ;
40 layer{layerNumber} = arcs(theArc,1) ;
41 otherNode = arcs(theArc,2) ;
42 theLoopOrientation = [-1] ;
45 treatedNodes = zeros(nnodes,1) ;
46 treatedNodes(layer{layerNumber}) = 1 ;
49 currentlayer = zeros(0) ;
50 for i = layer{layerNumber}
53 if (circ(arcij) > 0) && (treatedNodes(arcs(arcij,2)) == 0)
54 % printf("Forward arc: %d -> %d\n",arcs(arcij,1),arcs(arcij,2));
55 currentlayer(end + 1) = arcs(arcij,2) ;
57 treatedNodes(arcs(arcij,2)) = 1 ;
63 if (circ(arcij) < 0) && (treatedNodes(arcs(arcij,1)) == 0)
64 % printf("Backward arc: %d -> %d\n",arcs(arcij,1),arcs(arcij,2));
65 currentlayer(end + 1) = arcs(arcij,1) ;
67 treatedNodes(arcs(arcij,1)) = 1 ;
73 layer(layerNumber) = currentlayer;
76 % Here, layer contains the nodes organized into layers
77 % Search for the node to close the loop
79 for key = 1:length(layer)
80 theLayer = layer{key} ;
81 % printf(
"Look for %d in \n",otherNode)
82 idx = find(theLayer == otherNode) ;
84 theLayerOfOtherNode = key ;
85 theIndexOfOtherNode = idx ;
89 theLoop = [otherNode] ;
91 for k = theLayerOfOtherNode-1:-1:1
92 % Find an arc carrying positive flow. It must exist by construction
94 theForwardArc = adj(ell,theLoop(end)) ;
95 if (theForwardArc != 0)
96 forwardflow = circ(theForwardArc) ;
98 theLoop(end+1) = ell ;
99 theLoopArcs(end+1) = theForwardArc ;
100 theLoopOrientation(end+1) = 1 ;
101 if (forwardflow < flowMin)
102 flowMin = forwardflow ;
107 theBackwardArc = adj(theLoop(end),ell) ;
108 if (theBackwardArc != 0)
109 backwardflow = circ(theBackwardArc) ;
110 if (backwardflow < 0)
111 theLoop(end+1) = ell ;
112 theLoopArcs(end+1) = theBackwardArc ;
113 theLoopOrientation(end+1) = -1 ;
114 if (-backwardflow < flowMin)
115 flowMin = -backwardflow;
122 for i = length(theLoopArcs):-1:1
124 printf(
"%d",theLoop(i))
125 if (theLoopOrientation(i) == 1)
131 cycleflow(theLoopArcs(i)) = flowMin * theLoopOrientation(i) ;
134 printf(
"%d [%f]\n",theLoop(end),flowMin)
function circulationDecomposition(in adj, in circ)
Decompose a circulation into simple cycle flows.
function simpleCycle(in adj, in circ, in printlevel)
Extract a simple cycle flow from a circulation.
function prepareNetwork(in adj)
Identifies the upstream and downstream nodes of each arc.