- Gamma distributon random variable.
- Test file for random.h
3 * This file is a part of LEMON, a generic C++ optimization library
5 * Copyright (C) 2003-2007
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
19 #ifndef LEMON_SSP_MIN_COST_FLOW_H
20 #define LEMON_SSP_MIN_COST_FLOW_H
22 ///\ingroup min_cost_flow
25 ///\brief An algorithm for finding a flow of value \c k (for
26 ///small values of \c k) having minimal total cost
29 #include <lemon/dijkstra.h>
30 #include <lemon/graph_adaptor.h>
31 #include <lemon/maps.h>
36 /// \addtogroup min_cost_flow
39 /// \brief Implementation of an algorithm for finding a flow of
40 /// value \c k (for small values of \c k) having minimal total cost
44 /// The \ref lemon::SspMinCostFlow "Successive Shortest Path Minimum
45 /// Cost Flow" implements an algorithm for finding a flow of value
46 /// \c k having minimal total cost from a given source node to a
47 /// given target node in a directed graph with a cost function on
48 /// the edges. To this end, the edge-capacities and edge-costs have
49 /// to be nonnegative. The edge-capacities should be integers, but
50 /// the edge-costs can be integers, reals or of other comparable
51 /// numeric type. This algorithm is intended to be used only for
52 /// small values of \c k, since it is only polynomial in k, not in
53 /// the length of k (which is log k): in order to find the minimum
54 /// cost flow of value \c k it finds the minimum cost flow of value
55 /// \c i for every \c i between 0 and \c k.
57 ///\param Graph The directed graph type the algorithm runs on.
58 ///\param LengthMap The type of the length map.
59 ///\param CapacityMap The capacity map type.
61 ///\author Attila Bernath
62 template <typename Graph, typename LengthMap, typename CapacityMap>
63 class SspMinCostFlow {
65 typedef typename LengthMap::Value Length;
67 //Warning: this should be integer type
68 typedef typename CapacityMap::Value Capacity;
70 typedef typename Graph::Node Node;
71 typedef typename Graph::NodeIt NodeIt;
72 typedef typename Graph::Edge Edge;
73 typedef typename Graph::OutEdgeIt OutEdgeIt;
74 typedef typename Graph::template EdgeMap<int> EdgeIntMap;
76 typedef ResGraphAdaptor<const Graph,int,CapacityMap,EdgeIntMap> ResGW;
77 typedef typename ResGW::Edge ResGraphEdge;
82 const LengthMap& length;
83 const CapacityMap& capacity;
86 typedef typename Graph::template NodeMap<Length> PotentialMap;
87 PotentialMap potential;
95 typedef typename Graph::template NodeMap<Length> NodeMap;
97 const LengthMap &length;
100 typedef typename LengthMap::Key Key;
101 typedef typename LengthMap::Value Value;
103 ModLengthMap(const ResGW& _g,
104 const LengthMap &_length, const NodeMap &_pot) :
105 g(_g), /*rev(_rev),*/ length(_length), pot(_pot) { }
107 Value operator[](typename ResGW::Edge e) const {
109 return length[e]-(pot[g.target(e)]-pot[g.source(e)]);
111 return -length[e]-(pot[g.target(e)]-pot[g.source(e)]);
117 ModLengthMap mod_length;
118 Dijkstra<ResGW, ModLengthMap> dijkstra;
122 /// \brief The constructor of the class.
124 /// \param _g The directed graph the algorithm runs on.
125 /// \param _length The length (cost) of the edges.
126 /// \param _cap The capacity of the edges.
127 /// \param _s Source node.
128 /// \param _t Target node.
129 SspMinCostFlow(Graph& _g, LengthMap& _length, CapacityMap& _cap,
131 g(_g), length(_length), capacity(_cap), flow(_g), potential(_g),
133 res_graph(g, capacity, flow),
134 mod_length(res_graph, length, potential),
135 dijkstra(res_graph, mod_length) {
139 /// \brief Tries to augment the flow between s and t by 1. The
140 /// return value shows if the augmentation is successful.
143 if (!dijkstra.reached(t)) {
145 //Unsuccessful augmentation.
149 //We have to change the potential
150 for(typename ResGW::NodeIt n(res_graph); n!=INVALID; ++n)
151 potential.set(n, potential[n]+dijkstra.distMap()[n]);
153 //Augmenting on the shortest path
157 e = dijkstra.predEdge(n);
158 n = dijkstra.predNode(n);
159 res_graph.augment(e,1);
160 //Let's update the total length
161 if (res_graph.forward(e))
162 total_length += length[e];
164 total_length -= length[e];
171 /// \brief Runs the algorithm.
173 /// Runs the algorithm.
174 /// Returns k if there is a flow of value at least k from s to t.
175 /// Otherwise it returns the maximum value of a flow from s to t.
177 /// \param k The value of the flow we are looking for.
179 /// \todo May be it does make sense to be able to start with a
180 /// nonzero feasible primal-dual solution pair as well.
182 /// \todo If the current flow value is bigger than k, then
183 /// everything is cleared and the algorithm starts from zero
184 /// flow. Is it a good approach?
186 if (flowValue()>k) reset();
187 while (flowValue()<k && augment()) { }
191 /// \brief The class is reset to zero flow and potential.
194 for (typename Graph::EdgeIt e(g); e!=INVALID; ++e) flow.set(e, 0);
195 for (typename Graph::NodeIt n(g); n!=INVALID; ++n) potential.set(n, 0);
198 /// \brief Returns the value of the current flow.
199 int flowValue() const {
201 for (typename Graph::OutEdgeIt e(g, s); e!=INVALID; ++e) i+=flow[e];
202 for (typename Graph::InEdgeIt e(g, s); e!=INVALID; ++e) i-=flow[e];
206 /// \brief Total cost of the found flow.
208 /// This function gives back the total cost of the found flow.
209 Length totalLength(){
213 /// \brief Returns a const reference to the EdgeMap \c flow.
215 /// Returns a const reference to the EdgeMap \c flow.
216 const EdgeIntMap &getFlow() const { return flow;}
218 /// \brief Returns a const reference to the NodeMap \c potential
219 /// (the dual solution).
221 /// Returns a const reference to the NodeMap \c potential (the
223 const PotentialMap &getPotential() const { return potential;}
225 /// \brief Checking the complementary slackness optimality criteria.
227 /// This function checks, whether the given flow and potential
228 /// satisfy the complementary slackness conditions (i.e. these are optimal).
229 /// This function only checks optimality, doesn't bother with feasibility.
230 /// For testing purpose.
231 bool checkComplementarySlackness(){
234 for(typename Graph::EdgeIt e(g); e!=INVALID; ++e) {
236 mod_pot = length[e]-potential[g.target(e)]+potential[g.source(e)];
238 if (0<fl_e && fl_e<capacity[e]) {
239 /// \todo better comparison is needed for real types, moreover,
240 /// this comparison here is superfluous.
245 if (mod_pot > 0 && fl_e != 0)
247 if (mod_pot < 0 && fl_e != capacity[e])
254 }; //class SspMinCostFlow
260 #endif //LEMON_SSP_MIN_COST_FLOW_H