summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/Complete.hpp8
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp12
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/Log.hpp6
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/MaxStrategy.hpp82
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/Operator.hpp82
-rw-r--r--clang/include/clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp34
-rw-r--r--clang/lib/Analysis/Interval.cpp496
-rwxr-xr-xclang/lib/Analysis/MCFClass.h2436
-rw-r--r--clang/lib/Analysis/MCFSimplex.cpp68
-rw-r--r--clang/lib/Analysis/MCFSimplex.h1086
-rwxr-xr-xclang/lib/Analysis/OPTUtils.h475
-rw-r--r--impl/Operator.hpp4
12 files changed, 326 insertions, 4463 deletions
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/Complete.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/Complete.hpp
index 664d71f..73c8f23 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/Complete.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/Complete.hpp
@@ -4,9 +4,15 @@
#include <cassert>
#include <ostream>
#include <istream>
+#include <limits>
template<typename T>
-T infinity() { }
+T infinity();
+
+template<>
+double infinity() {
+ return std::numeric_limits<double>::infinity();
+}
template<typename T>
struct Complete {
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp
index 5ee5405..3342cc7 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp
@@ -76,14 +76,9 @@ struct EquationSystem {
}
Constant<Domain>& constant(const Domain& value) {
- if (_constants.find(value) == _constants.end()) {
- Constant<Domain>* constant = new Constant<Domain>(value);
- _expressions.insert(constant);
- _constants[value] = constant;
- return *constant;
- } else {
- return *_constants[value];
- }
+ Constant<Domain>* constant = new Constant<Domain>(value);
+ _expressions.insert(constant);
+ return *constant;
}
MaxExpression<Domain>* operator[](const Variable<Domain>& var) const {
@@ -132,7 +127,6 @@ struct EquationSystem {
private:
std::set<Operator<Domain>*> _operators;
std::set<Expression<Domain>*> _expressions;
- std::map<Domain,Constant<Domain>*> _constants;
std::vector<Variable<Domain>*> _variables;
std::map<std::string, Variable<Domain>*> _variable_names;
IdMap<MaxExpression<Domain>, Variable<Domain>*>* _expr_to_var;
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/Log.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/Log.hpp
index 90ab7e5..7502ac8 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/Log.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/Log.hpp
@@ -1,17 +1,13 @@
#ifndef LOG_HPP
#define LOG_HPP
-// could not be hackier, but C++ is annoying
-#define protected public
-#include <streambuf>
-#undef protected
-
#include <string>
#include <iostream>
#include <map>
#include <cstdio>
namespace solver_log {
+
struct Logger : public std::ostream {
Logger(std::streambuf* buffer, const std::string& name)
: std::ostream(NULL) { }
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/MaxStrategy.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/MaxStrategy.hpp
index 9ae7394..da05dd8 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/MaxStrategy.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/MaxStrategy.hpp
@@ -39,22 +39,13 @@ struct DynamicMaxStrategy : public MaxStrategy<Domain> {
_values(system.maxExpressionCount(), 0),
_stable(system.maxExpressionCount()),
_influence(system.maxExpressionCount(),
- IdSet<MaxExpression<Domain> >(system.maxExpressionCount())),
- _changed(false)
+ IdSet<MaxExpression<Domain> >(system.maxExpressionCount()))
{}
void setRho(DynamicVariableAssignment<Domain>& rho) {
_rho = &rho;
}
- void hasChanged(bool c) {
- _changed = c;
- }
-
- bool hasChanged() const {
- return _changed;
- }
-
unsigned int get(const MaxExpression<Domain>& e) const {
solver_log::strategy << indent() << "Look up " << e << std::endl;
return _values[e];
@@ -66,19 +57,6 @@ struct DynamicMaxStrategy : public MaxStrategy<Domain> {
return _values[e];
}
- void invalidate(const Variable<Domain>& v) {
- solver_log::strategy << indent() << "Invalidating " << v << " - " << *_system[v] << std::endl;
-
- IdSet<MaxExpression<Domain> > infl = _influence[*_system[v]];
- for (typename IdSet<MaxExpression<Domain> >::iterator
- it = infl.begin(),
- end = infl.end();
- it != end;
- ++it) {
- invalidate(_system.maxExpression(*it));
- }
- }
-
void invalidate(const MaxExpression<Domain>& v) {
if (_stable.contains(v)) {
solver_log::strategy << indent() << "Invalidating " << v << std::endl;
@@ -111,20 +89,18 @@ private:
if (val != _values[x]) {
solver_log::strategy << x << " => " << *x.arguments()[val] << std::endl;
- IdSet<MaxExpression<Domain> > oldInfluence = _influence[x];
- //_influence[x].clear();
_values[x] = val;
- _changed = true;
_rho->invalidate(*_system.varFromExpr(x));
//_rho->stabilise();
- _stable.filter(oldInfluence);
-
+ IdSet<MaxExpression<Domain> > infl = _influence[x];
+ _stable.filter(infl);
+
for (typename IdSet<MaxExpression<Domain> >::iterator
- it = oldInfluence.begin(),
- end = oldInfluence.end();
+ it = infl.begin(),
+ end = infl.end();
it != end;
++it) {
solve(_system.maxExpression(*it));
@@ -141,7 +117,7 @@ private:
struct DependencyAssignment : public VariableAssignment<Domain>{
DependencyAssignment(DynamicMaxStrategy& strat,
- VariableAssignment<Domain>& rho,
+ DynamicVariableAssignment<Domain>& rho,
const MaxExpression<Domain>& expr)
: _strat(strat),
_rho(rho),
@@ -149,14 +125,25 @@ private:
_current(strat._system.variableCount()) {
}
const Domain operator[](const Variable<Domain>& var) {
- // solve the strategy for this variable, too
- _strat.solve(*_strat._system[var]);
+ // force evaluation to get touched variables
+ Domain value = _rho[var];
+
_strat._influence[*_strat._system[var]].insert(_expr);
- return _rho[var];
+
+ // invalidate touched variables
+ IdSet<Variable<Domain> > changed = _rho.get_changed();
+ for (typename IdSet<Variable<Domain> >::iterator
+ it = changed.begin(),
+ ei = changed.end();
+ it != ei;
+ ++it) {
+ _strat.invalidate(*_strat._system[_strat._system.variable(*it)]);
+ }
+ return value;
}
private:
DynamicMaxStrategy& _strat;
- VariableAssignment<Domain>& _rho;
+ DynamicVariableAssignment<Domain>& _rho;
const MaxExpression<Domain>& _expr;
IdSet<Variable<Domain> > _current;
};
@@ -186,7 +173,6 @@ private:
IdMap<MaxExpression<Domain>,unsigned int> _values;
IdSet<MaxExpression<Domain> > _stable;
IdMap<MaxExpression<Domain>,IdSet<MaxExpression<Domain> > > _influence;
- bool _changed;
};
@@ -200,28 +186,12 @@ struct Solver {
}
Domain solve(Variable<Domain>& var) {
- MaxExpression<Domain>& rhs = *_system[var];
-
- // _rho automatically keeps up to date with changes made to the
- // strategy, so you don't need to stabilise it
-
+ MaxExpression<Domain>& rhs = *_system[var];
+ // this will automatically work sufficiently to get the final
+ // strategy for this variable's RHS
_strategy.get(rhs);
-
- /*
- do {
- _strategy.hasChanged(false);
-
- solver_log::debug << "Stabilise assignment (toplevel)" << std::endl;
- _rho.stabilise();
-
- solver_log::debug << "Improve strategy (toplevel)" << std::endl;
- // improve strategy
- _strategy.get(rhs);
- solver_log::debug << (_strategy.hasChanged() ? "Strategy has changed - loop again" : "Strategy has not changed - terminate") << std::endl;
- } while (_strategy.hasChanged());
- */
-
+ // this will automatically solve for the value of the RHS, if required
return _rho[var];
}
private:
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/Operator.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/Operator.hpp
index 64ef096..0b08866 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/Operator.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/Operator.hpp
@@ -135,6 +135,88 @@ struct Guard : public Operator<Domain> {
}
};
+#include "MCFSimplex.h"
+
+template<typename Domain>
+struct MinCostFlow : public Operator<Domain> {
+ MinCostFlow(const std::vector<Domain>& supplies, const std::vector<std::pair<int,int> > arcs)
+ : _supplies(supplies), _arcs(arcs), _solver(0,0) {
+ MCFClass::FNumber* deficits = new MCFClass::FNumber[supplies.size()];
+ MCFClass::Index* starts = new MCFClass::Index[arcs.size()];
+ MCFClass::Index* ends = new MCFClass::Index[arcs.size()];
+
+ for (int i = 0, size = supplies.size(); i < size; ++i) {
+ deficits[i] = -supplies[i].template as<MCFClass::FNumber>();
+ }
+ for (int i = 0, size = arcs.size(); i < size; ++i) {
+ starts[i] = arcs[i].first;
+ ends[i] = arcs[i].second;
+ }
+
+ _solver.LoadNet(supplies.size(), arcs.size(), // max nodes/arcs
+ supplies.size(), arcs.size(), // current nodes/arcs
+ NULL, NULL, // arcs have inf cap, zero cost (to begin)
+ deficits, // deficits for each node
+ starts, ends); // start/end of each edge
+
+ delete[] deficits;
+ delete[] starts;
+ delete[] ends;
+ }
+ Domain eval (const std::vector<Domain>& costs) const {
+ assert(costs.size() == _arcs.size());
+ if (costs.size() == 0)
+ return 0;
+ for (int i = 0, size = costs.size(); i < size; ++i) {
+ _solver.ChgCost(i, costs[i].template as<MCFClass::CNumber>());
+ }
+ _solver.SolveMCF();
+ if (_solver.MCFGetStatus() == MCFClass::kUnfeasible){
+ return -infinity<Domain>();
+ } else if (_solver.MCFGetFO() == MCFClass::Inf<MCFClass::FONumber>()) {
+ return infinity<Domain>();
+ } else if (_solver.MCFGetFO() == -MCFClass::Inf<MCFClass::FONumber>()) {
+ return -infinity<Domain>();
+ } else {
+ return _solver.MCFGetFO();
+ }
+ }
+ void print(std::ostream& cout) const {
+ std::string supplyString = "[";
+ for (int i = 0, size = _supplies.size(); i < size; ++i) {
+ if (i > 0)
+ supplyString += ",";
+ std::stringstream stream;
+ stream << _supplies[i];
+ supplyString += stream.str();
+ }
+ supplyString += ']';
+
+ std::string arcString = "[";
+ for (int i = 0, size = _arcs.size(); i < size; ++i) {
+ if (i > 0)
+ arcString += ",";
+ {
+ std::stringstream stream;
+ stream << _arcs[i].first;
+ arcString += stream.str() + ":";
+ }
+ {
+ std::stringstream stream;
+ stream << _arcs[i].second;
+ arcString += stream.str();
+ }
+ }
+ arcString += ']';
+
+ cout << "MCF<" << supplyString << "," << arcString << ">";
+ }
+private:
+ std::vector<Domain> _supplies;
+ std::vector<std::pair<int,int> > _arcs;
+ mutable MCFSimplex _solver;
+};
+
/*#include "TemplateConstraintMatrix.hpp"
template<typename Domain>
diff --git a/clang/include/clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp b/clang/include/clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp
index 746e5ef..4c4519a 100644
--- a/clang/include/clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp
+++ b/clang/include/clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp
@@ -24,9 +24,11 @@ struct DynamicVariableAssignment : public VariableAssignment<Domain> {
) : _system(system),
_strategy(strat),
_values(system.variableCount(), value),
+ _old_values(system.variableCount(), value),
_unstable(system.variableCount()),
_influence(system.variableCount(),
- IdSet<Variable<Domain> >(system.variableCount()))
+ IdSet<Variable<Domain> >(system.variableCount())),
+ _touched(system.variableCount())
{ }
const Domain operator[](const Variable<Domain>& var) {
@@ -34,18 +36,13 @@ struct DynamicVariableAssignment : public VariableAssignment<Domain> {
return _values[var];
}
- /*void stabilise() {
- if (!_unstable.empty()) {
- Variable<Domain>& var = _system.variable(*_unstable.begin());
- solve(var);
- }
- }*/
-
void invalidate(const Variable<Domain>& x) {
if (!_unstable.contains(x)) {
solver_log::fixpoint << indent() << "Invalidating " << x << std::endl;
_unstable.insert(x);
+ _touched.insert(x);
+ _old_values[x] = _values[x];
_values[x] = infinity<Domain>();
IdSet<Variable<Domain> > infl = _influence[x];
@@ -61,6 +58,23 @@ struct DynamicVariableAssignment : public VariableAssignment<Domain> {
}
}
+ IdSet<Variable<Domain> > get_changed() {
+ IdSet<Variable<Domain> > changed;
+ for (typename IdSet<Variable<Domain> >::iterator
+ it = _touched.begin(),
+ ei = _touched.end();
+ it != ei;
+ ++it) {
+ Variable<Domain>& var = _system.variable(*it);
+ if (!_unstable.contains(var) && _old_values[var] != _values[var]) {
+ changed.insert(var);
+ _touched.remove(var);
+ }
+ }
+ //_touched.clear();
+ return changed;
+ }
+
private:
void solve(const Variable<Domain>& x) {
@@ -79,8 +93,6 @@ private:
if (val != _values[x]) {
solver_log::fixpoint << x << " = " << val << std::endl;
- _strategy.invalidate(x);
-
IdSet<Variable<Domain> > oldInfluence = _influence[x];
_influence[x].clear();
_values[x] = val;
@@ -118,8 +130,10 @@ private:
const EquationSystem<Domain>& _system;
DynamicMaxStrategy<Domain>& _strategy;
IdMap<Variable<Domain>, Domain> _values;
+ IdMap<Variable<Domain>, Domain> _old_values;
IdSet<Variable<Domain> > _unstable;
IdMap<Variable<Domain>,IdSet<Variable<Domain> > > _influence;
+ IdSet<Variable<Domain> > _touched;
};
#endif
diff --git a/clang/lib/Analysis/Interval.cpp b/clang/lib/Analysis/Interval.cpp
index e643f5e..4b274ed 100644
--- a/clang/lib/Analysis/Interval.cpp
+++ b/clang/lib/Analysis/Interval.cpp
@@ -9,8 +9,6 @@
#include "clang/Analysis/Analyses/IntervalSolver/VariableAssignment.hpp"
#include "clang/Analysis/Analyses/IntervalSolver/EquationSystem.hpp"
-#include "MCFSimplex.h"
-
#include "llvm/ADT/PostOrderIterator.h"
#include "llvm/ADT/DenseMap.h"
#include "llvm/Support/Process.h"
@@ -22,8 +20,14 @@
#include <set>
using namespace clang;
+using namespace std;
+
+template<typename T>
+T neg(const T& t) {
+ return -t;
+}
-std::string negate_string(const std::string& str) {
+string operator-(const string& str) {
if (str[0] == '-')
return str.substr(1);
return '-' + str;
@@ -31,28 +35,28 @@ std::string negate_string(const std::string& str) {
#include <sstream>
template<typename T>
-std::string toString(const T& obj) {
- std::stringstream stream;
+string toString(const T& obj) {
+ stringstream stream;
stream << obj;
return stream.str();
}
#include <ostream>
template<typename K,typename V>
-std::ostream& operator<<(std::ostream& cout, const std::pair<K,V>& v) {
+ostream& operator<<(ostream& cout, const pair<K,V>& v) {
cout << "(" << v.first << ", " << v.second << ")";
return cout;
}
template<typename V>
-std::ostream& operator<<(std::ostream& cout, const std::pair<Variable<Complete<int64_t> >*, V>& v) {
+ostream& operator<<(ostream& cout, const pair<Variable<Complete<int64_t> >*, V>& v) {
cout << "(" << v.first->name() << ", " << v.second << ")";
return cout;
}
template<typename V>
-std::ostream& operator<<(std::ostream& cout, const std::vector<V>& v) {
+ostream& operator<<(ostream& cout, const vector<V>& v) {
cout << "[";
- for(typename std::vector<V>::const_iterator it = v.begin(), ei = v.end();
+ for(typename vector<V>::const_iterator it = v.begin(), ei = v.end();
it != ei;
++it) {
if (it != v.begin())
@@ -64,9 +68,9 @@ std::ostream& operator<<(std::ostream& cout, const std::vector<V>& v) {
}
template<typename K,typename V>
-std::ostream& operator<<(std::ostream& cout, const std::map<K,V>& v) {
+ostream& operator<<(ostream& cout, const map<K,V>& v) {
cout << "{";
- for (typename std::map<K,V>::const_iterator it = v.begin(), ei = v.end();
+ for (typename map<K,V>::const_iterator it = v.begin(), ei = v.end();
it != ei;
++it) {
if (it != v.begin())
@@ -91,219 +95,89 @@ template<>
ZBar infinity() {
return ZBar(1, true);
}
-//typedef std::map<std::string, ZBar> Vector;
-struct Vector : public std::map<std::string, ZBar> {
+
+
+
+
+struct Vector : public map<string, ZBar> {
Vector(const ZBar& val=0)
: _val(val) { }
- ZBar operator[](const std::string& key) const {
- if (this->find(key) != this->end())
- return this->find(key)->second;
+ ZBar operator[](const string& key) const {
+ const_iterator it = this->find(key);
+ if (it != this->end())
+ return it->second;
return _val;
}
- ZBar& operator[](const std::string& key) {
- if (this->find(key) != this->end())
- return this->find(key)->second;
- std::pair<iterator,bool> p = this->insert(std::pair<const std::string, ZBar>(key, _val));
+ ZBar& operator[](const string& key) {
+ iterator it = this->find(key);
+ if (it != this->end())
+ return it->second;
+ pair<iterator,bool> p = this->insert(pair<const string, ZBar>(key, _val));
return p.first->second;
}
ZBar _val;
};
-Vector negate_vector(const Vector& v) {
+Vector operator-(const Vector& v) {
Vector result;
for (Vector::const_iterator it = v.begin(),
ei = v.end();
it != ei;
++it) {
- result[negate_string(it->first)] = it->second;
+ result[-it->first] = it->second;
}
return result;
}
-typedef std::pair<Vector, ZBar> Result; // a "slice" of an equation
-
-Result negate_result(const Result& r) {
- return Result(negate_vector(r.first), -r.second);
-}
-
-//typedef std::map<std::string, Result> LinearEquation; // one `Result` per variable
-struct LinearEquation : public std::map<std::string, Result> {
- Result operator[](const std::string& key) const {
- if (this->find(key) != this->end())
- return this->find(key)->second;
- Result r;
- r.first[key] = 1;
- r.second = 0;
- return r;
- }
- Result& operator[](const std::string& key) {
- if (this->find(key) != this->end())
- return this->find(key)->second;
- Result r;
- r.first[key] = 1;
- r.second = 0;
- std::pair<iterator,bool> p = this->insert(std::pair<const std::string, Result>(key, r));
- return p.first->second;
- }
-};
-
-typedef Vector Condition;
-
-typedef EquationSystem<ZBar> EqnSys;
-typedef Expression<ZBar> EqnExpr;
-typedef Variable<ZBar> EqnVar;
-
-
-struct LinearOperator : public Operator<Vector> {
- LinearOperator(const LinearEquation* result)
- : _values(result) {}
-
- Vector eval(const std::vector<Vector>& vector) const {
- assert(vector.size() == 1);
- const Vector& v = vector[0];
- Vector result = v;
- for (LinearEquation::const_iterator it = _values->begin(),
- ei = _values->end();
- it != ei;
- ++it) {
- ZBar subresult = 0;
- for (Vector::const_iterator jt = it->second.first.begin(),
- ej = it->second.first.end();
- jt != ej;
- ++jt) {
- subresult += jt->second * v[jt->first];
+Vector operator+(const Vector& left, const Vector& right) {
+ Vector::const_iterator
+ left_iter = left.begin(),
+ left_end = left.end(),
+ right_iter = right.begin(),
+ right_end = right.end();
+
+ Vector result(left._val + right._val);
+ while (left_iter != left_end && right_iter != right_end) {
+ if (left_iter->first == right_iter->first) {
+ result[left_iter->first] = left_iter->second + right_iter->second;
+ left_iter++;
+ right_iter++;
+ } else {
+ if (left_iter->first < right_iter->first) {
+ result[left_iter->first] = left_iter->second;
+ left_iter++;
+ } else {
+ result[right_iter->first] = right_iter->second;
+ right_iter++;
}
- subresult += it->second.second;
- result[it->first] = subresult;
}
- return result;
}
+ Vector::const_iterator it = (right_iter == right_end ? left_iter : right_iter);
+ Vector::const_iterator end = (right_iter == right_end ? left_end : right_end);
+ for (; it != end; ++it)
+ result[it->first] = it->second;
- void print(std::ostream& cout) const {
- cout << "linear[" << *_values << "]";
- }
-
- const LinearEquation* _values;
-};
-
-
-
-template<class F, class M>
-void transform_values(const F& f, M& map) {
- for (typename M::iterator it = map.begin(),
- ei = map.end();
- it != ei;
- ++it) {
- it->second = f(it->second);
- }
-}
-
-template<class M, class F>
-M merge_maps_with(const F& f, const M& left, const M& right) {
- M result;
- typename M::const_iterator first1 = left.begin(), last1 = left.end(),
- first2 = right.begin(), last2 = right.end();
- for (; first1 != last1 && first2 != last2;) {
- if (first2->first < first1->first) {
- result[first2->first] = first2->second;
- ++first2;
- } else if (first1->first == first2->first) {
- result[first1->first] = f(first1->second, first2->second);
- ++first1;
- ++first2;
- } else {
- result[first1->first] = first1->second;
- ++first1;
- }
- }
- while (first1 != last1) {
- result[first1->first] = first1->second;
- ++first1;
- }
- while (first2 != last2) {
- result[first2->first] = first2->second;
- ++first2;
- }
return result;
}
+Vector operator*(const ZBar& left, const Vector& right) {
+ Vector result(left * right._val);
-template<class T>
-T max(const T& l, const T& r) {
- return (l < r ? l : r);
-}
-template<class T>
-T negate(const T& v) {
- return -v;
-}
-template<class T>
-T addValues(const T& l, const T& r) {
- return l + r;
-}
-
-Vector operator-(const Vector& vector) {
- Vector result(-vector._val);
- for (Vector::const_iterator it = vector.begin(),
- ei = vector.end();
- it != ei;
+ for (Vector::const_iterator
+ it = right.begin(),
+ end = right.end();
+ it != end;
++it) {
- result[it->first] = -it->second;
+ result[it->first] = left * it->second;
}
- return result;
-}
-
-Vector operator+(const Vector& left, const Vector& right) {
- return merge_maps_with(addValues<ZBar>, left, right);
-}
-
-Vector operator-(const Vector& left, const Vector& right) {
- return merge_maps_with(addValues<ZBar>, left, -right);
-}
-Vector operator*(const Vector& left, const ZBar& right) {
- Vector result;
- for (Vector::const_iterator it = left.begin(),
- ei = left.end();
- it != ei;
- ++it) {
- result[it->first] = (it->second * right);
- }
return result;
}
-Vector operator*(const ZBar& left, const Vector& right) {
+Vector operator*(const Vector& left, const ZBar& right) {
return right * left;
}
-bool operator<(const Vector& left, const Vector& right) {
- bool equal = true;
- for (Vector::const_iterator it = left.begin(),
- ei = left.end();
- it != ei;
- ++it) {
- if (it->second < right[it->first]) {
- equal = false;
- } else if (it->second > right[it->first]) {
- return false;
- }
- }
- for (Vector::const_iterator it = right.begin(),
- ei = right.end();
- it != ei;
- ++it) {
- if (left[it->first] < it->second) {
- equal = false;
- } else if (left[it->first] > it->second) {
- return false;
- }
- }
- return equal ? left._val < right._val : true;
-}
-
-template<>
-Vector infinity<Vector>() {
- return Vector(infinity<ZBar>());
-}
-std::ostream& operator<<(std::ostream& cout, const Vector& v) {
+ostream& operator<<(ostream& cout, const Vector& v) {
cout << "{";
for (Vector::const_iterator it = v.begin(), ei = v.end();
it != ei;
@@ -317,6 +191,22 @@ std::ostream& operator<<(std::ostream& cout, const Vector& v) {
+typedef pair<Vector, ZBar> Result; // a "slice" of an equation
+Result operator-(const Result& r) {
+ return Result(-r.first, -r.second);
+}
+
+typedef Vector Condition;
+
+typedef EquationSystem<ZBar> EqnSys;
+typedef Expression<ZBar> EqnExpr;
+typedef Variable<ZBar> EqnVar;
+
+
+
+
+
+
@@ -329,7 +219,7 @@ Result fromInteger(const IntegerLiteral* expr) {
}
Result fromDeclExpr(const DeclRefExpr* expr) {
- Vector val;
+ Vector val(0);
val[expr->getNameInfo().getAsString()] = 1;
return Result(val, 0);
}
@@ -356,14 +246,13 @@ Result fromBinary(const BinaryOperator* op) {
case BO_Assign:
return right;
case BO_Sub:
- right = negate_result(right);
+ right = -right;
//transform_values(negate<ZBar>, right.first);
//right.second *= -1;
case BO_Add:
{
Result result;
- result.first = merge_maps_with(addValues<ZBar>,
- left.first, right.first);
+ result.first = left.first + right.first;
result.second = left.second + right.second;
return result;
}
@@ -384,7 +273,7 @@ Result fromBinary(const BinaryOperator* op) {
if (scalar >= 0) {
return scalar * value;
} else {
- return -scalar * negate_result(value);
+ return scalar * -value;
}
}
case BO_LT:
@@ -405,15 +294,13 @@ Result fromExpr(const Expr* stmt) {
return fromInteger(static_cast<const IntegerLiteral*>(stmt));
case Stmt::DeclRefExprClass:
return fromDeclExpr(static_cast<const DeclRefExpr*>(stmt));
- case Stmt::UnaryOperatorClass:
- return fromUnary(static_cast<const UnaryOperator*>(stmt));
case Stmt::BinaryOperatorClass:
return fromBinary(static_cast<const BinaryOperator*>(stmt));
}
const Expr* expr = stmt->IgnoreParenCasts();
if (stmt != expr)
return fromExpr(expr);
- llvm::errs() << "we shouldn't get here...\n";
+ assert(false);
return Result();
}
@@ -433,7 +320,7 @@ Condition fromComparison(const BinaryOperator* op, bool negate) {
const Expr* right = op->getRHS()->IgnoreParenCasts();
bool flip = false;
- std::string name;
+ string name;
int64_t value;
if (left->getStmtClass() == Stmt::DeclRefExprClass) {
name = static_cast<const DeclRefExpr*>(left)->getNameInfo().getAsString();
@@ -465,13 +352,13 @@ Condition fromComparison(const BinaryOperator* op, bool negate) {
switch (operation) {
case BO_LT:
if (negate)
- cond[negate_string(name)] = -value;
+ cond[-name] = -value;
else
cond[name] = value - 1;
break;
case BO_LE:
if (negate)
- cond[negate_string(name)] = -(value + 1);
+ cond[-name] = -(value + 1);
else
cond[name] = value;
break;
@@ -479,13 +366,13 @@ Condition fromComparison(const BinaryOperator* op, bool negate) {
if (negate)
cond[name] = value - 1;
else
- cond[negate_string(name)] = -value;
+ cond[-name] = -value;
break;
case BO_GT:
if (negate)
cond[name] = value;
else
- cond['-' + name] = -(value + 1);
+ cond[-name] = -(value + 1);
break;
}
}
@@ -493,97 +380,20 @@ Condition fromComparison(const BinaryOperator* op, bool negate) {
}
-struct MinCostFlow : public Operator<ZBar> {
- MinCostFlow(const std::vector<std::pair<EqnVar*, ZBar> >& supplies)
- : solver(0, 0) {
- assert(supplies.size() % 2 == 0);
-
- if (supplies.size() == 0)
- return;
-
- MCFClass::FNumber* deficits = new MCFClass::FNumber[supplies.size()];
- MCFClass::Index* starts = new MCFClass::Index[supplies.size()];
- MCFClass::Index* ends = new MCFClass::Index[supplies.size()];
-
- for (int i = 0, size = supplies.size(); i < size; i += 2) {
- deficits[i] = -supplies[i].second.as<MCFClass::FNumber>();
- deficits[i+1] = -supplies[i+1].second.as<MCFClass::FNumber>();
-
- starts[i] = i+1;
- ends[i] = i+2;
-
- starts[i+1] = i+2;
- ends[i+1] = i+1;
-
- printableSupply.push_back(supplies[i].second);
- printableSupply.push_back(supplies[i+1].second);
- }
-
- solver.LoadNet(supplies.size(), supplies.size(), // max nodes/arcs
- supplies.size(), supplies.size(), // current nodes/arcs
- NULL, NULL, // arcs have inf cap, zero cost (to begin)
- deficits, // deficits for each node
- starts, ends); // start/end of each edge
-
- delete[] deficits;
- delete[] starts;
- delete[] ends;
- }
-
- ZBar eval(const std::vector<ZBar>& args) const {
- if (args.size() == 0)
- return 0;
- for (int i = 0, size = args.size(); i < size; ++i) {
- //if (args[i] == infinity<ZBar>()) {
- // if (!solver.IsClosedArc(i))
- // solver.CloseArc(i);
- //} else {
- // if (solver.IsClosedArc(i))
- // solver.OpenArc(i);
- solver.ChgCost(i, args[i].as<MCFClass::CNumber>());
- //}
- }
- std::cerr << "MCF" << this << "<" << printableSupply << ">(" << args << ")" << " = ";
- solver.SolveMCF();
- if (solver.MCFGetStatus() == MCFClass::kUnfeasible){
- std::cerr << -infinity<ZBar>() << std::endl;
- return -infinity<ZBar>();
- } else {
- if (solver.MCFGetFO() == MCFClass::Inf<MCFClass::FONumber>()) {
- std::cerr << infinity<ZBar>() << std::endl;
- return infinity<ZBar>();
- } else if (solver.MCFGetFO() == -MCFClass::Inf<MCFClass::FONumber>()) {
- std::cerr << -infinity<ZBar>() << std::endl;
- return -infinity<ZBar>();
- } else {
- ZBar result = solver.MCFGetFO();
- std::cerr << result << std::endl;
- return result;
- }
- }
- }
-
- virtual void print(std::ostream& cout) const {
- cout << "MCF" << this << "<" << printableSupply << ">";
- }
-
- mutable MCFSimplex solver;
- std::vector<ZBar> printableSupply;
-};
/* Blocks */
-typedef std::map<std::string, unsigned int> Counters;
-typedef std::map<std::string, EqnVar*> VarMap;
-typedef std::map<const CFGBlock*, std::set<std::string> > BlockVars;
+typedef map<string, unsigned int> Counters;
+typedef map<string, EqnVar*> VarMap;
+typedef map<const CFGBlock*, set<string> > BlockVars;
void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
Counters counters;
- std::string block_id = toString(block->getBlockID());
+ string block_id = toString(block->getBlockID());
VarMap vars;
- for (std::set<std::string>::iterator it = block_vars[block].begin(),
+ for (set<string>::iterator it = block_vars[block].begin(),
ei = block_vars[block].end();
it != ei;
++it) {
@@ -597,7 +407,7 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
const CFGStmt* cfg_stmt = it->getAs<CFGStmt>();
const Stmt* stmt = cfg_stmt->getStmt();
- std::string name = "";
+ string name = "";
Result result;
if (stmt->getStmtClass() == Stmt::BinaryOperatorClass) {
const BinaryOperator* binop = static_cast<const BinaryOperator*>(stmt);
@@ -630,41 +440,55 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
if (name == "")
continue;
- std::string count = toString(counters[name]);
+ string count = toString(counters[name]);
EqnVar* var = &system.variable(name + '-' + block_id + '[' + count + ']');
- EqnVar* negative_var = &system.variable(negate_string(name) + '-' + block_id + '[' + count + ']');
+ EqnVar* negative_var = &system.variable(-name + '-' + block_id + '[' + count + ']');
counters[name]++;
for (int negative = 0; negative < 2; ++negative) { // one loop for positive, the other for negative
if (negative) {
- result = negate_result(result);
+ result = -result;
}
EqnExpr* expression;
if (result.first.size() > 0) {
- // set up the min-cost stuff
- std::vector<EqnExpr*> minCostArgs;
- std::vector<std::pair<EqnVar*, ZBar> > varCosts;
- for (Vector::iterator it = result.first.begin(),
+ // make sure all our variables exist in vars
+ for (Vector::iterator
+ it = result.first.begin(),
ei = result.first.end();
it != ei;
++it) {
if (!vars[it->first])
vars[it->first] = &system.variable(it->first + '-' + block_id + "-pre");
+ }
- EqnVar& var = *vars[it->first];
- EqnVar& negVar = *vars[negate_string(it->first)];
-
- varCosts.push_back(std::pair<EqnVar*, ZBar>(&var, it->second));
- varCosts.push_back(std::pair<EqnVar*, ZBar>(&negVar, -it->second));
-
- minCostArgs.push_back(&var);
- minCostArgs.push_back(&negVar);
+ // set up the min-cost-flow operator
+ vector<ZBar> supplies;
+ vector<pair<int,int> > arcs;
+ vector<EqnExpr*> minCostArgs;
+ ZBar dummy_value = 0;
+ supplies.push_back(dummy_value); // dummy node to suck up flow
+ int index = 1; // the solver uses 1-indexing, for some reason
+ for (map<std::string,EqnVar*>::iterator
+ it = vars.begin(),
+ ei = vars.end();
+ it != ei;
+ it++) {
+ index++;
+ supplies.push_back(result.first[it->first]);
+ dummy_value -= result.first[it->first];
+ if (it->first[0] == '-')
+ arcs.push_back(pair<int,int>(1,index));
+ else
+ arcs.push_back(pair<int,int>(index,1));
+ minCostArgs.push_back(vars[it->first]);
}
- EqnExpr* minCostExpr = &system.expression(new MinCostFlow(varCosts), minCostArgs);
+ supplies[0] = dummy_value;
+
+ EqnExpr* minCostExpr = &system.expression(new MinCostFlow<ZBar>(supplies, arcs), minCostArgs);
// add the constant factor to the min-cost bit
- std::vector<EqnExpr*> additionArgs;
+ vector<EqnExpr*> additionArgs;
additionArgs.push_back(&system.constant(result.second));
additionArgs.push_back(minCostExpr);
expression = &system.expression(new Addition<ZBar>(), additionArgs);
@@ -673,7 +497,7 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
}
// max(-inf, expr), so strategy iteration will work
- std::vector<EqnExpr*> maxArgs;
+ vector<EqnExpr*> maxArgs;
maxArgs.push_back(&system.constant(-infinity<ZBar>()));
maxArgs.push_back(expression);
if (negative)
@@ -682,9 +506,9 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
system[*var] = &system.maxExpression(maxArgs);
}
vars[name] = var;
- vars[negate_string(name)] = negative_var;
+ vars[-name] = negative_var;
block_vars[block].insert(name);
- block_vars[block].insert(negate_string(name));
+ block_vars[block].insert(-name);
}
// add to our successor entry values
@@ -704,7 +528,7 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
ZBar val = cond[jt->first];
EqnVar* var = &system.variable(jt->first + '-' + toString((*it)->getBlockID()) + "-pre");
if (system[*var] == NULL) {
- std::vector<EqnExpr*> maxArgs;
+ vector<EqnExpr*> maxArgs;
maxArgs.push_back(&system.constant(-infinity<ZBar>()));
system[*var] = &system.maxExpression(maxArgs);
}
@@ -718,51 +542,13 @@ void runOnBlock(const CFGBlock* block, EqnSys& system, BlockVars& block_vars) {
expr = jt->second;
} else {
// need a min here
- std::vector<EqnExpr*> minArgs;
+ vector<EqnExpr*> minArgs;
minArgs.push_back(&system.constant(val));
minArgs.push_back(jt->second);
expr = &system.expression(new Minimum<ZBar>(), minArgs);
}
system[*var]->arguments().push_back(expr);
-
- /*if (expr) {
- std::set<std::string> ignore;
- for (VarMap::iterator
- variables = vars.begin(),
- variables_end = vars.end();
- variables != variables_end;
- ++variables) {
- if (ignore.find(variables->first) != ignore.end())
- continue;
- ignore.insert(negate_string(variables->first));
-
- std::vector<EqnExpr*> plusArgs;
- for (int negate = 0; negate < 2; ++negate) {
- std::string var_name = negate ? negate_string(variables->first) : variables->first;
- if (cond[var_name] == infinity<ZBar>()) {
- // min(x, inf) = x
- plusArgs.push_back(vars[var_name]);
- } else if (cond[var_name] == -infinity<ZBar>()) {
- // min(x, -inf) = -inf
- plusArgs.push_back(&system.constant(-infinity<ZBar>()));
- } else {
- // min(x, y) = ?
- std::vector<EqnExpr*> minArgs;
- minArgs.push_back(vars[var_name]);
- minArgs.push_back(&system.constant(cond[var_name]));
- plusArgs.push_back(&system.expression(new Minimum<ZBar>(), minArgs));
- }
- }
-
- std::vector<EqnExpr*> guard_args;
- guard_args.push_back(&system.expression(new Addition<ZBar>(), plusArgs)); // value
- guard_args.push_back(&system.constant(0)); // lower-bound (so value must be >= this)
- guard_args.push_back(expr); // result
- expr = &system.expression(new Guard<ZBar>(), guard_args);
- }
- system[*var]->arguments().push_back(expr);
- }*/
}
}
}
@@ -788,27 +574,27 @@ void IntervalAnalysis::runOnAllBlocks(const Decl& decl) {
EqnSys system;
BlockVars block_vars;
- std::vector<EqnExpr*> infArg;
+ vector<EqnExpr*> infArg;
infArg.push_back(&system.constant(-infinity<ZBar>())); // left-most argument has to be -infinity
infArg.push_back(&system.constant(infinity<ZBar>()));
- std::set<std::string>& function_arguments = block_vars[&cfg->getEntry()];
- std::string block_id = toString(cfg->getEntry().getBlockID());
+ set<string>& function_arguments = block_vars[&cfg->getEntry()];
+ string block_id = toString(cfg->getEntry().getBlockID());
if (const FunctionDecl* func = dyn_cast<const FunctionDecl>(&decl)) {
for (unsigned int i = func->getNumParams(); i > 0; i--) {
- std::string name = func->getParamDecl(i-1)->getNameAsString();
+ string name = func->getParamDecl(i-1)->getNameAsString();
// add the variables to the first block
function_arguments.insert(name);
- function_arguments.insert(negate_string(name));
+ function_arguments.insert(neg(name));
// set the vars to infinity (unconstrained)
system[system.variable(name + '-' + block_id + "-pre")] = &system.maxExpression(infArg);
- system[system.variable(negate_string(name) + '-' + block_id + "-pre")] = &system.maxExpression(infArg);
+ system[system.variable(neg(name) + '-' + block_id + "-pre")] = &system.maxExpression(infArg);
}
}
- std::set<const CFGBlock*> seen;
- std::deque<const CFGBlock*> todo;
+ set<const CFGBlock*> seen;
+ deque<const CFGBlock*> todo;
todo.push_back(&cfg->getEntry());
while (!todo.empty()) {
@@ -835,7 +621,7 @@ void IntervalAnalysis::runOnAllBlocks(const Decl& decl) {
Solver<ZBar> solver(system);
for (unsigned int i = 0, size = system.variableCount(); i < size; ++i) {
- EqnVar& var = system.variable(size - i - 1);
+ EqnVar& var = system.variable(i);
llvm::errs() << toString(var.name()) << " = " << toString(solver.solve(var)) << "\n";
}
}
diff --git a/clang/lib/Analysis/MCFClass.h b/clang/lib/Analysis/MCFClass.h
deleted file mode 100755
index 41adb82..0000000
--- a/clang/lib/Analysis/MCFClass.h
+++ /dev/null
@@ -1,2436 +0,0 @@
-/*--------------------------------------------------------------------------*/
-/*-------------------------- File MCFClass.h -------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @file
- * Header file for the abstract base class MCFClass, which defines a standard
- * interface for (linear or convex quadratic separable) Min Cost Flow Problem
- * solvers, to be implemented as derived classes.
- *
- * \version 3.01
- *
- * \date 30 - 09 - 2011
- *
- * \author Alessandro Bertolini \n
- * Operations Research Group \n
- * Dipartimento di Informatica \n
- * Universita' di Pisa \n
- *
- * \author Antonio Frangioni \n
- * Operations Research Group \n
- * Dipartimento di Informatica \n
- * Universita' di Pisa \n
- *
- * \author Claudio Gentile \n
- * Istituto di Analisi di Sistemi e Informatica \n
- * Consiglio Nazionale delle Ricerche \n
- *
- * Copyright &copy 1996 - 2011 by Antonio Frangioni, Claudio Gentile
- */
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*----------------------------- DEFINITIONS --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#ifndef __MCFClass
- #define __MCFClass /* self-identification: #endif at the end of the file */
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------------- MACROS ---------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup MCFCLASS_MACROS Compile-time switches in MCFClass.h
- These macros control some important details of the class interface.
- Although using macros for activating features of the interface is not
- very C++, switching off some unused features may allow some
- implementation to be more efficient in running time or memory.
- @{ */
-
-/*-------------------------------- USENAME0 --------------------------------*/
-
-#define USENAME0 1
-
-/**< Decides if 0 or 1 is the "name" of the first node.
- If USENAME0 == 1, (warning: it has to be *exactly* 1), then the node
- names go from 0 to n - 1, otherwise from 1 to n. Note that this does not
- affect the position of the deficit in the deficit vectors, i.e., the
- deficit of the i-th node - be its "name" `i' or `i - 1' - is always in
- the i-th position of the vector. */
-
-/*@} end( group( MCFCLASS_MACROS ) ) */
-/*--------------------------------------------------------------------------*/
-/*------------------------------ INCLUDES ----------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#include "OPTUtils.h"
-
-/* OPTtypes.h defines standard interfaces for timing and random routines, as
- well as the namespace OPTtypes_di_unipi_it and the macro
- OPT_USE_NAMESPACES, useful for switching off all namespaces in one blow
- for those strange cases where they create problems. */
-
-#include <iomanip>
-#include <sstream>
-#include <limits>
-#include <cassert>
-
-// QUICK HACK
-#define throw(x) assert(false)
-
-/*--------------------------------------------------------------------------*/
-/*------------------------- NAMESPACE and USING ----------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
-namespace MCFClass_di_unipi_it
-{
- /** @namespace MCFClass_di_unipi_it
- The namespace MCFClass_di_unipi_it is defined to hold the MCFClass
- class and all the relative stuff. It comprises the namespace
- OPTtypes_di_unipi_it. */
-
- using namespace OPTtypes_di_unipi_it;
-#endif
-
-/*@} end( group( MCFCLASS_CONSTANTS ) ) */
-/*--------------------------------------------------------------------------*/
-/*-------------------------- CLASS MCFClass --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup MCFCLASS_CLASSES Classes in MCFClass.h
- @{ */
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------- GENERAL NOTES --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** This abstract base class defines a standard interface for (linear or
- convex quadartic separable) Min Cost Flow (MCF) problem solvers.
-
- The data of the problem consist of a (directed) graph G = ( N , A ) with
- n = |N| nodes and m = |A| (directed) arcs. Each node `i' has a deficit
- b[ i ], i.e., the amount of flow that is produced/consumed by the node:
- source nodes (which produce flow) have negative deficits and sink nodes
- (which consume flow) have positive deficits. Each arc `(i, j)' has an
- upper capacity U[ i , j ], a linear cost coefficient C[ i , j ] and a
- (non negative) quadratic cost coefficient Q[ i , j ]. Flow variables
- X[ i , j ] represents the amount of flow to be sent on arc (i, j).
- Parallel arcs, i.e., multiple copies of the same arc `(i, j)' (with
- possibily different costs and/or capacities) are in general allowed.
- The formulation of the problem is therefore:
- \f[
- \min \sum_{ (i, j) \in A } C[ i , j ] X[ i, j ] +
- Q[ i , j ] X[ i, j ]^2 / 2
- \f]
- \f[
- (1) \sum_{ (j, i) \in A } X[ j , i ] -
- \sum_{ (i, j) \in A } X[ i , j ] = b[ i ]
- \hspace{1cm} i \in N
- \f]
- \f[
- (2) 0 \leq X[ i , j ] \leq U[ i , j ]
- \hspace{1cm} (i, j) \in A
- \f]
- The n equations (1) are the flow conservation constraints and the 2m
- inequalities (2) are the flow nonnegativity and capacity constraints.
- At least one of the flow conservation constraints is redundant, as the
- demands must be balanced (\f$\sum_{ i \in N } b[ i ] = 0\f$); indeed,
- exactly n - ConnectedComponents( G ) flow conservation constraints are
- redundant, as demands must be balanced in each connected component of G.
- Let us denote by QA and LA the disjoint subsets of A containing,
- respectively, "quadratic" arcs (with Q[ i , j ] > 0) and "linear" arcs
- (with Q[ i , j ] = 0); the (MCF) problem is linear if QA is empty, and
- nonlinear (convex quadratic) if QA is nonempty.
-
- The dual of the problem is:
- \f[
- \max \sum_{ i \in N } Pi[ i ] b[ i ] -
- \sum_{ (i, j) \in A } W[ i , j ] U[ i , j ] -
- \sum_{ (i, j) \in AQ } V[ i , j ]^2 / ( 2 * Q[ i , j ] )
- \f]
- \f[
- (3.a) C[ i , j ] - Pi[ j ] + Pi[ i ] + W[ i , j ] - Z[ i , j ] = 0
- \hspace{1cm} (i, j) \in AL
- \f]
- \f[
- (3.b) C[ i , j ] - Pi[ j ] + Pi[ i ] + W[ i , j ] - Z[ i , j ] =
- V[ i , j ]
- \hspace{1cm} (i, j) \in AQ
- \f]
- \f[
- (4.a) W[ i , j ] \geq 0 \hspace{1cm} (i, j) \in A
- \f]
- \f[
- (4.b) Z[ i , j ] \geq 0 \hspace{1cm} (i, j) \in A
- \f]
-
- Pi[] is said the vector of node potentials for the problem, W[] are
- bound variables and Z[] are slack variables. Given Pi[], the quantities
- \f[
- RC[ i , j ] = C[ i , j ] + Q[ i , j ] * X[ i , j ] - Pi[ j ] + Pi[ i ]
- \f]
- are said the "reduced costs" of arcs.
-
- A primal and dual feasible solution pair is optimal if and only if the
- complementary slackness conditions
- \f[
- RC[ i , j ] > 0 \Rightarrow X[ i , j ] = 0
- \f]
- \f[
- RC[ i , j ] < 0 \Rightarrow X[ i , j ] = U[ i , j ]
- \f]
- are satisfied for all arcs (i, j) of A.
-
- The MCFClass class provides an interface with methods for managing and
- solving problems of this kind. Actually, the class can also be used as
- an interface for more general NonLinear MCF problems, where the cost
- function either nonseparable ( C( X ) ) or arc-separable
- ( \f$\sum_{ (i, j) \in A } C_{i,j}( X[ i, j ] )\f$ ). However, solvers
- for NonLinear MCF problems are typically objective-function-specific,
- and there is no standard way for inputting a nonlinear function different
- from a separable convex quadratic one, so only the simplest form is dealt
- with in the interface, leaving more complex NonLinear parts to the
- interface of derived classes. */
-
-class MCFClass {
-
-/*--------------------------------------------------------------------------*/
-/*----------------------- PUBLIC PART OF THE CLASS -------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-- --*/
-/*-- The following methods and data are the actual interface of the --*/
-/*-- class: the standard user should use these methods and data only. --*/
-/*-- --*/
-/*--------------------------------------------------------------------------*/
-
- public:
-
-/*--------------------------------------------------------------------------*/
-/*---------------------------- PUBLIC TYPES --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Public types
- The MCFClass defines four main public types:
-
- - Index, the type of arc and node indices;
-
- - FNumber, the type of flow variables, arc capacities, and node deficits;
-
- - CNumber, the type of flow costs, node potentials, and arc reduced costs;
-
- - FONumber, the type of objective function value.
-
- By re-defining the types in this section, most MCFSolver should be made
- to work with any reasonable choice of data type (= one that is capable of
- properly representing the data of the instances to be solved). This may
- be relevant due to an important property of MCF problems: *if all arc
- capacities and node deficits are integer, then there exists an integral
- optimal primal solution*, and *if all arc costs are integer, then there
- exists an integral optimal dual solution*. Even more importantly, *many
- solution algorithms will in fact produce an integral primal/dual
- solution for free*, because *every primal/dual solution they generate
- during the solution process is naturally integral*. Therefore, one can
- use integer data types to represent everything connected with flows and/or
- costs if the corresponding data is integer in all instances one needs to
- solve. This directly translates in significant memory savings and/or speed
- improvements.
-
- *It is the user's responsibility to ensure that these types are set to
- reasonable values*. So, the experienced user may want to experiment with
- setting this data properly if memory footprint and/or speed is a primary
- concern. Note, however, that *not all solution algorithms will happily
- accept integer data*; one example are Interior-Point approaches, which
- require both flow and cost variables to be continuous (float). So, the
- viability of setting integer data (as well as its impact on performances)
- is strictly related to the specific kind of algorithm used. Since these
- types are common to all derived classes, they have to be set taking into
- account the needs of all the solvers that are going to be used, and
- adapting to the "worst case"; of course, FNumber == CNumber == double is
- going to always be an acceptable "worst case" setting. MCFClass may in a
- future be defined as a template class, with these as template parameters,
- but this is currently deemed overkill and avoided.
-
- Finally, note that the above integrality property only holds for *linear*
- MCF problems. If any arc has a nonzero quadratic cost coefficient, optimal
- flows and potentials may be fractional even if all the data of the problem
- (comprised quadratic cost coefficients) is integer. Hence, for *quadratic*
- MCF solvers, a setting like FNumber == CNumber == double is actually
- *mandatory*, for any reasonable algorithm will typically misbehave
- otherwise.
- @{ */
-
-/*--------------------------------------------------------------------------*/
-
- typedef unsigned int Index; ///< index of a node or arc ( >= 0 )
- typedef Index *Index_Set; ///< set (array) of indices
- typedef const Index cIndex; ///< a read-only index
- typedef cIndex *cIndex_Set; ///< read-only index array
-
-/*--------------------------------------------------------------------------*/
-
- typedef int SIndex; ///< index of a node or arc
- typedef SIndex *SIndex_Set; ///< set (array) of indices
- typedef const SIndex cSIndex; ///< a read-only index
- typedef cSIndex *cSIndex_Set; ///< read-only index array
-
-/*--------------------------------------------------------------------------*/
-
- typedef double FNumber; ///< type of arc flow
- typedef FNumber *FRow; ///< vector of flows
- typedef const FNumber cFNumber; ///< a read-only flow
- typedef cFNumber *cFRow; ///< read-only flow array
-
-/*--------------------------------------------------------------------------*/
-
- typedef double CNumber; ///< type of arc flow cost
- typedef CNumber *CRow; ///< vector of costs
- typedef const CNumber cCNumber; ///< a read-only cost
- typedef cCNumber *cCRow; ///< read-only cost array
-
-/*--------------------------------------------------------------------------*/
-
- typedef double FONumber;
- /**< type of the objective function: has to hold sums of products of
- FNumber(s) by CNumber(s) */
-
- typedef const FONumber cFONumber; ///< a read-only o.f. value
-
-/*--------------------------------------------------------------------------*/
-/** Very small class to simplify extracting the "+ infinity" value for a
- basic type (FNumber, CNumber, Index); just use Inf<type>(). */
-
- template <typename T>
- class Inf {
- public:
- Inf() {}
- operator T() { return( std::numeric_limits<T>::max() ); }
- };
-
-/*--------------------------------------------------------------------------*/
-/** Very small class to simplify extracting the "machine epsilon" for a
- basic type (FNumber, CNumber); just use Eps<type>(). */
-
- template <typename T>
- class Eps {
- public:
- Eps() {}
- operator T() { return( std::numeric_limits<T>::epsilon() ); }
- };
-
-/*--------------------------------------------------------------------------*/
-/** Small class for exceptions. Derives from std::exception implementing the
- virtual method what() -- and since what is virtual, remember to always
- catch it by reference (catch exception &e) if you want the thing to work.
- MCFException class are thought to be of the "fatal" type, i.e., problems
- for which no solutions exists apart from aborting the program. Other kinds
- of exceptions can be properly handled by defining derived classes with
- more information. */
-
- class MCFException : public exception {
- public:
- MCFException( const char *const msg = 0 ) { errmsg = msg; }
-
- // wow, such a hack
-#undef throw
- const char* what( void ) const throw () { return( errmsg ); }
-#define throw(x) assert(false)
- private:
- const char *errmsg;
- };
-
-/*--------------------------------------------------------------------------*/
-/** Public enum describing the possible parameters of the MCF solver, to be
- used with the methods SetPar() and GetPar(). */
-
- enum MCFParam { kMaxTime = 0 , ///< max time
- kMaxIter , ///< max number of iteration
- kEpsFlw , ///< tolerance for flows
- kEpsDfct , ///< tolerance for deficits
- kEpsCst , ///< tolerance for costs
- kReopt , ///< whether or not to reoptimize
- kLastParam /**< dummy parameter: this is used to
- allow derived classes to "extend"
- the set of parameters. */
- };
-
-/*--------------------------------------------------------------------------*/
-/** Public enum describing the possible status of the MCF solver. */
-
- enum MCFStatus { kUnSolved = -1 , ///< no solution available
- kOK = 0 , ///< optimal solution found
-
- kStopped , ///< optimization stopped
- kUnfeasible , ///< problem is unfeasible
- kUnbounded , ///< problem is unbounded
- kError ///< error in the solver
- };
-
-/*--------------------------------------------------------------------------*/
-/** Public enum describing the possible reoptimization status of the MCF
- solver. */
-
- enum MCFAnswer { kNo = 0 , ///< no
- kYes ///< yes
- };
-
-/*--------------------------------------------------------------------------*/
-/** Public enum describing the possible file formats in WriteMCF(). */
-
- enum MCFFlFrmt { kDimacs = 0 , ///< DIMACS file format for MCF
- kQDimacs , ///< quadratic DIMACS file format for MCF
- kMPS , ///< MPS file format for LP
- kFWMPS ///< "Fixed Width" MPS format
- };
-
-/*--------------------------------------------------------------------------*/
-/** Base class for representing the internal state of the MCF algorithm. */
-
- class MCFState {
- public:
- MCFState( void ) {};
- virtual ~MCFState() {};
- };
-
- typedef MCFState *MCFStatePtr; ///< pointer to a MCFState
-
-/*@} -----------------------------------------------------------------------*/
-/*--------------------------- PUBLIC METHODS -------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*---------------------------- CONSTRUCTOR ---------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Constructors
- @{ */
-
- MCFClass( cIndex nmx = 0 , cIndex mmx = 0 )
- {
- nmax = nmx;
- mmax = mmx;
- n = m = 0;
-
- status = kUnSolved;
- Senstv = true;
-
- EpsFlw = Eps<FNumber>() * 100;
- EpsCst = Eps<CNumber>() * 100;
- EpsDfct = EpsFlw * ( nmax ? nmax : 100 );
-
- MaxTime = 0;
- MaxIter = 0;
-
- MCFt = NULL;
- }
-
-/**< Constructor of the class.
-
- nmx and mmx, if provided, are taken to be respectively the maximum number
- of nodes and arcs in the network. If nonzero values are passed, memory
- allocation can be anticipated in the constructor, which is sometimes
- desirable. The maximum values are stored in the protected fields nmax and
- mmax, and can be changed with LoadNet() [see below]; however, changing
- them typically requires memory allocation/deallocation, which is
- sometimes undesirable outside the constructor.
-
- After that an object has been constructed, no problem is loaded; this has
- to be done with LoadNet() [see below]. Thus, it is an error to invoke any
- method which requires the presence of a problem (typicall all except those
- in the initializations part). The base class provides two protected fields
- n and m for the current number of nodes and arcs, respectively, that are
- set to 0 in the constructor precisely to indicate that no instance is
- currently loaded. */
-
-/*@} -----------------------------------------------------------------------*/
-/*-------------------------- OTHER INITIALIZATIONS -------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Other initializations
- @{ */
-
- virtual void LoadNet( cIndex nmx = 0 , cIndex mmx = 0 , cIndex pn = 0 ,
- cIndex pm = 0 , cFRow pU = NULL , cCRow pC = NULL ,
- cFRow pDfct = NULL , cIndex_Set pSn = NULL ,
- cIndex_Set pEn = NULL ) = 0;
-
-/**< Inputs a new network.
-
- The parameters nmx and mmx are the new max number of nodes and arcs,
- possibly overriding those set in the constructor [see above], altough at
- the likely cost of memory allocation and deallocation. Passing nmx ==
- mmx == 0 is intended as a signal to the solver to deallocate everything
- and wait for new orders; in this case, all the other parameters are ignored.
-
- Otherwise, in principle all the other parameters have to be provided.
- Actually, some of them may not be needed for special classes of MCF
- problems (e.g., costs in a MaxFlow problem, or start/end nodes in a
- problem defined over a graph with fixed topology, such as a complete
- graph). Also, passing NULL is allowed to set default values.
-
- The meaning of the parameters is the following:
-
- - pn is the current number of nodes of the network (<= nmax).
- - pm is the number of arcs of the network (<= mmax).
-
- - pU is the m-vector of the arc upper capacities; capacities must be
- nonnegative, but can in principle be infinite (== F_INF); passing
- pU == NULL means that all capacities are infinite;
-
- - pC is the m-vector of the arc costs; costs must be finite (< C_INF);
- passing pC == NULL means that all costs must be 0.
-
- - pDfct is the n-vector of the node deficits; source nodes have negative
- deficits and sink nodes have positive deficits; passing pDfct ==
- NULL means that all deficits must be 0 (a circulation problem);
-
- - pSn is the m-vector of the arc starting nodes; pSn == NULL is in
- principle not allowed, unless the topology of the graph is fixed;
- - pEn is the m-vector of the arc ending nodes; same comments as for pSn.
-
- Note that node "names" in the arrays pSn and pEn must go from 1 to pn if
- the macro USANAME0 [see above] is set to 0, while they must go from 0 to
- pn - 1 if USANAME0 is set to 1. In both cases, however, the deficit of the
- first node is read from the first (0-th) position of pDfct, that is if
- USANAME0 == 0 then the deficit of the node with name `i' is read from
- pDfct[ i - 1 ].
-
- The data passed to LoadNet() can be used to specify that the arc `i' must
- not "exist" in the problem. This is done by passing pC[ i ] == C_INF;
- solvers which don't read costs are forced to read them in order to check
- this, unless they provide alternative solver-specific ways to accomplish
- the same tasks. These arcs are "closed", as for the effect of CloseArc()
- [see below]. "invalid" costs (== C_INF) are set to 0 in order to being
- subsequently capable of "opening" them back with OpenArc() [see below].
- The way in which these non-existent arcs are phisically dealt with is
- solver-specific; in some solvers, for instance, this could be obtained by
- simply putting their capacity to zero. Details about these issues should
- be found in the interface of derived classes.
-
- Note that the quadratic part of the objective function, if any, is not
- dealt with in LoadNet(); it can only be separately provided with
- ChgQCoef() [see below]. By default, the problem is linear, i.e., all
- coefficients of the second-order terms in the objective function are
- assumed to be zero. */
-
-/*--------------------------------------------------------------------------*/
-
- virtual inline void LoadDMX( istream &DMXs , bool IsQuad = false );
-
-/**< Read a MCF instance in DIMACS standard format from the istream. The
- format is the following. The first line must be
-
- p min <number of nodes> <number of arcs>
-
- Then the node definition lines must be found, in the form
-
- n <node number> <node supply>
-
- Not all nodes need have a node definition line; these are given zero
- supply, i.e., they are transhipment nodes (supplies are the inverse of
- deficits, i.e., a node with positive supply is a source node). Finally,
- the arc definition lines must be found, in the form
-
- a <start node> <end node> <lower bound> <upper bound> <flow cost>
-
- There must be exactly <number of arcs> arc definition lines in the file.
-
- This method is *not* pure virtual because an implementation is provided by
- the base class, using the LoadNet() method (which *is* pure virtual).
- However, the method *is* virtual to allow derived classes to implement
- more efficient versions, should they have any reason to do so.
-
- \note Actually, the file format accepted by LoadDMX (at least in the
- base class implementation) is more general than the DIMACS standard
- format, in that it is allowed to mix node and arc definitions in
- any order, while the DIMACS file requires all node information to
- appear before all arc information.
-
- \note Other than for the above, this method is assumed to allow for
- *quadratic* Dimacs files, encoding for convex quadratic separable
- Min Cost Flow instances. This is a simple extension where each arc
- descriptor has a sixth field, <quadratic cost>. The provided
- istream is assumed to be quadratic Dimacs file if IsQuad is true,
- and a regular linear Dimacs file otherwise. */
-
-/*--------------------------------------------------------------------------*/
-
- virtual void PreProcess( void ) {}
-
-/**< Extract a smaller/easier equivalent MCF problem. The data of the instance
- is changed and the easier one is solved instead of the original one. In the
- MCF case, preprocessing may involve reducing bounds, identifying
- disconnected components of the graph etc. However, proprocessing is
- solver-specific.
-
- This method can be implemented by derived classes in their solver-specific
- way. Preprocessing may reveal unboundedness or unfeasibility of the
- problem; if that happens, PreProcess() should properly set the `status'
- field, that can then be read with MCFGetStatus() [see below].
-
- Note that preprocessing may destroy all the solution information. Also, it
- may be allowed to change the data of the problem, such as costs/capacities
- of the arcs.
-
- A valid preprocessing is doing nothing, and that's what the default
- implementation of this method (that is *not* pure virtual) does. */
-
-/*--------------------------------------------------------------------------*/
-
- virtual inline void SetPar( int par , int val );
-
-/**< Set integer parameters of the algorithm.
-
- @param par is the parameter to be set; the enum MCFParam can be used, but
- 'par' is an int (every enum is an int) so that the method can
- be extended by derived classes for the setting of their
- parameters
-
- @param value is the value to assign to the parameter.
-
- The base class implementation handles these parameters:
-
- - kMaxIter: the max number of iterations in which the MCF Solver can find
- an optimal solution (default 0, which means no limit)
-
- - kReopt: tells the solver if it has to reoptimize. The implementation in
- the base class sets a flag, the protected \c bool field \c
- Senstv; if true (default) this field instructs the MCF solver to
- to try to exploit the information about the latest optimal
- solution to speedup the optimization of the current problem,
- while if the field is false the MCF solver should restart the
- optimization "from scratch" discarding any previous information.
- Usually reoptimization speeds up the computation considerably,
- but this is not always true, especially if the data of the
- problem changes a lot.*/
-
-/*--------------------------------------------------------------------------*/
-
- virtual inline void SetPar( int par , double val );
-
-/**< Set float parameters of the algorithm.
-
- @param par is the parameter to be set; the enum MCFParam can be used, but
- 'par' is an int (every enum is an int) so that the method can
- be extended by derived classes for the setting of their
- parameters
-
- @param value is the value to assign to the parameter.
-
- The base class implementation handles these parameters:
-
- - kEpsFlw: sets the tolerance for controlling if the flow on an arc is zero
- to val. This also sets the tolerance for controlling if a node
- deficit is zero (see kEpsDfct) to val * < max number of nodes >;
- this value should be safe for graphs in which any node has less
- than < max number of nodes > adjacent nodes, i.e., for all graphs
- but for very dense ones with "parallel arcs"
-
- - kEpsDfct: sets the tolerance for controlling if a node deficit is zero to
- val, in case a better value than that autmatically set by
- kEpsFlw (see above) is available (e.g., val * k would be good
- if no node has more than k neighbours)
-
- - kEpsCst: sets the tolerance for controlling if the reduced cost of an arc
- is zero to val. A feasible solution satisfying eps-complementary
- slackness, i.e., such that
- \f[
- RC[ i , j ] < - eps \Rightarrow X[ i , j ] = U[ ij ]
- \f]
- and
- \f[
- RC[ i , j ] > eps \Rightarrow X[ i , j ] == 0 ,
- \f]
- is known to be ( eps * n )-optimal.
-
- - kMaxTime: sets the max time (in seconds) in which the MCF Solver can find
- an optimal solution (default 0, which means no limit). */
-
-/*--------------------------------------------------------------------------*/
-
- virtual inline void GetPar( int par , int &val );
-
-/**< This method returns one of the integer parameter of the algorithm.
-
- @param par is the parameter to return [see SetPar( int ) for comments];
-
- @param val upon return, it will contain the value of the parameter.
-
- The base class implementation handles the parameters kMaxIter and kReopt.
- */
-
-/*--------------------------------------------------------------------------*/
-
- virtual inline void GetPar( int par , double &val );
-
-/**< This method returns one of the integer parameter of the algorithm.
-
- @param par is the parameter to return [see SetPar( double ) for comments];
-
- @param val upon return, it will contain the value of the parameter.
-
- The base class implementation handles the parameters kEpsFlw, kEpsDfct,
- kEpsCst, and kMaxTime. */
-
-/*--------------------------------------------------------------------------*/
-
- virtual void SetMCFTime( bool TimeIt = true )
- {
- if( TimeIt )
- if( MCFt )
- MCFt->ReSet();
- else
- MCFt = new OPTtimers();
- else {
- delete MCFt;
- MCFt = NULL;
- }
- }
-
-/**< Allocate an OPTtimers object [see OPTtypes.h] to be used for timing the
- methods of the class. The time can be read with TimeMCF() [see below]. By
- default, or if SetMCFTime( false ) is called, no timing is done. Note that,
- since all the relevant methods ot the class are pure virtual, MCFClass can
- only manage the OPTtimers object, but it is due to derived classes to
- actually implement the timing.
-
- @note time accumulates over the calls: calling SetMCFTime(), however,
- resets the counters, allowing to time specific groups of calls.
-
- @note of course, setting kMaxTime [see SetPar() above] to any nonzero
- value has no effect unless SetMCFTime( true ) has been called. */
-
-/*@} -----------------------------------------------------------------------*/
-/*---------------------- METHODS FOR SOLVING THE PROBLEM -------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Solving the problem
- @{ */
-
- virtual void SolveMCF( void ) = 0;
-
-/**< Solver of the Min Cost Flow Problem. Attempts to solve the MCF instance
- currently loaded in the object. */
-
-/*--------------------------------------------------------------------------*/
-
- inline int MCFGetStatus( void )
- {
- return( status );
- }
-
-/**< Returns an int describing the current status of the MCF solver. Possible
- return values are:
-
- - kUnSolved SolveMCF() has not been called yet, or the data of the
- problem has been changed since the last call;
-
- - kOK optimization has been carried out succesfully;
-
- - kStopped optimization have been stopped before that the stopping
- conditions of the solver applied, e.g. because of the
- maximum allowed number of "iterations" [see SetPar( int )]
- or the maximum allowed time [see SetPar( double )] has been
- reached; this is not necessarily an error, as it might just
- be required to re-call SolveMCF() giving it more "resources"
- in order to solve the problem;
-
- - kUnfeasible if the current MCF instance is (primal) unfeasible;
-
- - kUnbounded if the current MCF instance is (primal) unbounded (this can
- only happen if the solver actually allows F_INF capacities,
- which is nonstandard in the interface);
-
- - kError if there was an error during the optimization; this
- typically indicates that computation cannot be resumed,
- although solver-dependent ways of dealing with
- solver-dependent errors may exist.
-
- MCFClass has a protected \c int \c member \c status \c that can be used by
- derived classes to hold status information and that is returned by the
- standard implementation of this method. Note that \c status \c is an
- \c int \c and not an \c enum \c, and that an \c int \c is returned by this
- method, in order to allow the derived classes to extend the set of return
- values if they need to do so. */
-
-/*@} -----------------------------------------------------------------------*/
-/*---------------------- METHODS FOR READING RESULTS -----------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Reading flow solution
- @{ */
-
- virtual void MCFGetX( FRow F , Index_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Write the optimal flow solution in the vector F[]. If nms == NULL, F[]
- will be in "dense" format, i.e., the flow relative to arc `i'
- (i in 0 .. m - 1) is written in F[ i ]. If nms != NULL, F[] will be in
- "sparse" format, i.e., the indices of the nonzero elements in the flow
- solution are written in nms (that is then Inf<Index>()-terminated) and
- the flow value of arc nms[ i ] is written in F[ i ]. Note that nms is
- *not* guaranteed to be ordered. Also, note that, unlike MCFGetRC() and
- MCFGetPi() [see below], nms is an *output* of the method.
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cFRow MCFGetX( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal data structure containing the
- flow solution in "dense" format. Since this may *not always be available*,
- depending on the implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual bool HaveNewX( void )
- {
- return( false );
- }
-
-/**< Return true if a different (approximately) optimal primal solution is
- available. If the method returns true, then any subsequent call to (any
- form of) MCFGetX() will return a different primal solution w.r.t. the one
- that was being returned *before* the call to HaveNewX(). This solution need
- not be optimal (although, ideally, it has to be "good); this can be checked
- by comparing its objective function value, that will be returned by a call
- to MCFGetFO() [see below].
-
- Any subsequent call of HaveNewX() that returns true produces a new
- solution, until the first that returns false; from then on, no new
- solutions will be generated until something changes in the problem's
- data.
-
- Note that a default implementation of HaveNewX() is provided which is
- good for those solvers that only produce one optimal primal solution. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading potentials
- @{ */
-
- virtual void MCFGetPi( CRow P , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Writes the optimal node potentials in the vector P[]. If nms == NULL,
- the node potential of node `i' (i in 0 .. n - 1) is written in P[ i ]
- (note that here node names always start from zero, regardless to the value
- of USENAME0). If nms != NULL, it must point to a vector of indices in
- 0 .. n - 1 (ordered in increasing sense and Inf<Index>()-terminated), and
- the node potential of nms[ i ] is written in P[ i ]. Note that, unlike
- MCFGetX() above, nms is an *input* of the method.
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the nodes `i' with strt <= i < min( MCFn() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to nodes which are *both* in nms[] and whose index is
- in the correct range are returned. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cCRow MCFGetPi( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal data structure containing the
- node potentials. Since this may *not always be available*, depending on
- the implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual bool HaveNewPi( void )
- {
- return( false );
- }
-
-/**< Return true if a different (approximately) optimal dual solution is
- available. If the method returns true, then any subsequent call to (any
- form of) MCFGetPi() will return a different dual solution, and MCFGetRC()
- [see below] will return the corresponding reduced costs. The new solution
- need not be optimal (although, ideally, it has to be "good); this can be
- checked by comparing its objective function value, that will be returned
- by a call to MCFGetDFO() [see below].
-
- Any subsequent call of HaveNewPi() that returns true produces a new
- solution, until the first that returns false; from then on, no new
- solutions will be generated until something changes in the problem's
- data.
-
- Note that a default implementation of HaveNewPi() is provided which is
- good for those solvers that only produce one optimal dual solution. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading reduced costs
- @{ */
-
- virtual void MCFGetRC( CRow CR , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Write the reduced costs corresponding to the current dual solution in
- RC[]. If nms == NULL, the reduced cost of arc `i' (i in 0 .. m - 1) is
- written in RC[ i ]; if nms != NULL, it must point to a vector of indices
- in 0 .. m - 1 (ordered in increasing sense and Inf<Index>()-terminated),
- and the reduced cost of arc nms[ i ] is written in RC[ i ]. Note that,
- unlike MCFGetX() above, nms is an *input* of the method.
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to arcs which are *both* in nms[] and whose index is
- in the correct range are returned.
-
- @note the output of MCFGetRC() will change after any call to HaveNewPi()
- [see above] which returns true. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cCRow MCFGetRC( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal data structure containing the
- reduced costs. Since this may *not always be available*, depending on the
- implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method.
-
- @note the output of MCFGetRC() will change after any call to HaveNewPi()
- [see above] which returns true. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual CNumber MCFGetRC( cIndex i ) = 0;
-
-/**< Return the reduced cost of the i-th arc. This information should be
- cheapily available in most implementations.
-
- @note the output of MCFGetRC() will change after any call to HaveNewPi()
- [see above] which returns true. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading the objective function value
- @{ */
-
- virtual FONumber MCFGetFO( void ) = 0;
-
-/**< Return the objective function value of the primal solution currently
- returned by MCFGetX().
-
- If MCFGetStatus() == kOK, this is guaranteed to be the optimal objective
- function value of the problem (to within the optimality tolerances), but
- only prior to any call to HaveNewX() that returns true. MCFGetFO()
- typically returns Inf<FONumber>() if MCFGetStatus() == kUnfeasible and
- - Inf<FONumber>() if MCFGetStatus() == kUnbounded. If MCFGetStatus() ==
- kStopped and MCFGetFO() returns a finite value, it must be an upper bound
- on the optimal objective function value (typically, the objective function
- value of one primal feasible solution). */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual FONumber MCFGetDFO( void )
- {
- switch( MCFGetStatus() ) {
- case( kUnSolved ):
- case( kStopped ):
- case( kError ): return( -Inf<FONumber>() );
- default: return( MCFGetFO() );
- }
- }
-
-/**< Return the objective function value of the dual solution currently
- returned by MCFGetPi() / MCFGetRC(). This value (possibly) changes after
- any call to HaveNewPi() that returns true. The relations between
- MCFGetStatus() and MCFGetDFO() are analogous to these of MCFGetFO(),
- except that a finite value corresponding to kStopped must be a lower
- bound on the optimal objective function value (typically, the objective
- function value one dual feasible solution).
-
- A default implementation is provided for MCFGetDFO(), which is good for
- MCF solvers where the primal and dual optimal solution values always
- are identical (except if the problem is unfeasible/unbounded). */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Getting unfeasibility certificate
- @{ */
-
- virtual FNumber MCFGetUnfCut( Index_Set Cut )
- {
- *Cut = Inf<Index>();
- return( 0 );
- }
-
-/**< Return an unfeasibility certificate. In an unfeasible MCF problem,
- unfeasibility can always be reduced to the existence of a cut (subset of
- nodes of the graph) such as either:
-
- - the inverse of the deficit of the cut (the sum of all the deficits of the
- nodes in the cut) is larger than the forward capacity of the cut (sum of
- the capacities of forward arcs in the cut); that is, the nodes in the
- cut globally produce more flow than can be routed to sinks outside the
- cut;
-
- - the deficit of the cut is larger than the backward capacity of the cut
- (sum of the capacities of backward arcs in the cut); that is, the nodes
- in the cut globally require more flow than can be routed to them from
- sources outside the cut.
-
- When detecting unfeasibility, MCF solvers are typically capable to provide
- one such cut. This information can be useful - typically, the only way to
- make the problem feasible is to increase the capacity of at least one of
- the forward/backward arcs of the cut -, and this method is provided for
- getting it. It can be called only if MCFGetStatus() == kUnfeasible, and
- should write in Cut the set of names of nodes in the unfeasible cut (note
- that node names depend on USENAME0), Inf<Index>()-terminated, returning the
- deficit of the cut (which allows to distinguish which of the two cases
- above hold). In general, no special properties can be expected from the
- returned cut, but solvers should be able to provide e.g. "small" cuts.
-
- However, not all solvers may be (easily) capable of providing this
- information; thus, returning 0 (no cut) is allowed, as in the base class
- implementation, to signify that this information is not available. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Getting unboundedness certificate
- @{ */
-
- virtual Index MCFGetUnbCycl( Index_Set Pred , Index_Set ArcPred )
- {
- return( Inf<Index>() );
- }
-
-/**< Return an unboundedness certificate. In an unbounded MCF problem,
- unboundedness can always be reduced to the existence of a directed cycle
- with negative cost and all arcs having infinite capacity.
- When detecting unboundedness, MCF solvers are typically capable to provide
- one such cycle. This information can be useful, and this method is provided
- for getting it. It can be called only if MCFGetStatus() == kUnbounded, and
- writes in Pred[] and ArcPred[], respectively, the node and arc predecessor
- function of the cycle. That is, if node `i' belongs to the cycle then
- `Pred[ i ]' contains the name of the predecessor of `j' of `i' in the cycle
- (note that node names depend on USENAME0), and `ArcPred[ i ]' contains the
- index of the arc joining the two (note that in general there may be
- multiple copies of each arc). Entries of the vectors for nodes not
- belonging to the cycle are in principle undefined, and the name of one node
- belonging to the cycle is returned by the method.
- Note that if there are multiple cycles with negative costs this method
- will return just one of them (finding the cycle with most negative cost
- is an NO-hard problem), although solvers should be able to produce cycles
- with "large negative" cost.
-
- However, not all solvers may be (easily) capable of providing this
- information; thus, returning Inf<Index>() is allowed, as in the base class
- implementation, to signify that this information is not available. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Saving/restoring the state of the solver
- @{ */
-
- virtual MCFStatePtr MCFGetState( void )
- {
- return( NULL );
- }
-
- /**< Save the state of the MCF solver. The MCFClass interface supports the
- notion of saving and restoring the state of the MCF solver, such as the
- current/optimal basis in a simplex solver. The "empty" class MCFState is
- defined as a placeholder for state descriptions.
-
- MCFGetState() creates and returns a pointer to an object of (a proper
- derived class of) class MCFState which describes the current state of the
- MCF solver. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void MCFPutState( MCFStatePtr S ) {}
-
-/**< Restore the solver to the state in which it was when the state `S' was
- created with MCFGetState() [see above].
-
- The typical use of this method is the following: a MCF problem is solved
- and the "optimal state" is set aside. Then the problem changes and it is
- re-solved. Then, the problem has to be changed again to a form that is
- close to the original one but rather different from the second one (think
- of a long backtracking in a Branch & Bound) to which the current "state"
- referes. Then, the old optimal state can be expected to provide a better
- starting point for reoptimization [see ReOptimize() below].
-
- Note, however, that the state is only relative to the optimization process,
- i.e., this operation is meaningless if the data of the problem has changed
- in the meantime. So, if a state has to be used for speeding up
- reoptimization, the following has to be done:
-
- - first, the data of the solver is brought back to *exactly* the same as
- it was at the moment where the state `S' was created (prior than this
- operation a ReOptimize( false ) call is probably advisable);
-
- - then, MCFPutState() is called (and ReOptimize( true ) is called);
-
- - only afterwards the data of the problem is changed to the final state
- and the problem is solved.
-
- A "put state" operation does not "deplete" the state, which can therefore
- be used more than once. Indeed, a state is constructed inside the solver
- for each call to MCFGetState(), but the solver never deletes statuses;
- this has to be done on the outside when they are no longer needed (the
- solver must be "resistent" to deletion of the state at any moment).
-
- Since not all the MCF solvers reoptimize (efficiently enough to make these
- operations worth), an "empty" implementation that does nothing is provided
- by the base class. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Time the code
- @{ */
-
- void TimeMCF( double &t_us , double &t_ss )
- {
- t_us = t_ss = 0;
- if( MCFt )
- MCFt->Read( t_us , t_ss );
- }
-
-/**< Time the code. If called within any of the methods of the class that are
- "actively timed" (this depends on the subclasses), this method returns the
- user and sistem time (in seconds) since the start of that method. If
- methods that are actively timed call other methods that are actively
- timed, TimeMCF() returns the (...) time since the beginning of the
- *outer* actively timed method. If called outside of any actively timed
- method, this method returns the (...) time spent in all the previous
- executions of all the actively timed methods of the class.
-
- Implementing the proper calls to MCFt->Start() and MCFt->Stop() is due to
- derived classes; these should at least be placed at the beginning and at
- the end, respectively, of SolveMCF() and presumably the Chg***() methods,
- that is, at least these methods should be "actively timed". */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- double TimeMCF( void )
- {
- return( MCFt ? MCFt->Read() : 0 );
- }
-
-/**< Like TimeMCF(double,double) [see above], but returns the total time. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Check the solutions
- @{ */
-
- inline void CheckPSol( void );
-
-/**< Check that the primal solution returned by the solver is primal feasible.
- (to within the tolerances set by SetPar(kEps****) [see above], if any). Also,
- check that the objective function value is correct.
-
- This method is implemented by the base class, using the above methods
- for collecting the solutions and the methods of the next section for
- reading the data of the problem; as such, they will work for any derived
- class that properly implements all these methods. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- inline void CheckDSol( void );
-
-/**< Check that the dual solution returned by the solver is dual feasible.
- (to within the tolerances set by SetPar(kEps****) [see above], if any). Also,
- check that the objective function value is correct.
-
- This method is implemented by the base class, using the above methods
- for collecting the solutions and the methods of the next section for
- reading the data of the problem; as such, they will work for any derived
- class that properly implements all these methods. */
-
-/*@} -----------------------------------------------------------------------*/
-/*-------------- METHODS FOR READING THE DATA OF THE PROBLEM ---------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Reading graph size
- @{ */
-
- inline Index MCFnmax( void )
- {
- return( nmax );
- }
-
-/**< Return the maximum number of nodes for this instance of MCFClass.
- The implementation of the method in the base class returns the protected
- fields \c nmax, which is provided for derived classes to hold this
- information. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- inline Index MCFmmax( void )
- {
- return( mmax );
- }
-
-/**< Return the maximum number of arcs for this instance of MCFClass.
- The implementation of the method in the base class returns the protected
- fields \c mmax, which is provided for derived classes to hold this
- information. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- inline Index MCFn( void )
- {
- return( n );
- }
-
-/**< Return the number of nodes in the current graph.
- The implementation of the method in the base class returns the protected
- fields \c n, which is provided for derived classes to hold this
- information. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- inline Index MCFm( void )
- {
- return( m );
- }
-
-/**< Return the number of arcs in the current graph.
- The implementation of the method in the base class returns the protected
- fields \c m, which is provided for derived classes to hold this
- information. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading graph topology
- @{ */
-
- virtual void MCFArcs( Index_Set Startv , Index_Set Endv ,
- cIndex_Set nms = NULL , cIndex strt = 0 ,
- Index stp = Inf<Index>() ) = 0;
-
-/**< Write the starting (tail) and ending (head) nodes of the arcs in Startv[]
- and Endv[]. If nms == NULL, then the information relative to all arcs is
- written into Startv[] and Endv[], otherwise Startv[ i ] and Endv[ i ]
- contain the information relative to arc nms[ i ] (nms[] must be
- Inf<Index>()-terminated).
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to arcs which are *both* in nms and whose index is in
- the correct range are returned.
-
- Startv or Endv can be NULL, meaning that only the other information is
- required.
-
- @note If USENAME0 == 0 then the returned node names will be in the range
- 1 .. n, while if USENAME0 == 1 the returned node names will be in
- the range 0 .. n - 1.
-
- @note If the graph is "dynamic", be careful to use MCFn() e MCFm() to
- properly choose the dimension of nodes and arcs arrays. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual Index MCFSNde( cIndex i ) = 0;
-
-/**< Return the starting (tail) node of the arc `i'.
-
- @note If USENAME0 == 0 then the returned node names will be in the range
- 1 .. n, while if USENAME0 == 1 the returned node names will be in
- the range 0 .. n - 1. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual Index MCFENde( cIndex i ) = 0;
-
-/**< Return the ending (head) node of the arc `i'.
-
- @note If USENAME0 == 0 then the returned node names will be in the range
- 1 .. n, while if USENAME0 == 1 the returned node names will be in
- the range 0 .. n - 1. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cIndex_Set MCFSNdes( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the starting
- (tail) nodes for each arc. Since this may *not always be available*,
- depending on the implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cIndex_Set MCFENdes( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the ending
- (head) nodes for each arc. Since this may *not always be available*,
- depending on the implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading arc costs
- @{ */
-
- virtual void MCFCosts( CRow Costv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Write the arc costs into Costv[]. If nms == NULL, then all the costs are
- written, otherwise Costv[ i ] contains the information relative to arc
- nms[ i ] (nms[] must be Inf<Index>()-terminated).
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to arcs which are *both* in nms and whose index is in
- the correct range are returned. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual CNumber MCFCost( cIndex i ) = 0;
-
-/**< Return the cost of the i-th arc. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cCRow MCFCosts( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the arc
- costs. Since this may *not always be available*, depending on the
- implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void MCFQCoef( CRow Qv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() )
-
-/**< Write the quadratic coefficients of the arc costs into Qv[]. If
- nms == NULL, then all the costs are written, otherwise Costv[ i ] contains
- the information relative to arc nms[ i ] (nms[] must be
- Inf<Index>()-terminated).
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to arcs which are *both* in nms and whose index is in
- the correct range are returned.
-
- Note that the method is *not* pure virtual: an implementation is provided
- for "pure linear" MCF solvers that only work with all zero quadratic
- coefficients. */
- {
- if( nms ) {
- while( *nms < strt )
- nms++;
-
- for( Index h ; ( h = *(nms++) ) < stp ; )
- *(Qv++) = 0;
- }
- else {
- if( stp > m )
- stp = m;
-
- while( stp-- > strt )
- *(Qv++) = 0;
- }
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual CNumber MCFQCoef( cIndex i )
-
-/**< Return the quadratic coefficients of the cost of the i-th arc. Note that
- the method is *not* pure virtual: an implementation is provided for "pure
- linear" MCF solvers that only work with all zero quadratic coefficients. */
- {
- return( 0 );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cCRow MCFQCoef( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the arc
- costs. Since this may *not always be available*, depending on the
- implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready (such as "pure linear" MCF solvers that only
- work with all zero quadratic coefficients) does not need to implement the
- method. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading arc capacities
- @{ */
-
- virtual void MCFUCaps( FRow UCapv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Write the arc capacities into UCapv[]. If nms == NULL, then all the
- capacities are written, otherwise UCapv[ i ] contains the information
- relative to arc nms[ i ] (nms[] must be Inf<Index>()-terminated).
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the arcs `i' with strt <= i < min( MCFm() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to arcs which are *both* in nms and whose index is in
- the correct range are returned. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual FNumber MCFUCap( cIndex i ) = 0;
-
-/**< Return the capacity of the i-th arc. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cFRow MCFUCaps( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the arc
- capacities. Since this may *not always be available*, depending on the
- implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Reading node deficits
- @{ */
-
- virtual void MCFDfcts( FRow Dfctv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Write the node deficits into Dfctv[]. If nms == NULL, then all the
- defcits are written, otherwise Dfctvv[ i ] contains the information
- relative to node nms[ i ] (nms[] must be Inf<Index>()-terminated).
-
- The parameters `strt' and `stp' allow to restrict the output of the method
- to all and only the nodes `i' with strt <= i < min( MCFn() , stp ). `strt'
- and `stp' work in "&&" with nms; that is, if nms != NULL then only the
- values corresponding to nodes which are *both* in nms and whose index is in
- the correct range are returned.
-
- @note Here node "names" (strt and stp, those contained in nms[] or `i'
- in MCFDfct()) go from 0 to n - 1, regardless to the value of
- USENAME0; hence, if USENAME0 == 0 then the first node is "named 1",
- but its deficit is returned by MCFDfcts( Dfctv , NULL , 0 , 1 ). */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual FNumber MCFDfct( cIndex i ) = 0;
-
-/**< Return the deficit of the i-th node.
-
- @note Here node "names" go from 0 to n - 1, regardless to the value of
- USENAME0; hence, if USENAME0 == 0 then the first node is "named 1",
- but its deficit is returned by MCFDfct( 0 ). */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual cFRow MCFDfcts( void )
- {
- return( NULL );
- }
-
-/**< Return a read-only pointer to an internal vector containing the node
- deficits. Since this may *not always be available*, depending on the
- implementation, this method can (uniformly) return NULL.
- This is done by the base class already, so a derived class that does not
- have the information ready does not need to implement the method.
-
- @note Here node "names" go from 0 to n - 1, regardless to the value of
- USENAME0; hence, if USENAME0 == 0 then the first node is "named 1",
- but its deficit is contained in MCFDfcts()[ 0 ]. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Write problem to file
- @{ */
-
- virtual inline void WriteMCF( ostream &oStrm , int frmt = 0 );
-
-/**< Write the current MCF problem to an ostream. This may be useful e.g. for
- debugging purposes.
-
- The base MCFClass class provides output in two different formats, depending
- on the value of the parameter frmt:
-
- - kDimacs the problem is written in DIMACS standard format, read by most
- MCF codes available;
-
- - kMPS the problem is written in the "modern version" (tab-separated)
- of the MPS format, read by most LP/MIP solvers;
-
- - kFWMPS the problem is written in the "old version" (fixed width
- fields) of the MPS format; this is read by most LP/MIP
- solvers, but some codes still require the old format.
-
- The implementation of WriteMCF() in the base class uses all the above
- methods for reading the data; as such it will work for any derived class
- that properly implements this part of the interface, but it may not be
- very efficient. Thus, the method is virtual to allow the derived classes
- to either re-implement WriteMCF() for the above two formats in a more
- efficient way, and/or to extend it to support other solver-specific
- formats.
-
- @note None of the above two formats supports quadratic MCFs, so if
- nonzero quadratic coefficients are present, they are just ignored.
- */
-
-/*@} -----------------------------------------------------------------------*/
-/*------------- METHODS FOR ADDING / REMOVING / CHANGING DATA --------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Changing the costs
- @{ */
-
- virtual void ChgCosts( cCRow NCost , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Change the arc costs. In particular, change the costs that are:
-
- - listed in into the vector of indices `nms' (ordered in increasing
- sense and Inf<Index>()-terminated),
-
- - *and* whose name belongs to the interval [`strt', `stp').
-
- That is, if strt <= nms[ i ] < stp, then the nms[ i ]-th cost will be
- changed to NCost[ i ]. If nms == NULL (as the default), *all* the
- entries in the given range will be changed; if stp > MCFm(), then the
- smaller bound is used.
-
- @note changing the costs of arcs that *do not exist* is *not allowed*;
- only arcs which have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChgCost( Index arc , cCNumber NCost ) = 0;
-
-/**< Change the cost of the i-th arc.
-
- @note changing the costs of arcs that *do not exist* is *not allowed*;
- only arcs which have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChgQCoef( cCRow NQCoef = NULL , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() )
-
-/**< Change the quadratic coefficients of the arc costs. In particular,
- change the coefficients that are:
-
- - listed in into the vector of indices `nms' (ordered in increasing
- sense and Inf<Index>()-terminated),
-
- - *and* whose name belongs to the interval [`strt', `stp').
-
- That is, if strt <= nms[ i ] < stp, then the nms[ i ]-th cost will be
- changed to NCost[ i ]. If nms == NULL (as the default), *all* the
- entries in the given range will be changed; if stp > MCFm(), then the
- smaller bound is used. If NQCoef == NULL, all the specified coefficients
- are set to zero.
-
- Note that the method is *not* pure virtual: an implementation is provided
- for "pure linear" MCF solvers that only work with all zero quadratic
- coefficients.
-
- @note changing the costs of arcs that *do not exist* is *not allowed*;
- only arcs which have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
- {
- assert(NQCoef);
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChgQCoef( Index arc , cCNumber NQCoef )
-
-/**< Change the quadratic coefficient of the cost of the i-th arc.
-
- Note that the method is *not* pure virtual: an implementation is provided
- for "pure linear" MCF solvers that only work with all zero quadratic
- coefficients.
-
- @note changing the costs of arcs that *do not exist* is *not allowed*;
- only arcs which have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
- {
- assert(NQCoef);
- }
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Changing the capacities
- @{ */
-
- virtual void ChgUCaps( cFRow NCap , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Change the arc capacities. In particular, change the capacities that are:
-
- - listed in into the vector of indices `nms' (ordered in increasing sense
- and Inf<Index>()-terminated),
-
- - *and* whose name belongs to the interval [`strt', `stp').
-
- That is, if strt <= nms[ i ] < stp, then the nms[ i ]-th capacity will be
- changed to NCap[ i ]. If nms == NULL (as the default), *all* the
- entries in the given range will be changed; if stp > MCFm(), then the
- smaller bound is used.
-
- @note changing the capacities of arcs that *do not exist* is *not allowed*;
- only arcs that have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChgUCap( Index arc , cFNumber NCap ) = 0;
-
-/**< Change the capacity of the i-th arc.
-
- @note changing the capacities of arcs that *do not exist* is *not allowed*;
- only arcs that have not been closed/deleted [see CloseArc() /
- DelArc() below and LoadNet() above about C_INF costs] can be
- touched with these methods. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Changing the deficits
- @{ */
-
- virtual void ChgDfcts( cFRow NDfct , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) = 0;
-
-/**< Change the node deficits. In particular, change the deficits that are:
-
- - listed in into the vector of indices `nms' (ordered in increasing sense
- and Inf<Index>()-terminated),
-
- - *and* whose name belongs to the interval [`strt', `stp').
-
- That is, if strt <= nms[ i ] < stp, then the nms[ i ]-th deficit will be
- changed to NDfct[ i ]. If nms == NULL (as the default), *all* the
- entries in the given range will be changed; if stp > MCFn(), then the
- smaller bound is used.
-
- Note that, in ChgDfcts(), node "names" (strt, stp or those contained
- in nms[]) go from 0 to n - 1, regardless to the value of USENAME0; hence,
- if USENAME0 == 0 then the first node is "named 1", but its deficit can be
- changed e.g. with ChgDfcts( &new_deficit , NULL , 0 , 1 ).
-
- @note changing the capacities of nodes that *do not exist* is *not
- allowed*; only nodes that have not been deleted [see DelNode()
- below] can be touched with these methods. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChgDfct( Index node , cFNumber NDfct ) = 0;
-
-/**< Change the deficit of the i-th node.
-
- Note that, in ChgDfct[s](), node "names" (i, strt/ stp or those contained
- in nms[]) go from 0 to n - 1, regardless to the value of USENAME0; hence,
- if USENAME0 == 0 then the first node is "named 1", but its deficit can be
- changed e.g. with ChgDfcts( 0 , new_deficit ).
-
- @note changing the capacities of nodes that *do not exist* is *not
- allowed*; only nodes that have not been deleted [see DelNode()
- below] can be touched with these methods. */
-
-/*@} -----------------------------------------------------------------------*/
-/** @name Changing graph topology
- @{ */
-
- virtual void CloseArc( cIndex name ) = 0;
-
-/**< "Close" the arc `name'. Although all the associated information (name,
- cost, capacity, end and start node) is kept, the arc is removed from the
- problem until OpenArc( i ) [see below] is called.
-
- "closed" arcs always have 0 flow, but are otherwise counted as any other
- arc; for instance, MCFm() does *not* decrease as an effect of a call to
- CloseArc(). How this closure is implemented is solver-specific. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual bool IsClosedArc( cIndex name ) = 0;
-
-/**< IsClosedArc() returns true if and only if the arc `name' is closed. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void DelNode( cIndex name ) = 0;
-
-/**< Delete the node `name'.
-
- For any value of `name', all incident arcs to that node are closed [see
- CloseArc() above] (*not* Deleted, see DelArc() below) and the deficit is
- set to zero.
-
- Il furthermore `name' is the last node, the number of nodes as reported by
- MCFn() is reduced by at least one, until the n-th node is not a deleted
- one. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void OpenArc( cIndex name ) = 0;
-
-/**< Restore the previously closed arc `name'. It is an error to open an arc
- that has not been previously closed. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual Index AddNode( cFNumber aDfct ) = 0;
-
-/**< Add a new node with deficit aDfct, returning its name. Inf<Index>() is
- returned if there is no room for a new node. Remember that the node names
- are either { 0 .. nmax - 1 } or { 1 .. nmax }, depending on the value of
- USENAME0. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void ChangeArc( cIndex name , cIndex nSN = Inf<Index>() ,
- cIndex nEN = Inf<Index>() ) = 0;
-
-/**< Change the starting and/or ending node of arc `name' to nSN and nEN.
- Each parameter being Inf<Index>() means to leave the previous starting or
- ending node untouched. When this method is called `name' can be either the
- name of a "normal" arc or that of a "closed" arc [see CloseArc() above]:
- in the latter case, at the end of ChangeArc() the arc is *still closed*,
- and it remains so until OpenArc( name ) [see above] is called. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual void DelArc( cIndex name ) = 0;
-
-/**< Delete the arc `name'. Unlike "closed" arcs, all the information
- associated with a deleted arc is lost and `name' is made available as a
- name for new arcs to be created with AddArc() [see below].
-
- Il furthermore `name' is the last arc, the number of arcs as reported by
- MCFm() is reduced by at least one, until the m-th arc is not a deleted
- one. Otherwise, the flow on the arc is always ensured to be 0. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual bool IsDeletedArc( cIndex name ) = 0;
-
-/**< Return true if and only if the arc `name' is deleted. It should only
- be called with name < MCFm(), as every other arc is deleted by
- definition. */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
-
- virtual Index AddArc( cIndex Start , cIndex End , cFNumber aU ,
- cCNumber aC ) = 0;
-
-/**< Add the new arc ( Start , End ) with cost aC and capacity aU,
- returning its name. Inf<Index>() is returned if there is no room for a new
- arc. Remember that arc names go from 0 to mmax - 1. */
-
-/*@} -----------------------------------------------------------------------*/
-/*------------------------------ DESTRUCTOR --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Destructor
- @{ */
-
- virtual ~MCFClass()
- {
- delete MCFt;
- }
-
-/**< Destructor of the class. The implementation in the base class only
- deletes the MCFt field. It is virtual, as it should be. */
-
-/*@} -----------------------------------------------------------------------*/
-/*-------------------- PROTECTED PART OF THE CLASS -------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-- --*/
-/*-- The following fields are believed to be general enough to make it --*/
-/*-- worth adding them to the abstract base class. --*/
-/*-- --*/
-/*--------------------------------------------------------------------------*/
-
- protected:
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------- PROTECTED METHODS ----------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-------------------------- MANAGING COMPARISONS --------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @name Managing comparisons.
- The following methods are provided for making it easier to perform
- comparisons, with and without tolerances.
- @{ */
-
-/*--------------------------------------------------------------------------*/
-/** true if flow x is equal to zero (possibly considering tolerances). */
-
- template<class T>
- inline bool ETZ( T x , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x == 0 );
- else
- return( (x <= eps ) && ( x >= -eps ) );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is greater than zero (possibly considering tolerances). */
-
- template<class T>
- inline bool GTZ( T x , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x > 0 );
- else
- return( x > eps );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is greater than or equal to zero (possibly considering
- tolerances). */
-
- template<class T>
- inline bool GEZ( T x , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x >= 0 );
- else
- return( x >= - eps );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is less than zero (possibly considering tolerances). */
-
- template<class T>
- inline bool LTZ( T x , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x < 0 );
- else
- return( x < - eps );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is less than or equal to zero (possibly considering
- tolerances). */
-
- template<class T>
- inline bool LEZ( T x , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x <= 0 );
- else
- return( x <= eps );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is greater than flow y (possibly considering tolerances).
- */
-
- template<class T>
- inline bool GT( T x , T y , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x > y );
- else
- return( x > y + eps );
- }
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-/** true if flow x is less than flow y (possibly considering tolerances). */
-
- template<class T>
- inline bool LT( T x , T y , const T eps )
- {
- if( numeric_limits<T>::is_integer )
- return( x < y );
- else
- return( x < y - eps );
- }
-
-/*@} -----------------------------------------------------------------------*/
-/*---------------------- PROTECTED DATA STRUCTURES -------------------------*/
-/*--------------------------------------------------------------------------*/
-
- Index n; ///< total number of nodes
- Index nmax; ///< maximum number of nodes
-
- Index m; ///< total number of arcs
- Index mmax; ///< maximum number of arcs
-
- int status; ///< return status, see the comments to MCFGetStatus() above.
- ///< Note that the variable is defined int to allow derived
- ///< classes to return their own specialized status codes
- bool Senstv; ///< true <=> the latest optimal solution should be exploited
-
- OPTtimers *MCFt; ///< timer for performances evaluation
-
- FNumber EpsFlw; ///< precision for comparing arc flows / capacities
- FNumber EpsDfct; ///< precision for comparing node deficits
- CNumber EpsCst; ///< precision for comparing arc costs
-
- double MaxTime; ///< max time (in seconds) in which MCF Solver can find
- ///< an optimal solution (0 = no limits)
- int MaxIter; ///< max number of iterations in which MCF Solver can find
- ///< an optimal solution (0 = no limits)
-
-/*--------------------------------------------------------------------------*/
-
- }; // end( class MCFClass )
-
-/*@} -----------------------------------------------------------------------*/
-/*------------------- inline methods implementation ------------------------*/
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::LoadDMX( istream &DMXs , bool IsQuad )
-{
- // read first non-comment line - - - - - - - - - - - - - - - - - - - - - - -
- // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- char c;
- for(;;) {
- if( ! ( DMXs >> c ) )
- throw( MCFException( "LoadDMX: error reading the input stream" ) );
-
- if( c != 'c' ) // if it's not a comment
- break;
-
- do // skip the entire line
- if( ! DMXs.get( c ) )
- throw( MCFException( "LoadDMX: error reading the input stream" ) );
- while( c != '\n' );
- }
-
- if( c != 'p' )
- throw( MCFException( "LoadDMX: format error in the input stream" ) );
-
- char buf[ 3 ]; // free space
- DMXs >> buf; // skip (in has to be "min")
-
- Index tn;
- if( ! ( DMXs >> tn ) )
- throw( MCFException( "LoadDMX: error reading number of nodes" ) );
-
- Index tm;
- if( ! ( DMXs >> tm ) )
- throw( MCFException( "LoadDMX: error reading number of arcs" ) );
-
- // allocate memory - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
- // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- FRow tU = new FNumber[ tm ]; // arc upper capacities
- FRow tDfct = new FNumber[ tn ]; // node deficits
- Index_Set tStartn = new Index[ tm ]; // arc start nodes
- Index_Set tEndn = new Index[ tm ]; // arc end nodes
- CRow tC = new CNumber[ tm ]; // arc costs
- CRow tQ = NULL;
- if( IsQuad )
- tQ = new CNumber[ tm ]; // arc quadratic costs
-
-
- for( Index i = 0 ; i < tn ; ) // all deficits are 0
- tDfct[ i++ ] = 0; // unless otherwise stated
-
- // read problem data - - - - - - - - - - - - - - - - - - - - - - - - - - - -
- // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- Index i = 0; // arc counter
- for(;;) {
- if( ! ( DMXs >> c ) )
- break;
-
- switch( c ) {
- case( 'c' ): // comment line- - - - - - - - - - - - - - - - - - - - - - -
- do
- if( ! DMXs.get( c ) )
- throw( MCFException( "LoadDMX: error reading the input stream" ) );
- while( c != '\n' );
- break;
-
- case( 'n' ): // description of a node - - - - - - - - - - - - - - - - - -
- Index j;
- if( ! ( DMXs >> j ) )
- throw( MCFException( "LoadDMX: error reading node name" ) );
-
- if( ( j < 1 ) || ( j > tn ) )
- throw( MCFException( "LoadDMX: invalid node name" ) );
-
- FNumber Dfctj;
- if( ! ( DMXs >> Dfctj ) )
- throw( MCFException( "LoadDMX: error reading deficit" ) );
-
- tDfct[ j - 1 ] -= Dfctj;
- break;
-
- case( 'a' ): // description of an arc - - - - - - - - - - - - - - - - - -
- if( i == tm )
- throw( MCFException( "LoadDMX: too many arc descriptors" ) );
-
- if( ! ( DMXs >> tStartn[ i ] ) )
- throw( MCFException( "LoadDMX: error reading start node" ) );
-
- if( ( tStartn[ i ] < 1 ) || ( tStartn[ i ] > tn ) )
- throw( MCFException( "LoadDMX: invalid start node" ) );
-
- if( ! ( DMXs >> tEndn[ i ] ) )
- throw( MCFException( "LoadDMX: error reading end node" ) );
-
- if( ( tEndn[ i ] < 1 ) || ( tEndn[ i ] > tn ) )
- throw( MCFException( "LoadDMX: invalid end node" ) );
-
- if( tStartn[ i ] == tEndn[ i ] )
- throw( MCFException( "LoadDMX: self-loops not permitted" ) );
-
- FNumber LB;
- if( ! ( DMXs >> LB ) )
- throw( MCFException( "LoadDMX: error reading lower bound" ) );
-
- if( ! ( DMXs >> tU[ i ] ) )
- throw( MCFException( "LoadDMX: error reading upper bound" ) );
-
- if( ! ( DMXs >> tC[ i ] ) )
- throw( MCFException( "LoadDMX: error reading arc cost" ) );
-
- if( tQ ) {
- if( ! ( DMXs >> tQ[ i ] ) )
- throw( MCFException( "LoadDMX: error reading arc quadratic cost" ) );
-
- if( tQ[ i ] < 0 )
- throw( MCFException( "LoadDMX: negative arc quadratic cost" ) );
- }
-
- if( tU[ i ] < LB )
- throw( MCFException( "LoadDMX: lower bound > upper bound" ) );
-
- tU[ i ] -= LB;
- tDfct[ tStartn[ i ] - 1 ] += LB;
- tDfct[ tEndn[ i ] - 1 ] -= LB;
- #if( USENAME0 )
- tStartn[ i ]--; // in the DIMACS format, node names start from 1
- tEndn[ i ]--;
- #endif
- i++;
- break;
-
- default: // invalid code- - - - - - - - - - - - - - - - - - - - - - - - -
- throw( MCFException( "LoadDMX: invalid code" ) );
-
- } // end( switch( c ) )
- } // end( for( ever ) )
-
- if( i < tm )
- throw( MCFException( "LoadDMX: too few arc descriptors" ) );
-
- // call LoadNet- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
- // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- LoadNet( tn , tm , tn , tm , tU , tC , tDfct , tStartn , tEndn );
-
- // then pass quadratic costs, if any
-
- if( tQ )
- ChgQCoef( tQ );
-
- // delete the original data structures - - - - - - - - - - - - - - - - - - -
- // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- delete[] tQ;
- delete[] tC;
- delete[] tEndn;
- delete[] tStartn;
- delete[] tDfct;
- delete[] tU;
-
- } // end( MCFClass::LoadDMX )
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::SetPar( int par , int val )
-{
- switch( par ) {
- case( kMaxIter ):
- MaxIter = val;
- break;
-
- case( kReopt ):
- Senstv = (val == kYes);
- break;
-
- default:
- throw( MCFException( "Error using SetPar: unknown parameter" ) );
- }
- }
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::SetPar( int par, double val )
-{
- switch( par ) {
- case( kEpsFlw ):
- if( EpsFlw != FNumber( val ) ) {
- EpsFlw = FNumber( val );
- EpsDfct = EpsFlw * ( nmax ? nmax : 100 );
- status = kUnSolved;
- }
- break;
-
- case( kEpsDfct ):
- if( EpsDfct != FNumber( val ) ) {
- EpsDfct = FNumber( val );
- status = kUnSolved;
- }
- break;
-
- case( kEpsCst ):
- if( EpsCst != CNumber( val ) ) {
- EpsCst = CNumber( val );
- status = kUnSolved;
- }
- break;
-
- case( kMaxTime ):
- MaxTime = val;
- break;
-
- default:
- throw( MCFException( "Error using SetPar: unknown parameter" ) );
- }
- }
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::GetPar( int par, int &val )
-{
- switch( par ) {
- case( kMaxIter ):
- val = MaxIter;
- break;
-
- case( kReopt ):
- if( Senstv ) val = kYes;
- else val = kNo;
- break;
-
- default:
- throw( MCFException( "Error using GetPar: unknown parameter" ) );
- }
- }
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::GetPar( int par , double &val )
-{
- switch( par ) {
- case( kEpsFlw ):
- val = double( EpsFlw );
- break;
-
- case( kEpsDfct ):
- val = double( EpsDfct );
- break;
-
- case( kEpsCst ):
- val = double( EpsCst );
- break;
-
- case( kMaxTime ):
- val = MaxTime;
- break;
-
- default:
- throw( MCFException( "Error using GetPar: unknown parameter" ) );
- }
- }
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::CheckPSol( void )
-{
- FRow tB = new FNumber[ MCFn() ];
- MCFDfcts( tB );
- #if( ! USENAME0 )
- tB--;
- #endif
-
- FRow tX = new FNumber[ MCFm() ];
- MCFGetX( tX );
- FONumber CX = 0;
- for( Index i = 0 ; i < MCFm() ; i++ )
- if( GTZ( tX[ i ] , EpsFlw ) ) {
- if( IsClosedArc( i ) )
- throw( MCFException( "Closed arc with nonzero flow (CheckPSol)" ) );
-
- if( GT( tX[ i ] , MCFUCap( i ) , EpsFlw ) )
- throw( MCFException( "Arc flow exceeds capacity (CheckPSol)" ) );
-
- if( LTZ( tX[ i ] , EpsFlw ) )
- throw( MCFException( "Arc flow is negative (CheckPSol)" ) );
-
- CX += ( FONumber( MCFCost( i ) ) +
- FONumber( MCFQCoef( i ) ) * FONumber( tX[ i ] ) / 2 )
- * FONumber( tX[ i ] );
- tB[ MCFSNde( i ) ] += tX[ i ];
- tB[ MCFENde( i ) ] -= tX[ i ];
- }
- delete[] tX;
- #if( ! USENAME0 )
- tB++;
- #endif
-
- for( Index i = 0 ; i < MCFn() ; i++ )
- if( ! ETZ( tB[ i ] , EpsDfct ) )
- throw( MCFException( "Node is not balanced (CheckPSol)" ) );
-
- delete[] tB;
- CX -= MCFGetFO();
- if( ( CX >= 0 ? CX : - CX ) > EpsCst * MCFn() )
- throw( MCFException( "Objective function value is wrong (CheckPSol)" ) );
-
- } // end( CheckPSol )
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::CheckDSol( void )
-{
- CRow tPi = new CNumber[ MCFn() ];
- MCFGetPi( tPi );
-
- FONumber BY = 0;
- for( Index i = 0 ; i < MCFn() ; i++ )
- BY += FONumber( tPi[ i ] ) * FONumber( MCFDfct( i ) );
-
- CRow tRC = new CNumber[ MCFm() ];
- MCFGetRC( tRC );
-
- FRow tX = new FNumber[ MCFm() ];
- MCFGetX( tX );
-
- #if( ! USENAME0 )
- tPi--;
- #endif
-
- FONumber QdrtcCmpnt = 0; // the quadratic part of the objective
- for( Index i = 0 ; i < MCFm() ; i++ ) {
- if( IsClosedArc( i ) )
- continue;
-
- cFONumber tXi = FONumber( tX[ i ] );
- cCNumber QCi = CNumber( MCFQCoef( i ) * tXi );
- QdrtcCmpnt += QCi * tXi;
- CNumber RCi = MCFCost( i ) + QCi
- + tPi[ MCFSNde( i ) ] - tPi[ MCFENde( i ) ];
- RCi -= tRC[ i ];
-
- if( ! ETZ( RCi , EpsCst ) )
- throw( MCFException( "Reduced Cost value is wrong (CheckDSol)" ) );
-
- if( LTZ( tRC[ i ] , EpsCst ) ) {
- BY += FONumber( tRC[ i ] ) * FONumber( MCFUCap( i ) );
-
- if( LT( tX[ i ] , MCFUCap( i ) , EpsFlw ) )
- throw( MCFException( "Csomplementary Slackness violated (CheckDSol)" ) );
- }
- else
- if( GTZ( tRC[ i ] , EpsCst ) && GTZ( tX[ i ] , EpsFlw ) )
- throw( MCFException( "Complementary Slackness violated (CheckDSol)" ) );
-
- } // end( for( i ) )
-
- delete[] tX;
- delete[] tRC;
- #if( ! USENAME0 )
- tPi++;
- #endif
- delete[] tPi;
-
- BY -= ( MCFGetFO() + QdrtcCmpnt / 2 );
- if( ( BY >= 0 ? BY : - BY ) > EpsCst * MCFn() )
- throw( MCFException( "Objective function value is wrong (CheckDSol)" ) );
-
- } // end( CheckDSol )
-
-/*--------------------------------------------------------------------------*/
-
-inline void MCFClass::WriteMCF( ostream &oStrm , int frmt )
-{
- if( ( ! numeric_limits<FNumber>::is_integer ) ||
- ( ! numeric_limits<CNumber>::is_integer ) )
- oStrm.precision( 12 );
-
- switch( frmt ) {
- case( kDimacs ): // DIMACS standard format - - - - - - - - - - - - - - - -
- //- - - - - - - - - - - - - - - - - - - - - - - - - - - -
- case( kQDimacs ): // ... linear or quadratic
-
- // compute and output preamble and size- - - - - - - - - - - - - - - - - -
-
- oStrm << "p min " << MCFn() << " ";
- {
- Index tm = MCFm();
- for( Index i = MCFm() ; i-- ; )
- if( IsClosedArc( i ) || IsDeletedArc( i ) )
- tm--;
-
- oStrm << tm << endl;
- }
-
- // output arc information- - - - - - - - - - - - - - - - - - - - - - - - -
-
- for( Index i = 0 ; i < MCFm() ; i++ )
- if( ( ! IsClosedArc( i ) ) && ( ! IsDeletedArc( i ) ) ) {
- oStrm << "a\t";
- #if( USENAME0 )
- oStrm << MCFSNde( i ) + 1 << "\t" << MCFENde( i ) + 1 << "\t";
- #else
- oStrm << MCFSNde( i ) << "\t" << MCFENde( i ) << "\t";
- #endif
- oStrm << "0\t" << MCFUCap( i ) << "\t" << MCFCost( i );
-
- if( frmt == kQDimacs )
- oStrm << "\t" << MCFQCoef( i );
-
- oStrm << endl;
- }
-
- // output node information - - - - - - - - - - - - - - - - - - - - - - - -
-
- for( Index i = 0 ; i < MCFn() ; ) {
- cFNumber Dfcti = MCFDfct( i++ );
- if( Dfcti )
- oStrm << "n\t" << i << "\t" << - Dfcti << endl;
- }
-
- break;
-
- case( kMPS ): // MPS format(s)- - - - - - - - - - - - - - - - - - - - -
- case( kFWMPS ): //- - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- // writing problem name- - - - - - - - - - - - - - - - - - - - - - - - - -
-
- oStrm << "NAME MCF" << endl;
-
- // writing ROWS section: flow balance constraints- - - - - - - - - - - - -
-
- oStrm << "ROWS" << endl << " N obj" << endl;
-
- for( Index i = 0 ; i < MCFn() ; )
- oStrm << " E c" << i++ << endl;
-
- // writing COLUMNS section - - - - - - - - - - - - - - - - - - - - - - - -
-
- oStrm << "COLUMNS" << endl;
-
- if( frmt == kMPS ) { // "modern" MPS
- for( Index i = 0 ; i < MCFm() ; i++ )
- if( ( ! IsClosedArc( i ) ) && ( ! IsDeletedArc( i ) ) ) {
- oStrm << " x" << i << "\tobj\t" << MCFCost( i ) << "\tc";
- #if( USENAME0 )
- oStrm << MCFSNde( i ) << "\t-1" << endl << " x" << i << "\tc"
- << MCFENde( i ) << "\t1" << endl;
- #else
- oStrm << MCFSNde( i ) - 1 << "\t-1" << endl << " x" << i << "\tc"
- << MCFENde( i ) - 1 << "\t1" << endl;
- #endif
- }
- }
- else // "old" MPS
- for( Index i = 0 ; i < MCFm() ; i++ ) {
- ostringstream myField;
- myField << "x" << i << ends;
- oStrm << " " << setw( 8 ) << left << myField.str()
- << " " << setw( 8 ) << left << "obj"
- << " " << setw( 12 ) << left << MCFCost( i );
-
- myField.seekp( 0 );
- myField << "c" << MCFSNde( i ) - ( 1 - USENAME0 ) << ends;
-
- oStrm << " " << setw( 8 ) << left << myField.str()
- << " " << setw( 12 ) << left << -1 << endl;
-
- myField.seekp( 0 );
- myField << "x" << i << ends;
-
- oStrm << " " << setw( 8 ) << left << myField.str();
-
- myField.seekp( 0 );
- myField << "c" << MCFENde( i ) - ( 1 - USENAME0 ) << ends;
-
- oStrm << " " << setw( 8 ) << left << myField.str() << endl;
- }
-
- // writing RHS section: flow balance constraints - - - - - - - - - - - - -
-
- oStrm << "RHS" << endl;
-
- if( frmt == kMPS ) // "modern" MPS
- for( Index i = 0 ; i < MCFn() ; i++ )
- oStrm << " rhs\tc" << i << "\t" << MCFDfct( i ) << endl;
- else // "old" MPS
- for( Index i = 0 ; i < MCFn() ; i++ ) {
- ostringstream myField;
- oStrm << " " << setw( 8 ) << left << "rhs";
- myField << "c" << i << ends;
-
- oStrm << " " << setw( 8 ) << left << myField.str()
- << " " << setw( 12 ) << left << MCFDfct( i ) << endl;
- }
-
- // writing BOUNDS section- - - - - - - - - - - - - - - - - - - - - - - - -
-
- oStrm << "BOUNDS" << endl;
-
- if( frmt == kMPS ) { // "modern" MPS
- for( Index i = 0 ; i < MCFm() ; i++ )
- if( ( ! IsClosedArc( i ) ) && ( ! IsDeletedArc( i ) ) )
- oStrm << " UP BOUND\tx" << i << "\t" << MCFUCap( i ) << endl;
- }
- else // "old" MPS
- for( Index i = 0 ; i < MCFm() ; i++ )
- if( ( ! IsClosedArc( i ) ) && ( ! IsDeletedArc( i ) ) ) {
- ostringstream myField;
- oStrm << " " << setw( 2 ) << left << "UP"
- << " " << setw( 8 ) << left << "BOUND";
-
- myField << "x" << i << ends;
-
- oStrm << " " << setw( 8 ) << left << myField.str()
- << " " << setw( 12 ) << left << MCFUCap( i ) << endl;
- }
-
- oStrm << "ENDATA" << endl;
- break;
-
- default: // unknown format - - - - - - - - - - - - - - - - - - - -
- //- - - - - - - - - - - - - - - - - - - - - - - - - - - -
- oStrm << "Error: unknown format " << frmt << endl;
- }
- } // end( MCFClass::WriteMCF )
-
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
- }; // end( namespace MCFClass_di_unipi_it )
-#endif
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#endif /* MCFClass.h included */
-
-/*--------------------------------------------------------------------------*/
-/*------------------------- End File MCFClass.h ----------------------------*/
-/*--------------------------------------------------------------------------*/
diff --git a/clang/lib/Analysis/MCFSimplex.cpp b/clang/lib/Analysis/MCFSimplex.cpp
index 4b1b75a..3da685d 100644
--- a/clang/lib/Analysis/MCFSimplex.cpp
+++ b/clang/lib/Analysis/MCFSimplex.cpp
@@ -30,7 +30,7 @@
/*------------------------------ INCLUDES ----------------------------------*/
/*--------------------------------------------------------------------------*/
-#include "MCFSimplex.h"
+#include "clang/Analysis/Analyses/IntervalSolver/MCFSimplex.h"
#include <algorithm>
#include <iostream>
@@ -116,9 +116,10 @@ using namespace MCFClass_di_unipi_it;
optimize the update of the potential.
Unfortunately this doesn't work well: for this reason it is set to 0. */
-
+#ifndef throw
// QUICK HACK
#define throw(x) assert(false);
+#endif
/*--------------------------------------------------------------------------*/
/*--------------------------- FUNCTIONS ------------------------------------*/
@@ -410,7 +411,7 @@ void MCFSimplex::SetAlg( bool UsPrml , char WhchPrc )
MemDeAlloc(false);
if( Senstv && ( status != kUnSolved ) ) {
nodePType *node = dummyRootP;
- for( int i = 0 ; i < n ; i++ )
+ for( unsigned int i = 0 ; i < n ; i++ )
node = node->nextInT;
node->nextInT = NULL;
@@ -483,7 +484,7 @@ void MCFSimplex::SetAlg( bool UsPrml , char WhchPrc )
CreateAdditionalDualStructures();
MemDeAlloc(true);
nodeDType *node = dummyRootD;
- for( int i = 0 ; i < n ; i++ )
+ for( unsigned int i = 0 ; i < n ; i++ )
node = node->nextInT;
node->nextInT = NULL;
@@ -513,11 +514,12 @@ void MCFSimplex::SetAlg( bool UsPrml , char WhchPrc )
//#endif
}
#endif
- if( newPricingRule == kFirstEligibleArc )
+ if( newPricingRule == kFirstEligibleArc ) {
if( newUsePrimalSimplex )
arcToStartP = arcsP;
else
arcToStartD = arcsD;
+ }
if( ( nmax && mmax ) && ( newPricingRule == kCandidateListPivot ) )
MemAllocCandidateList();
@@ -941,13 +943,13 @@ void MCFSimplex::ChgCosts( cCRow NCost , cIndex_Set nms ,
for( arcDType *arc = arcsD + strt ; arc < (arcsD + stp) ; arc++ )
arc->cost = *(NCost++);
- if( Senstv && ( status != kUnSolved ) )
+ if( Senstv && ( status != kUnSolved ) ) {
if( usePrimalSimplex )
ComputePotential( dummyRootP );
else {
#if( QUADRATICCOST == 0 )
for( arcDType *arc = arcsD ; arc != stopArcsD ; arc++ )
- if( arc->ident > BASIC )
+ if( arc->ident > BASIC ) {
if( GTZ( ReductCost( arc ) , EpsCst ) ) {
arc->flow = 0;
arc->ident = AT_LOWER;
@@ -956,13 +958,15 @@ void MCFSimplex::ChgCosts( cCRow NCost , cIndex_Set nms ,
arc->flow = arc->upper;
arc->ident = AT_UPPER;
}
+ }
CreateInitialDModifiedBalanceVector();
PostDVisit( dummyRootD );
#endif
}
- else
+ } else {
status = kUnSolved;
+ }
} // end( MCFSimplex::ChgCosts )
@@ -993,8 +997,8 @@ void MCFSimplex::ChgCost( Index arc , cCNumber NCost )
node = ( arcsD + arc )->head;
ComputePotential( dummyRootD );
- for( arcDType *a = arcsD ; a != stopArcsD ; a++)
- if( a->ident > BASIC )
+ for( arcDType *a = arcsD ; a != stopArcsD ; a++)
+ if( a->ident > BASIC ) {
if( GTZ( ReductCost( a ) , EpsCst ) ) {
a->flow = 0;
a->ident = AT_LOWER;
@@ -1003,6 +1007,7 @@ void MCFSimplex::ChgCost( Index arc , cCNumber NCost )
a->flow = a->upper;
a->ident = AT_UPPER;
}
+ }
CreateInitialDModifiedBalanceVector();
PostDVisit( dummyRootD );
@@ -1016,8 +1021,8 @@ void MCFSimplex::ChgCost( Index arc , cCNumber NCost )
/*-------------------------------------------------------------------------*/
-void MCFSimplex::ChgQCoef( cCRow NQCoef , cIndex_Set nms ,
- cIndex strt , Index stp )
+void MCFSimplex::ChgQCoef( cCRow NQCoef , cIndex_Set ,
+ cIndex , Index stp )
{
if( stp > m )
stp = m;
@@ -1057,7 +1062,7 @@ void MCFSimplex::ChgQCoef( cCRow NQCoef , cIndex_Set nms ,
/*-------------------------------------------------------------------------*/
-void MCFSimplex::ChgQCoef( Index arc , cCNumber NQCoef )
+void MCFSimplex::ChgQCoef( Index , cCNumber NQCoef )
{
#if( QUADRATICCOST )
if( arc >= m )
@@ -1362,7 +1367,7 @@ void MCFSimplex::CloseArc( cIndex name )
ComputePotential( dummyRootD );
for( arcDType *a = arcsD ; a != stopArcsD ; a++ )
- if( a->ident > BASIC )
+ if( a->ident > BASIC ) {
if( GTZ( ReductCost( a ) , EpsCst ) ) {
a->flow = 0;
a->ident = AT_LOWER;
@@ -1371,6 +1376,7 @@ void MCFSimplex::CloseArc( cIndex name )
a->flow = a->upper;
a->ident = AT_UPPER;
}
+ }
}
CreateInitialDModifiedBalanceVector();
@@ -2564,11 +2570,12 @@ void MCFSimplex::DualSimplex( void )
while( fine == false ) {
/* If node is the root of subtree T2, Dual Simplex jumps to the node
(if exists) which follows the last node of T2 */
- if( node == h2 )
+ if( node == h2 ) {
if( lastNodeOfT2->nextInT )
node = lastNodeOfT2->nextInT;
else
break;
+ }
// Search arc in the Backward Star of nodes of T1
arcDType *arc = node->firstBs;
@@ -2760,7 +2767,6 @@ void MCFSimplex::DualSimplex( void )
theta = leavingArc->flow - leavingArc->upper;
// Initial value of theta is the infeasibility of the leaving arc
- FNumber t;
nodeDType *k1;
nodeDType *k2;
/* if entering arc is in U, Dual Simplex pushs flow in the cycle
@@ -2905,7 +2911,7 @@ void MCFSimplex::DualSimplex( void )
/*--------------------------------------------------------------------------*/
template<class N, class A>
-void MCFSimplex::UpdateT( A *h , A *k , N *h1 , N *h2 , N *k1 , N *k2 )
+void MCFSimplex::UpdateT( A * , A *k , N * , N *h2 , N *k1 , N *k2 )
{
/* In subtree T2 there is a path from node h2 (deepest node of the leaving
arc h and root of T2) to node k2 (deepest node of the leaving arc h and
@@ -3421,7 +3427,7 @@ MCFSimplex::arcDType* MCFSimplex::RuleDualCandidateListPivot( void )
// Search other arcs to fill the list
do {
nodeDType *node = nodesD + groupPos;
- for( node ; node < stopNodesD ; node += numGroup ) {
+ for( ; node < stopNodesD ; node += numGroup ) {
arcDType *arc = node->enteringTArc;
cFNumber flow = arc->flow;
if( LTZ( flow , EpsFlw ) ) {
@@ -3598,18 +3604,19 @@ void MCFSimplex::PostPVisit( nodePType *r )
// The method controls if "r" is a leaf in T
bool rLeaf = false;
int i = r - nodesP;
- if( r->nextInT )
+ if( r->nextInT ) {
if( ( r->nextInT )->subTreeLevel <= r->subTreeLevel )
rLeaf = true;
else
rLeaf = true;
+ }
- if( rLeaf ) // If "r" is a leaf
+ if( rLeaf ) { // If "r" is a leaf
if( ( r->enteringTArc)->head == r ) // If enteringTArc of "r" goes in "r"
( r->enteringTArc )->flow = modifiedBalance[ i ];
else // If enteringTArc of "r" goes out "r"
( r->enteringTArc )->flow = - modifiedBalance[ i ];
- else { // If "r" isn't a leaf
+ } else { // If "r" isn't a leaf
nodePType *desc = r->nextInT;
// Call PostPVisit for every child of "r"
while( ( desc ) && ( desc->subTreeLevel > r->subTreeLevel ) ) {
@@ -3624,12 +3631,13 @@ void MCFSimplex::PostPVisit( nodePType *r )
desc = desc->nextInT;
}
- if( r != dummyRootP )
+ if( r != dummyRootP ) {
if( ( r->enteringTArc )->head == r ) // If enteringTArc of "r" goes in "r"
( r->enteringTArc )->flow = modifiedBalance[ i ];
else // If enteringTArc of "r" goes out "r"
( r->enteringTArc )->flow = - modifiedBalance[ i ];
}
+ }
}
/*--------------------------------------------------------------------------*/
@@ -3650,11 +3658,12 @@ void MCFSimplex::BalanceFlow( nodePType *r )
else {
// The method controls if "r" is a leaf in T
bool rLeaf = false;
- if( r->nextInT )
+ if( r->nextInT ) {
if( ( r->nextInT )->subTreeLevel <= r->subTreeLevel )
rLeaf = true;
else
rLeaf = true;
+ }
if( rLeaf ) // If "r" is a leaf
AdjustFlow( r ); // The method controls if entering basic arc in "r" is
@@ -3680,7 +3689,7 @@ void MCFSimplex::BalanceFlow( nodePType *r )
void MCFSimplex::AdjustFlow( nodePType *r )
{
arcPType *arc = r->enteringTArc;
- if( arc >= dummyArcsP ) // If entering arc of "r" is a dummy arc
+ if( arc >= dummyArcsP ) { // If entering arc of "r" is a dummy arc
if( LTZ( arc->flow , EpsFlw ) ) {
// If this dummy arc has flow < 0, the algorithm overturns the arc
nodePType *temp = arc->tail;
@@ -3749,6 +3758,7 @@ void MCFSimplex::AdjustFlow( nodePType *r )
}
}
}
+ }
/*--------------------------------------------------------------------------*/
@@ -3792,11 +3802,12 @@ void MCFSimplex::PostDVisit( nodeDType *r )
// The method controls if "r" is a leaf in T
bool rLeaf = false;
int i = r - nodesD;
- if( r->nextInT )
+ if( r->nextInT ) {
if( ( r->nextInT )->subTreeLevel <= r->subTreeLevel )
rLeaf = true;
else
rLeaf = true;
+ }
if( rLeaf ) // If "r" is a leaf
if( ( r->enteringTArc)->head == r ) // If enteringTArc of "r" goes in "r"
@@ -3819,12 +3830,13 @@ void MCFSimplex::PostDVisit( nodeDType *r )
desc = desc->nextInT;
}
- if( r != dummyRootD )
+ if( r != dummyRootD ) {
if( ( r->enteringTArc )->head == r ) // If enteringTArc of "r" goes in "r"
( r->enteringTArc )->flow = modifiedBalance[ i ];
else // If enteringTArc of "r" goes out "r"
( r->enteringTArc )->flow = - modifiedBalance[ i ];
}
+ }
#endif
}
@@ -3961,7 +3973,7 @@ void MCFSimplex::PrintDArc( arcDType *arc )
MCFSimplex::nodePType* MCFSimplex::RecoverPNode( Index ind )
{
- if( ( ind < 0 ) || ( ind > n ) )
+ if( ( ind > n ) )
return( NULL );
if( ind )
return( nodesP + ind - 1 );
@@ -3993,7 +4005,7 @@ MCFSimplex::arcPType* MCFSimplex::RecoverPArc( nodePType *tail ,
MCFSimplex::nodeDType* MCFSimplex::RecoverDNode( Index ind )
{
- if( ( ind < 0 ) || ( ind > n ) )
+ if( ( ind > n ) )
return( NULL );
if( ind )
diff --git a/clang/lib/Analysis/MCFSimplex.h b/clang/lib/Analysis/MCFSimplex.h
deleted file mode 100644
index a8d84f6..0000000
--- a/clang/lib/Analysis/MCFSimplex.h
+++ /dev/null
@@ -1,1086 +0,0 @@
-/*--------------------------------------------------------------------------*/
-/*------------------------- File MCFSimplex.h ------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @file
- * Linear and Quadratic Min Cost Flow problems solver based on the (primal and
- * dual) simplex algorithm. Conforms to the standard MCF interface defined in
- * MCFClass.h.
- *
- * \Version 1.00
- *
- * \date 29 - 08 - 2011
- *
- * \author Alessandro Bertolini \n
- * Antonio Frangioni \n
- * Operations Research Group \n
- * Dipartimento di Informatica \n
- * Universita' di Pisa \n
- *
- * Copyright &copy 2008 - 2011 by Alessandro Bertolini, Antonio Frangioni
- */
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*----------------------------- DEFINITIONS --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#ifndef __MCFSimplex
- #define __MCFSimplex /* self-identification: #endif at the end of the file */
-
-/*--------------------------------------------------------------------------*/
-/*------------------------------ INCLUDES ----------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#include "MCFClass.h"
-
-/*@} -----------------------------------------------------------------------*/
-/*--------------------------- NAMESPACE ------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
-namespace MCFClass_di_unipi_it
-{
-#endif
-
-/*--------------------------------------------------------------------------*/
-/*-------------------------------- MACROS ----------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup MCFSimplex_MACROS Compile-time switches in MCFSimplex.h
- There is only one macro in MCFSimplex, but it is very important!
- @{ */
-
-#define QUADRATICCOST 0
-
-/**< Setting QUADRATICCOST == 1 the solver can solve problems with linear and
- quadratic costs too (but the latter only with the Primal Simplex).
- The reason for having a macro is that when quadratic costs are present the
- "arcType" struct has the additional field "quadraticCost" to hold it.
- Furthermore, the field "ident" is not created because the solver doesn't
- use the classical TLU tripartition. Instead, closed arcs and deleted arcs
- are characterized as follows:
- - closed arcs have the field "cost" to INFINITY (Inf<FNumber>());
- - deleted arcs have the field "upper" to INFINITY and the "tail" and
- "head" field are NULL.
- Furthermore, the solver needs the variables "ignoredEnteringArc" and
- "firstIgnoredEnteringArc", used to avoid nasty loops during the execution
- of the Quadratic Primal Simplex algorithm.
- If, instead, QUADRATICCOST == 0 then the solver can solve only problems
- with linear costs. Hence, the field "quadraticCost" is useless and it
- isn't created. Furthermore, Primal Simplex and Dual Simplex use the
- tripartition TLU to divide the arcs, so the solver creates the field
- "ident", which differentiates the set of the arcs in among deleted arcs,
- closed arcs, arcs in T, arcs in L, arcs in U.
- Thus, with QUADRATICCOST == 0 the solver cannot solve problems with
- quadratic costs, but it does solve problems with linear costs faster. */
-
-/*@} end( group( MCFCLASS_MACROS ) ) */
-
-/*--------------------------------------------------------------------------*/
-/*---------------------------- CLASSES -------------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup MCFSimplex_CLASSES Classes in MCFSimplex.h
- @{ */
-
-/** The MCFSimplex class derives from the abstract base class MCFClass, thus
- sharing its (standard) interface, and implements both the Primal and
- Dual network simplex algorithms for solving (Linear and Quadratic)
- Min Cost Flow problems */
-
-class MCFSimplex: public MCFClass
-{
-
-/*--------------------------------------------------------------------------*/
-/*----------------------- PUBLIC PART OF THE CLASS -------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-- --*/
-/*-- The following methods and data are the actual interface of the --*/
-/*-- class: the standard user should use these methods and data only. --*/
-/*-- --*/
-/*--------------------------------------------------------------------------*/
-
- public:
-
-/*--------------------------------------------------------------------------*/
-/*---------------------------- PUBLIC TYPES --------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-/** Public enum describing the parameters of MCFSimplex. */
-
- enum SimplexParam
- {
- kAlgPrimal = kLastParam , ///< parameter to set algorithm (Primal/Dual):
- kAlgPricing , ///< parameter to set algorithm of pricing
- kNumCandList , /**< parameter to set the number of candidate
- list for Candidate List Pivot method */
- kHotListSize , /**< parameter to set the size of Hot List
- for Candidate List Pivot method */
- kRecomputeFOLimits , /**< parameter to set the number of iterations
- in which quadratic Primal Simplex computes
- "manually" the f.o. value */
- kEpsOpt /**< parameter to set the precision of the
- objective function value for the
- quadratic Primal Simplex */
- };
-
-/** Public enum describing the pricing rules in MCFSimplex::SetAlg(). */
-
- enum enumPrcngRl
- {
- kDantzig = 0, ///< Dantzig's rule (most violated constraint)
- kFirstEligibleArc , ///< First eligible arc in round-robin
- kCandidateListPivot ///< Candidate List Pivot Rule
- };
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------- PUBLIC METHODS -------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*---------------------------- CONSTRUCTOR ---------------------------------*/
-/*--------------------------------------------------------------------------*/
-
- MCFSimplex( cIndex nmx = 0 , cIndex mmx = 0 );
-
-/**< Constructor of the class, as in MCFClass::MCFClass(). */
-
-/*--------------------------------------------------------------------------*/
-/*-------------------------- OTHER INITIALIZATIONS -------------------------*/
-/*--------------------------------------------------------------------------*/
-
- void LoadNet( cIndex nmx = 0 , cIndex mmx = 0 , cIndex pn = 0 ,
- cIndex pm = 0 , cFRow pU = NULL , cCRow pC = NULL ,
- cFRow pDfct = NULL , cIndex_Set pSn = NULL ,
- cIndex_Set pEn = NULL );
-
-/**< Inputs a new network, as in MCFClass::LoadNet(). */
-
-/*--------------------------------------------------------------------------*/
-
- void SetAlg( bool UsPrml , char WhchPrc );
-
-/**< Selects which algorithm (Primal vs Dual Network Simplex), and which
- pricing rule within the algorithm, is used.
-
- If UsPrml == TRUE then the Primal Network Simplex algorithm is used,
- otherwise the Dual Network Simplex is used.
-
- The allowed values for WhchPrc are:
-
- - kDantzig Dantzig's pricing rule, i.e., most violated dual
- constraint; this can only be used with the Primal
- Network Simplex
-
- - kFirstEligibleArcA "dumb" rule, first eligible arc in round-robin;
-
- - kCandidateListPivot Candidate List Pivot Rule
-
- If this method is *not* called, the Primal Network Simplex with the
- Candidate List Pivot Rule (the best setting on most instances) is
- used. */
-
-/*--------------------------------------------------------------------------*/
-
- void SetPar( int par , int val );
-
-/**< Set general integer parameters.
-
- @param par is the parameter to set; since this method accepts an int
- value, the enum SimplexParam can be used in addition to the
- enum MCFParam to specify the integer parameters (every enum
- is an int).
-
- @param value is the value to assign to the parameter. */
-
-/*-------------------------------------------------------------------------*/
-
- void SetPar( int par , double val );
-
-/**< Set general float parameters.
-
- @param par is the parameter to set; since this method accepts an int
- value, the enum SimplexParam can be used in addition to the
- enum MCFParam to specify the float parameters (every enum
- is an int).
-
- @param value is the value to assign to the parameter. */
-
-/*--------------------------------------------------------------------------*/
-/*-------------------- METHODS FOR SOLVING THE PROBLEM ---------------------*/
-/*--------------------------------------------------------------------------*/
-
- void SolveMCF( void );
-
-/*--------------------------------------------------------------------------*/
-/*---------------------- METHODS FOR READING RESULTS -----------------------*/
-/*--------------------------------------------------------------------------*/
-
- void MCFGetX( FRow F , Index_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFGetRC( CRow CR , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() ) ;
-
- CNumber MCFGetRC( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFGetPi( CRow P , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
-/*--------------------------------------------------------------------------*/
-
- FONumber MCFGetFO( void );
-
-/*--------------------------------------------------------------------------*/
-/*-------------- METHODS FOR READING THE DATA OF THE PROBLEM ---------------*/
-/*--------------------------------------------------------------------------*/
-
- virtual void MCFArcs( Index_Set Startv , Index_Set Endv ,
- cIndex_Set nms = NULL , cIndex strt = 0 ,
- Index stp = Inf<Index>() );
-
- inline Index MCFSNde( cIndex i );
-
- inline Index MCFENde( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFCosts( CRow Costv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- inline CNumber MCFCost( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFQCoef( CRow Qv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- inline CNumber MCFQCoef( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFUCaps( FRow UCapv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- inline FNumber MCFUCap( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-
- void MCFDfcts( FRow Dfctv , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- inline FNumber MCFDfct( cIndex i );
-
-/*--------------------------------------------------------------------------*/
-/*------------- METHODS FOR ADDING / REMOVING / CHANGING DATA --------------*/
-/*--------------------------------------------------------------------------*/
-/*------- Changing the costs, QCoef, deficits and upper capacities ---------*/
-/*--------------------------------------------------------------------------*/
-
- void ChgCosts( cCRow NCost , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- void ChgCost( Index arc , cCNumber NCost );
-
-/*--------------------------------------------------------------------------*/
-
- void ChgQCoef( cCRow NQCoef = NULL , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- void ChgQCoef( Index arc , cCNumber NQCoef );
-
-/*--------------------------------------------------------------------------*/
-
- void ChgDfcts( cFRow NDfct , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- void ChgDfct( Index nod , cFNumber NDfct );
-
-/*--------------------------------------------------------------------------*/
-
- void ChgUCaps( cFRow NCap , cIndex_Set nms = NULL ,
- cIndex strt = 0 , Index stp = Inf<Index>() );
-
- void ChgUCap( Index arc , cFNumber NCap );
-
-/*--------------------------------------------------------------------------*/
-/*--------------- Modifying the structure of the graph ---------------------*/
-/*--------------------------------------------------------------------------*/
-
- void CloseArc( cIndex name );
-
- void DelNode( cIndex name );
-
- bool IsClosedArc( cIndex name );
-
- void OpenArc( cIndex name ) ;
-
- Index AddNode( cFNumber aDfct );
-
- void ChangeArc( cIndex name , cIndex nSS = Inf<Index>() ,
- cIndex nEN = Inf<Index>() );
-
- void DelArc( cIndex name );
-
- Index AddArc( cIndex Start , cIndex End , cFNumber aU , cCNumber aC );
-
- bool IsDeletedArc( cIndex name );
-
-/*--------------------------------------------------------------------------*/
-/*------------------------------ DESTRUCTOR --------------------------------*/
-/*--------------------------------------------------------------------------*/
-
- ~MCFSimplex();
-
-/*--------------------------------------------------------------------------*/
-/*--------------------- PRIVATE PART OF THE CLASS --------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-- --*/
-/*-- Nobody should ever look at this part: everything that is under this --*/
-/*-- advice may be changed without notice in any new release of the code. --*/
-/*-- --*/
-/*--------------------------------------------------------------------------*/
-
- private:
-
-/*--------------------------------------------------------------------------*/
-/*-------------------------- PRIVATE DATA TYPES ----------------------------*/
-/*--------------------------------------------------------------------------*/
-/*-- --*/
-/*-- Let T \subseteq A be a spanning tree, and consider some node v \in V --*/
-/*-- \setminus { 0 }. There is an unique (undirected) path, denoted by --*/
-/*-- P(v), defined by T from v to the root node 0. The arc in P(v), which --*/
-/*-- is incident to v, is called the *basic arc* of v. The other terminal --*/
-/*-- node u of this basic arc is called the *father node* of v. The --*/
-/*-- spanning tree T is represented saving the basic arc of every node, --*/
-/*-- and maintaining the order of the nodes and the depth as to T root --*/
-/*-- after a Post-Visit of T. This order is saved in a bidirectional list --*/
-/*-- written in the node. --*/
-/*-- --*/
-/*-- The Primal Simplex uses a different data structure than the Dual --*/
-/*-- Simplex, because the Dual Simplex needs additional data (mainly the --*/
-/*-- Backward Star and Forward Star. Furthermore, the Primal Simplex uses --*/
-/*-- different data structures in the quadratic case. --*/
-/*-- --*/
-/*--------------------------------------------------------------------------*/
-
- struct arcPType; // pre-declaration of the arc structure (pointers to arcs
- // are contained in the node structure) for Primal Simplex
-
- struct arcDType; // pre-declaration of the arc structure (pointers to arcs
- // are contained in the node structure) for Dual Simplex
-
- typedef double iteratorType; // type for the iteration counter and array
- // "whenInT2"
-
- struct nodePType { // node structure for Primal Simplex - - - - - - - -
- nodePType *prevInT; // previous node in the order of the Post-Visit on T
-
- nodePType *nextInT; // next node in the order of the Post-Visit on T
-
- arcPType *enteringTArc; // entering basic arc of this node
-
- FNumber balance; // supply/demand of this node; a node is called a
- // supply node, a demand node, or a transshipment
- // node depending upon whether balance is larger
- // than, smaller than, or equal to zero
- #if( QUADRATICCOST )
- CNumber sumQuadratic; // the sum of the quadratic coefficients of the tree's arcs
- // from root of T to the node
-
- FONumber potential; // the node potential corresponding with the flow
- // conservation constrait of this node
- #else
- CNumber potential; // the node potential corresponding with the flow
- // conservation constrait of this node
- #endif
-
- int subTreeLevel; // the depth of the node in T as to T root
-
- }; // end( struct( nodePType ) )
-
- struct nodeDType { // node structure for Dual Simplex - - - - - - - - -
- nodeDType *prevInT; // previous node in the order of the Post-Visit on T
-
- nodeDType *nextInT; // next node in the order of the Post-Visit on T
-
- arcDType *enteringTArc; // entering basic arc of this node
-
- FNumber balance; // supply/demand of this node; a node is called a
- // supply node, a demand node, or a transshipment
- // node depending upon whether balance is larger
- // than, smaller than, or equal to zero
-
- #if( QUADRATICCOST )
- CNumber sumQuadratic; // the sum of the quadratic coefficients of the tree's arcs
- // from root of T to the node
-
- FONumber potential; // the node potential corresponding with the flow
- // conservation constrait of this node
- #else
- CNumber potential; // the node potential corresponding with the flow
- // conservation constrait of this node
- #endif
-
- int subTreeLevel; // the depth of the node in T as to T root
- iteratorType whenInT2; // the last iteration where a node is in subtree T2
-
- Index numArcs; // the number of the arcs which enter/exit from node
- arcDType *firstBs; // the first arc in the node's Backward Star
- arcDType *firstFs; // the first arc in the node's Forward Star
-
- }; // end( struct( nodeDType ) )
-
- struct arcPType { // arc structure for Primal Simplex - - - - - - - -
- nodePType *tail; // tail node
- nodePType *head; // head node
- FNumber flow; // arc flow
- CNumber cost; // arc linear cost
-
- #if( QUADRATICCOST )
- CNumber quadraticCost; // arc quadratic cost
- #else
- char ident; // tells if arc is deleted, closed, in T, L, or U
- #endif
-
- FNumber upper; // arc upper bound
- }; // end( struct( arcPType ) )
-
- struct arcDType { // arc structure for Primal Simplex - - - - - - - -
- nodeDType *tail; // tail node
- nodeDType *head; // head node
- FNumber flow; // arc flow
- CNumber cost; // arc linear cost
-
- #if( QUADRATICCOST )
- CNumber quadraticCost; // arc quadratic cost
- #else
- char ident; // indicates if arc is deleted, closed, in T, in L, or in U
- #endif
-
- FNumber upper; // arc upper bound
- arcDType *nextBs; // the next arc in the Backward Star of the arc's head
- arcDType *nextFs; // the next arc in the Forward Star of the arc's tail
-
- }; // end( struct( arcDType ) )
-
- struct primalCandidType { // Primal Candidate List- - - - - - - - - - - - -
- arcPType *arc; // pointer to the violating primal bound arc
-
- #if(QUADRATICCOST )
- FONumber absRC; // absolute value of the arc's reduced cost
- #else
- CNumber absRC; // absolute value of the arc's reduced cost
- #endif
- }; // end( struct( primalCandidateType ) )
-
- struct dualCandidType { // Dual Candidate List- - - - - - - - - - - - - -
- nodeDType *node; // deepest node violating the dual bound arc
- FNumber absInfeas; // absolute value of the arc's flow infeasibility
- }; // end( struct( dualCandidateType ) )
-
-/*--------------------------------------------------------------------------*/
-/*----------------------- PRIVATE DATA STRUCTURES -------------------------*/
-/*--------------------------------------------------------------------------*/
-
- bool usePrimalSimplex; // TRUE if the Primal Network Simplex is used
-
- char pricingRule; // which pricing rule is used
-
- nodePType *nodesP; // vector of nodes: points to the n + 1 node
- // structs (including the dummy root node)
- // where the first node is indexed by zero
- // and the last node is the dummy root node
-
- nodePType *dummyRootP; // the dummy root node
-
- nodePType *stopNodesP; // first infeasible node address = nodes + n
-
- arcPType *arcsP; // vector of arcs; this variable points to
- // the m arc structs.
-
- arcPType *dummyArcsP; // vector of artificial dummy arcs: points
- // to the artificial dummy arc variables and
- // contains n arc structs. The artificial
- // arcs are used to build artificial feasible
- // starting bases. For each node i, there is
- // exactly one dummy arc defined to connect
- // the node i with the dummy root node.
-
- arcPType *stopArcsP; // first infeasible arc address = arcs + m
-
- arcPType *stopDummyP; // first infeasible dummy arc address
- // = arcs + m + n
-
- arcPType *arcToStartP; // Dantzig Rule and First Eligible Arc Rule
- // start their search from this arc
-
- nodeDType *nodesD; // vector of nodes: points to the n + 1 node
- // structs (including the dummy root node)
- // where the first node is indexed by zero
- // and the last node is the dummy root node
-
- nodeDType *dummyRootD; // the dummy root node
-
- nodeDType *stopNodesD; // first infeasible node address = nodes + n
-
- arcDType *arcsD; // vector of arcs; this variable points to
- // the m arc structs.
-
- arcDType *dummyArcsD; // vector of artificial dummy arcs: points to
- // to the artificial dummy arc variables and
- // contains n arc structs. The artificial
- // arcs are used to build artificial feasible
- // starting bases. For each node i, there is
- // exactly one dummy arc defined to connect
- // the node i with the dummy root node.
-
- arcDType *stopArcsD; // first infeasible arc address = arcs + m
-
- arcDType *stopDummyD; // first infeasible dummy arc address
- // = arcs + m + n
-
- arcDType *arcToStartD; // Dantzig Rule and First Eligible Arc Rule
- // start their search from this arc
-
- iteratorType iterator; // the current number of iterations
-
- primalCandidType *candP; // every element points to an element of the
- // arcs vector which contains an arc violating
- // dual bound
-
- dualCandidType *candD; // every element points to an element of the
- // arcs vector which contains an arc violating
- // primal bond
-
- Index numGroup; // number of the candidate lists
-
- Index tempCandidateListSize; // hot list dimension (it is variable)
-
- Index groupPos; // contains the actual candidate list
-
- Index numCandidateList; // number of candidate lists
-
- Index hotListSize; // number of candidate lists and hot list dimension
-
- Index forcedNumCandidateList; // value used to force the number of candidate list
-
- Index forcedHotListSize; // value used to force the number of candidate list
- // and hot list dimension
-
- bool newSession; // true if algorithm is just started
-
- CNumber MAX_ART_COST; // large cost for artificial arcs
-
- FNumber *modifiedBalance; // vector of balance used by the PostVisit
-
- FONumber EpsOpt; // the precision of the objective function value
- // for the quadratic case of the Primal Simplex
-
- int recomputeFOLimits; // after how many iterations the quadratic Primal
- // Simplex computes "manually" the o.f. value
-
- #if( QUADRATICCOST )
- FONumber foValue; // the temporary objective function value
- #endif
-
-/*--------------------------------------------------------------------------*/
-/*-------------------------- PRIVATE METHODS -------------------------------*/
-/*--------------------------------------------------------------------------*/
-
- void MemAlloc( void );
-
-/**< Method to allocate memory for the main data structures. It creates the
- vector of nodes, the vector of arcs (real and dummy). If the algorithm is
- using the Dual Simplex, it creates also the vectors whenInT2, firstIn and
- nextIn, usefull to identify the next entering. */
-
-/*--------------------------------------------------------------------------*/
-
- void MemDeAlloc( bool whatDeAlloc );
-
-/**< Method to deallocate memory for the main data structures created in
- MemAlloc(). */
-
-/*--------------------------------------------------------------------------*/
-
- void MemAllocCandidateList( void );
-
-/**< Method to allocate memory for the data structure used by the Candidate
- List Pivot Rule. It creates the vector of candP (or candD in the Dual
- Simplex case), determining the size of this vector on the basis of number
- of arcs or nodes. */
-
-/*--------------------------------------------------------------------------*/
-
- void MemDeAllocCandidateList( void );
-
-/**< Method to deallocate memory for the data structures created in
- MemAllocCandidateList(). */
-
-/*--------------------------------------------------------------------------*/
-
- void CreateInitialPrimalBase( void );
-
-/**< Method to create an initial feasible primal base. Add one node (dummyRoot)
- to the network and connect this dummy root to each real node with n dummy
- arcs. For each source node i, a dummy arc (i, r) exists. For each sink or
- transit node j, a dummy arc (r, j) exists. The source nodes send their flows
- to the dummy root through their dummy arcs, so the dummy root send all to
- the sink nodes. The dummy root balance is 0, the costs of dummy arcs are
- fixed to "infinity". */
-
-/*--------------------------------------------------------------------------*/
-
- void CreateInitialDualBase( void );
-
-/**< Method to create an initial feasible dual base. Add one node (dummyRoot)
- to the network and connect this dummy root to each real node with n dummy
- arcs. The source nodes send their flows to the dummy root through their
- dummy arcs, so the dummy root send all to the sink nodes. The dummy root
- balance is 0, the costs of dummy arcs are fixed to "infinity". */
-
-/*--------------------------------------------------------------------------*/
-
- void CreateAdditionalDualStructures( void );
-
-/**< The Dual Simplex needs nodes' Backward and Forward Star to work. So when
- the Primal Simplex runs, these structures don't exist. When the Dual
- Simplex starts, these structure are created in this method. */
-
-/*--------------------------------------------------------------------------*/
-
- void PrimalSimplex( void );
-
-/**< Main method to implement the Primal Simplex algorithm. */
-
-/*--------------------------------------------------------------------------*/
-
- void DualSimplex( void );
-
-/**< Main method to implement the Dual Simplex algorithm. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N, class A>
- void UpdateT( A *h , A *k , N *h1 , N *h2 , N *k1 , N *k2 );
-
- /**< Method to update the spanning tree T every iteration.
- The spanning tree T is implemented with a bidirectional list stored
- in the node structure, which represents the nodes' order after a
- Post-Visit. The parameter of the method are the outgoing arc "h", the
- incoming arc "k", and the four nodes on the extremity of these arcs
- (for example h2 is the deepest node of the outgoing arc). Removing the
- arc "h" splits T in two subtrees: T1 (which contains the root of T) and
- T2, which will be re-connected by the incoming arc "k".
- T2 will be reordered; in fact, the node "k2" becomes the root of T2
- instead of "h2" and the hierarchy of T2 will be overturned. Then T2 will
- be moved; the root of T2 is changed, therefore the predecessor of the
- root will become the node "k1".
-
- This method uses the methods cutAndUpdateSubtree() and pasteSubtree().
- First it cuts the node "k2" and its subtree from T using
- cutAndUpdateSubtree(). Moreover the method cutAndUpdateSubtree()
- updates the field "subTreeLevel" of every subtree's nodes, since k2's
- subtree will be moved from the bottom to the top of T2. Then the method
- pasteSubtree() puts this subtree in the bidirectional list after the node
- "k1". The same operations will be applied to the old precedessor of "k2"
- (which will become one of the childs of "k2"). This second subtree will
- be cut, the subTreeLevel fields will be updated, and it will be inserted
- in the bidirectional list after the k2's subtree. This is iterated until
- the node "h2" is reached. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N>
- N* CutAndUpdateSubtree( N *root, int delta );
-
-/**< This method cuts a generic subtree from the spanning tree T. Then it
- updates the field "subTreeLevel" of every subtree's nodes adding the value
- "delta". This method returns the last node of the subtree. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N>
- void PasteSubtree( N *root , N *lastNode , N *previousNode );
-
-/**< This method inserts a generic subtree with root passed by parameter into
- the spanning tree T, between the nodes "previousNode" and "lastNode". */
-
-/*--------------------------------------------------------------------------*/
-
- arcPType* RuleDantzig( void );
-
-/**< This method returns an arc which violates the dual conditions. It searchs
- the arc with most violation of dual conditions in the entire set of real
- arcs. It can be used only by the Primal Simplex in the case of networks
- with linear costs. */
-
-/*--------------------------------------------------------------------------*/
-
- arcPType* PRuleFirstEligibleArc( void );
-
-/**< This method returns the first found arc which violates the dual conditions
- in the case of Primal Simplex, the primal condition in the case of Dual
- Simplex. It can be used only in the case of networks with linear costs. */
-
-/*--------------------------------------------------------------------------*/
-
- arcDType* DRuleFirstEligibleArc( void );
-
-/**< This method returns the first found arc which violates the dual conditions
- in the case of Primal Simplex, the primal condition in the case of Dual
- Simplex. It can be used only in the case of networks with linear costs. */
-
-/*--------------------------------------------------------------------------*/
-
- arcPType* RulePrimalCandidateListPivot( void );
-
-/**< This method returns an arc which violates the dual conditions. It searches
- the arc with most violation of dual conditions in a small set of candidate
- arcs. In every iteration the method rebuilds this set of arcs executing three
- phases:
-
- - in the first phase it analyzes the remaining arcs and delete the arcs which
- don't violate the dual condition any more;
-
- - in the second phase it tries to fill the set, so it searchs other arcs
- which violate the dual condition: the set of arcs is divided into "buckets"
- which are searched sequentially until the candidate list is full; the
- last visited bucket is retained, and the search is restarted from that
- one at later iterations
-
- - in the third phase the small set of candidate arcs is ordered according
- to the violation of dual condition by the method SortPrimalCandidateList()
- using an implementation of the algorithm "quicksort".
-
- At last the method returns the first arc in the ordered small set. If the
- arc doesn't exist (the set is empty), it returns NULL. */
-
-/*--------------------------------------------------------------------------*/
-
- inline void InitializePrimalCandidateList( void );
-
-/**< Method to initialize some important details for Primal Candidate List
- Rule. */
-
-/*--------------------------------------------------------------------------*/
-
- inline void SortPrimalCandidateList( Index min , Index max );
-
-/**< Method to order the little set of candidate arcs according to
- infeasibility of dual conditions. It implements the "quicksort"
- algorithm. */
-
-/*--------------------------------------------------------------------------*/
-
- arcDType* RuleDualCandidateListPivot( void );
-
-/**< Similar to RulePrimalCandidateListPivot() for the Dual Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- inline void InitializeDualCandidateList( void );
-
-/**< Method to initialize some important details for Dual Candidate List Rule.
- */
-
-/*--------------------------------------------------------------------------*/
-
- inline void SortDualCandidateList( Index min , Index max );
-
-/**< Similar to SortPrimalCandidateList() for the Dual Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N, class RCT>
- inline void AddPotential( N *r , RCT delta );
-
-/**< Method to quickly update the dual solutions. During the change of the
- base, the potential of node "k2" (deepest node in T of incoming arc "k",
- and new root of T2) changes according to the new structure of T. In fact,
- the precedessor of "k2" changes: now the predecessor of "k2"is "k1" (the
- other node of the incoming arc "k"). This change of predecessor causes a
- change of potential of "k2". The change of potential of "k2" causes the
- changes of potential of all nodes of T2. This method computes the change
- of potential of "k2" and applies it on all the nodes of T2. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N>
- inline void ComputePotential( N *r );
-
-/**< Method to update the dual solutions. It computes all the potential of the
- nodes of the subtree which has r as root. */
-
-/*--------------------------------------------------------------------------*/
-
- inline void ResetWhenInT2( void );
-
-/**< Method to order the small set of candidate arcs according to dual
- infeasibility. It implements the algorithm "quicksort". */
-
-/*--------------------------------------------------------------------------*/
-
- void CreateInitialPModifiedBalanceVector( void );
-
-/**< Method to initialize the vector of modified balance for the postvisit on
- T in the Primal Simplex data structure. */
-
-/*--------------------------------------------------------------------------*/
-
- void PostPVisit( nodePType *r );
-
-/**< Method to calculate the flow on the basic arcs with the Primal Simplex's
- data structure. It uses the set of the upper bound arcs, the construction
- of a modified balance vector and the postvisit on T. */
-
-/*--------------------------------------------------------------------------*/
-
- void BalanceFlow( nodePType *r );
-
-/**< This method works after the method PostPVisit (in a Primal context). It
- restores primal admissimibility on the r's subtree. It starts from the leaf
- of the subtree and goes up to the root, using the method AdjustFlow. */
-
-/*--------------------------------------------------------------------------*/
-
- void AdjustFlow( nodePType *r );
-
-/**< This method checks the primal admissibility of the r's basic entering arc.
- If it is out of bounds, the method removes it from the tree (and keeps the
- relative dummy arc) and push flow in the cycle (some tree's arc and the old
- entering arc) to restores the right balances of the node. */
-
-/*--------------------------------------------------------------------------*/
-
- void CreateInitialDModifiedBalanceVector( void );
-
-/**< Method to initialize the vector of modified balance for the postvisit on
- T in the Dual Simplex data structure. */
-
-/*--------------------------------------------------------------------------*/
-
- void PostDVisit( nodeDType *r );
-
-/**< Method to calculate the flow on the basic arcs with the Dual Simplex's data
- structure, using the set of the upper bound arcs, the construction of a
- modified balance vector and the postvisit on T. */
-
-/*--------------------------------------------------------------------------*/
-
- template<class N, class A>
- inline N* Father( N *n, A *a );
-
-/**< Method to find the predecessor of the node in the tree. */
-
-/*--------------------------------------------------------------------------*/
-
- #if(QUADRATICCOST)
- template<class A>
- inline FONumber ReductCost( A *a );
- #else
- template<class A>
- inline CNumber ReductCost( A *a );
- #endif
-
-/**< Method to calculate the reduct cost of the arc. */
-
-/*--------------------------------------------------------------------------*/
-
- inline FONumber GetFO( void );
-
-/**< Method to calculate the temporary (or the final) objective function
- value. */
-
-/*--------------------------------------------------------------------------*/
-
- void PrintPNode( nodePType *nodo );
-
-/**< Method to print the "name" of the node in the Primal Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- void PrintPArc( arcPType *arc );
-
-/**< Method to print the "name" of the arc in the Primal Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- void PrintDNode( nodeDType *nodo );
-
-/**< Method to print the "name" of the node in the Dual Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- void PrintDArc( arcDType *arc );
-
-/**< Method to print the "name" of the arc in the Dual Simplex. */
-
-/*--------------------------------------------------------------------------*/
-
- nodePType* RecoverPNode( Index ind );
-
-/**< Method to find a node (in the Primal Simplex) using its index. */
-
-/*--------------------------------------------------------------------------*/
-
- arcPType* RecoverPArc( nodePType *tail , nodePType *head );
-
-/**< Method to find an arc (in the Primal Simplex) using 2 pointers to tail
- node and head node. */
-
-/*--------------------------------------------------------------------------*/
-
- nodeDType* RecoverDNode( Index ind );
-
-/**< Method to find a node (in the Dual Simplex) using its index. */
-
-/*--------------------------------------------------------------------------*/
-
- arcDType* RecoverDArc( nodeDType *tail , nodeDType *head );
-
-/**< Method to find an arc (in the Dual Simplex) using 2 pointers to tail
- node and head node. */
-
-/*--------------------------------------------------------------------------*/
-
- void infoPNode( nodePType *node , int tab );
-
-/**< Method to print some information of the node (in the Primal Simplex). */
-
-/*--------------------------------------------------------------------------*/
-
- void infoPArc( arcPType *arc , int ind , int tab );
-
-/**< Method to print some information of the arc (in the Primal Simplex). */
-
-/*--------------------------------------------------------------------------*/
-
- void infoDNode( nodeDType *node , int tab );
-
-/**< Method to print some information of the node (in the Dual Simplex). */
-
-/*--------------------------------------------------------------------------*/
-
- void infoDArc( arcDType *arc , int ind , int tab );
-
-/**< Method to print some information of the arc (in the Dual Simplex). */
-
-/*--------------------------------------------------------------------------*/
-
- void ShowSituation( int tab );
-
-/**< Method to show the actual complete situation. */
-
-/*--------------------------------------------------------------------------*/
-
- }; // end( class MCFSimplex )
-
-/* @} end( group( MCFSimplex_CLASSES ) ) */
-
-#endif /* MCFSimplex.h included */
-
-/*-------------------------------------------------------------------------*/
-/*-------------------inline methods implementation-------------------------*/
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::Index MCFSimplex::MCFSNde( MCFClass::cIndex i )
-{
- if( usePrimalSimplex )
- return( Index( ( (arcsP + i)->tail - nodesP + 1 ) - USENAME0 ) );
- else
- return( Index( ( (arcsD + i)->tail - nodesD + 1 ) - USENAME0 ) );
- }
-
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::Index MCFSimplex::MCFENde( MCFClass::cIndex i )
-{
- if( usePrimalSimplex )
- return( Index( ( (arcsP + i)->head - nodesP + 1 ) - USENAME0 ) );
- else
- return( Index( ( (arcsD + i)->head - nodesD + 1 ) - USENAME0 ) );
- }
-
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::CNumber MCFSimplex::MCFCost( MCFClass::cIndex i )
-{
- if( usePrimalSimplex )
- return( (arcsP + i)->cost );
- else
- return( (arcsD + i)->cost );
- }
-
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::CNumber MCFSimplex::MCFQCoef( MCFClass::cIndex i )
-{
- #if( QUADRATICCOST )
- if( usePrimalSimplex )
- return( (arcsP + i)->quadraticCost );
- else
- return( (arcsD + i)->quadraticCost );
- #else
- return( 0 );
- #endif
- }
-
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::FNumber MCFSimplex::MCFUCap( MCFClass::cIndex i )
-{
- if( usePrimalSimplex )
- return( (arcsP + i)->upper );
- else
- return( (arcsD + i)->upper );
- }
-
-/*-------------------------------------------------------------------------*/
-
-inline MCFClass::FNumber MCFSimplex::MCFDfct( MCFClass::cIndex i )
-{
- if( usePrimalSimplex )
- return( (nodesP + i)->balance );
- else
- return( (nodesD + i)->balance );
- }
-
-/*-------------------------------------------------------------------------*/
-
-#if( QUADRATICCOST )
-
-template <class A>
-inline MCFSimplex::FONumber MCFSimplex::ReductCost( A *a )
-{
- FONumber redc = (a->tail)->potential - (a->head)->potential;
- redc = redc + a->cost;
- redc = redc + a->quadraticCost * a->flow;
- return( redc );
- }
-
-#else
-
-template <class A>
-inline MCFSimplex::CNumber MCFSimplex::ReductCost( A *a )
-{
- CNumber redc = (a->tail)->potential - (a->head)->potential;
- redc = redc + a->cost;
- return( redc );
- }
-
-#endif
-
-/*-------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
-}; // end( namespace MCFClass_di_unipi_it )
-#endif
-
-/*-------------------------------------------------------------------------*/
-/*-------------------------------------------------------------------------*/
-
-/*-------------------------------------------------------------------------*/
-/*---------------------- End File MCFSimplex.h ----------------------------*/
-/*-------------------------------------------------------------------------*/
diff --git a/clang/lib/Analysis/OPTUtils.h b/clang/lib/Analysis/OPTUtils.h
deleted file mode 100755
index 66ad4af..0000000
--- a/clang/lib/Analysis/OPTUtils.h
+++ /dev/null
@@ -1,475 +0,0 @@
-/*--------------------------------------------------------------------------*/
-/*--------------------------- File OPTutils.h ------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @file
- * Small classes are provided for:
- * - reading the time of a code;
- * - generating random numbers.
- *
- * The classes can be adapted to different environments setting a
- * compile-time switch in this file.
- *
- * Additionally, a function is provided for safely reading parameters
- * out of a stream.
- *
- * \version 1.00
- *
- * \date 04 - 10 - 2010
- *
- * \author Antonio Frangioni \n
- * Operations Research Group \n
- * Dipartimento di Informatica \n
- * Universita' di Pisa \n
- *
- * Copyright &copy 1994 - 2010 by Antonio Frangioni
- */
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#ifndef __OPTutils
- #define __OPTutils /* self-identification: #endif at the end of the file */
-
-/*--------------------------------------------------------------------------*/
-/*----------------------------- MACROS -------------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup OPTUTILS_MACROS Compile-time switches in OPTutils.h
- These macros control how the classes OPTTimers and OPTrand are
- implemented; choose the appropriate value for your environment,
- or program a new version if no value suits you.
- Also, namespaces can be eliminated if they create problems.
- @{ */
-
-
-/*----------------------- OPT_USE_NAMESPACES -------------------------------*/
-
-#define OPT_USE_NAMESPACES 0
-
-/**< Setting OPT_USE_NAMESPACES == 0 should instruct all codes that use
- OPTutils stuff to avoid using namespaces; to start with, the common
- namespace OPTutils_di_unipi_it, that contains all the types defined
- herein, is *not* defined. */
-
-/*---------------------------- OPT_TIMERS ----------------------------------*/
-
-#define OPT_TIMERS 5
-
-/**< The class OPTtimers is defined below to give an abstract interface to the
- different timing routines that are used in different platforms. This is
- needed since time-related functions are one of the less standard parts of
- the C[++] library. The value of the OPT_TIMERS constant selects among the
- different timing routines:
-
- - 1 = Use the Unix times() routine in sys/times.h
-
- - 2 = As 1 but uses sys/timeb.h (typical for Microsoft(TM) compilers)
-
- - 3 = Still use times() of sys/times.h, but return wallclock time
- rather than CPU time
-
- - 4 = As 3 but uses sys/timeb.h (typical for Microsoft(TM) compilers)
-
- - 5 = return the user time obtained with ANSI C clock() function; this
- may result in more accurate running times w.r.t. but may be limited
- to ~ 72 hours on systems where ints are 32bits.
-
- - 6 = Use the Unix gettimeofday() routine of sys/time.h.
-
- Any unsupported value would simply make the class to report constant
- zero as the time.
-
- The values 1 .. 4 rely on the constant CLK_TCK for converting between
- clock ticks and seconds; for the case where the constant is not defined by
- the compiler -- should not happen, but it does -- or it is defined in a
- wrong way, the constant is re-defined below. */
-
-/*---------------------------- OPT_RANDOM ---------------------------------*/
-
-#define OPT_RANDOM 1
-
-/**< The class OPTrand is defined below to give an abstract interface to the
- different random generators that are used in different platforms. This is
- needed since random generators are one of the less standard parts of the
- C[++] library. The value of the OPT_RANDOM constant selects among the
- different timing routines:
-
- - 0 = an hand-made implementation of a rather good random number generator
- is used; note that this assumes that long ints >= 32 bits
-
- - 1 = standard rand() / srand() pair, common to all C libreries but not
- very sophisticated
-
- - 2 = drand48() / srand48(), common on Unix architectures and pretty good.
-
- Any unsupported value would simply make the functions to report constant
- zero, which is not nice but useful to quickly fix problems if you don't
- use random numbers at all. */
-
-/*@} -----------------------------------------------------------------------*/
-/*------------------------------ INCLUDES ----------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_RANDOM )
- #include <stdlib.h>
-
- /* For random routines, see OPTrand() below. */
-#endif
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-
-#if( OPT_TIMERS <= 4 )
- #if( ( OPT_TIMERS == 1 ) || ( OPT_TIMERS == 3 ) )
- #include <sys/times.h>
- #else
- #include <sys/timeb.h>
- #endif
- #include <limits.h>
-#elif( OPT_TIMERS == 5 )
- #include <time.h>
-#elif( OPT_TIMERS == 6 )
- #include <sys/time.h>
-#endif
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-
-#include <iostream>
-#include <sstream>
-
-/* For istream and the >> operator, used in DfltdSfInpt(). */
-
-/*- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -*/
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------- NAMESPACE ------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
-namespace OPTtypes_di_unipi_it
-{
- /** @namespace OPTtypes_di_unipi_it
- The namespace OPTtypes_di_unipi_it is defined to hold all the data
- types, constants, classes and functions defined here. It also
- comprises the namespace std. */
-#endif
-
- using namespace std; // I know it's not elegant, but ...
-
-/*--------------------------------------------------------------------------*/
-/*----------------------------- NULL ---------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#ifndef NULL
- #define NULL 0
-#endif
-
-/* Curiously enough, some compilers do not always define NULL; for instance,
- sometimes it is defined in stdlib.h. */
-
-/*--------------------------------------------------------------------------*/
-/*---------------------------- CLK_TCK -------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#ifndef CLK_TCK
- #define CLK_TCK 100
-#endif
-
-/* CLK_TCK is the constant used to transform (integer) time values in seconds
- for times()-based routines. Its normal value is 100. */
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------- OPT_TIMERS -----------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup OPTtypes_CLASSES Classes in OPTutils.h
- @{ */
-
-#if( OPT_TIMERS )
-
-/** Provides a common interface to the different timing routines that are
- available in different platforms. */
-
-class OPTtimers {
-
- public: //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- /// constructor of the class
- OPTtimers( void )
- {
- ReSet();
- }
-
- /// start the timer
- void Start( void )
- {
- if( ! ticking ) {
- #if( ( OPT_TIMERS == 1 ) || ( OPT_TIMERS == 2 ) )
- times( &buff );
- t_u = buff.tms_utime;
- t_s = buff.tms_stime;
- #elif( ( OPT_TIMERS == 3 ) || ( OPT_TIMERS == 4 ) )
- t_u = times( &buff );
- #elif( OPT_TIMERS == 5 )
- t_u = clock();
- #elif( OPT_TIMERS == 6 )
- struct timeval t;
- gettimeofday( &t , NULL );
- t_u = double( t.tv_sec + t.tv_usec * 1e-6 );
- #endif
-
- ticking = 1;
- }
- }
-
- /// stop the timer
- void Stop( void )
- {
- if( ticking ) {
- Read( u , s );
- ticking = 0;
- }
- }
-
- /** Return the elapsed time. If the clock is ticking, return the *total*
- time since the last Start() without stopping the clock; otherwise,
- return the total elapsed time of all the past runs of the clock since
- the last ReSet() [see below]. */
-
- double Read( void )
- {
- double tu = 0;
- double ts = 0;
- Read( tu , ts );
- return( tu + ts );
- }
-
- /// As Read( void ) but *adds* user and system time to tu and ts.
- void Read( double &tu , double &ts )
- {
- if( ticking ) {
- #if( ( OPT_TIMERS == 1 ) || ( OPT_TIMERS == 2 ) )
- times( &buff );
- tu += ( double( buff.tms_utime - t_u ) ) / double( CLK_TCK );
- ts += ( double( buff.tms_stime - t_s ) ) / double( CLK_TCK );
- #elif( ( OPT_TIMERS == 3 ) || ( OPT_TIMERS == 4 ) )
- tu += ( double( times( &buff ) - t_u ) ) / double( CLK_TCK );
- #elif( OPT_TIMERS == 5 )
- tu += ( double( clock() - t_u ) ) / double( CLOCKS_PER_SEC );
- #elif( OPT_TIMERS == 6 )
- struct timeval t;
- gettimeofday( &t , NULL );
- tu += double( t.tv_sec + t.tv_usec * 1e-6 ) - t_u;
- #endif
- }
- else { tu += u; ts += s; }
- }
-
- /// reset the timer
- void ReSet( void )
- {
- u = s = 0; ticking = 0;
- }
-
- private: //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- double u; // elapsed *user* time, in seconds
- double s; // elapsed *system* time, in seconds
- char ticking; // 1 if the clock is ticking
-
- #if( ( OPT_TIMERS > 0 ) && ( OPT_TIMERS <= 5 ) )
- clock_t t_u;
-
- #if( OPT_TIMERS <= 4 )
- struct tms buff;
-
- #if( ( OPT_TIMERS == 1 ) || ( OPT_TIMERS == 3 ) )
- clock_t t_s;
- #endif
- #endif
- #elif( OPT_TIMERS == 6 )
- double t_u;
- #endif
-
- }; // end( class OPTtimers );
-
-#endif
-
-/*--------------------------------------------------------------------------*/
-/*------------------------------ OPTrand() ---------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-/** Provide a common interface to the different random generators that are
- available in different platforms. */
-
-class OPTrand {
-
- public: //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- /// constructor of the class
- OPTrand( void )
- {
- #if( OPT_RANDOM == 0 )
- A[ 0 ] = -1;
- srand( long( 123456789 ) );
- #else
- OPTrand::srand( long( 1 ) );
- #endif
- }
-
- /** Returns a random number uniformly distributed in [0, 1).
- \note each object of class OPTrand has its own sequence, so that
- multiple OPTrand objects being used within the same program do not
- interfere with each other (as opposed to what C random routines
- would do). */
-
- double rand( void )
- {
- #if( OPT_RANDOM == 0 )
- long nmbr = *(gb_fptr--);
- if( nmbr < 0 )
- nmbr = gb_flip_cycle();
- return( double( nmbr ) / double( (unsigned long)0x80000000 ) );
- #elif( OPT_RANDOM == 1 )
- ::srand( myseed );
- myseed = ::rand();
- return( double( myseed ) / double( RAND_MAX ) );
- #elif( OPT_RANDOM == 2 )
- return( erand48( myseed ) );
- #else
- return( 0 ); // random, eh?
- #endif
- }
-
- /// Seeds the random generator for this instance of OPTrand.
- void srand( long seed )
- {
- #if( OPT_RANDOM == 0 )
- long prev = seed , next = 1;
- seed = prev = mod_diff( prev , 0 );
- A[ 55 ] = prev;
-
- for( long i = 21 ; i ; i = ( i + 21 ) % 55 ) {
- A[ i ] = next;
- next = mod_diff( prev , next );
- if( seed & 1 )
- seed = 0x40000000 + ( seed >> 1 );
- else
- seed >>= 1;
-
- next = mod_diff( next , seed );
- prev = A[ i ];
- }
-
- gb_flip_cycle();
- gb_flip_cycle();
- gb_flip_cycle();
- gb_flip_cycle();
- gb_flip_cycle();
- #elif( OPT_RANDOM == 1 )
- myseed = int( seed );
- #elif( OPT_RANDOM == 2 )
- long *sp = (long*)( &myseed );
- *sp = seed; // copy higher 32 bits
- myseed[ 2 ] = 0x330E; // initialize lower 16 bits
- #endif
- }
-
- private: //- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
-
- #if( OPT_RANDOM == 0 )
- long A[ 56 ];
- long *gb_fptr;
-
- long mod_diff( long x , long y )
- {
- return( ( ( x ) - ( y ) ) & 0x7fffffff );
- }
-
- long gb_flip_cycle( void )
- {
- long *ii, *jj;
- for( ii = &A[ 1 ] , jj = &A[ 32 ] ; jj <= &A[ 55 ] ; ii++ , jj++ )
- *ii = mod_diff( *ii , *jj );
-
- for( jj = &A[ 1 ] ; ii <= &A[ 55 ] ; ii++ , jj++ )
- *ii = mod_diff( *ii , *jj );
-
- gb_fptr = &A[ 54 ];
-
- return A[ 55 ];
- }
- #elif( OPT_RANDOM == 1 )
- int myseed;
- #elif( OPT_RANDOM == 2 )
- unsigned short int myseed[ 3 ];
- #endif
-
- }; // end( class( OPTrand ) )
-
-/* @} end( group( OPTtypes_CLASSES ) ) */
-/*--------------------------------------------------------------------------*/
-/*--------------------------- DfltdSfInpt() --------------------------------*/
-/*--------------------------------------------------------------------------*/
-/** @defgroup OPTtypes_FUNCTIONS Functions in OPTutils.h
- @{ */
-
-/** Template function for reading parameters from a istream. The function is
- "safe" because it works also if the istream is not given, is not be long
- enough or contains erroneous things.
-
- Given a &istream (possibly NULL), DfltdSfInpt() attempts to read Param out
- of it, skipping any line that begins with the comment carachter (defaulted
- to '#'), any blank line and any line starting with anything that can not
- be interpreted as a `T'. If, for any reason, the read operation fails,
- then the parameter is given the default value `Dflt'. Otherwise, all the
- rest of the line up to the nearest newline ('\n') carachter is flushed. */
-
-template<class T>
-inline void DfltdSfInpt( istream *iStrm , T &Param , const T Dflt ,
- const char cmntc = '#' )
-{
- string comm( 1 , cmntc );
-
- // the "! ! stream" trick is there to force the compiler to apply the
- // stream -> bool conversion, in case it is too dumb to do it by itself
- if( iStrm && ( ! ( ! (*iStrm) ) ) ) {
- string buf;
-
- for( ;; ) {
- if( ! ( (*iStrm) >> ws ) ) // skip whitespace
- break;
-
- if( ! (*iStrm).good() ) // check stream is OK
- break;
-
- getline( *iStrm , buf );
-
- if( ! buf.empty() )
- if( buf.substr( 0 , 1 ).compare( comm ) ) {
- stringstream ss;
- ss << buf;
- if( ! ( ss >> Param ) )
- break;
-
- return;
- }
- }
- }
-
- Param = Dflt;
-
- } // end( DfltdSfInpt )
-
-/* @} end( group( OPTtypes_FUNCTIONS ) ) */
-/*--------------------------------------------------------------------------*/
-
-#if( OPT_USE_NAMESPACES )
- }; // end( namespace OPTtypes_di_unipi_it )
-#endif
-
-/*--------------------------------------------------------------------------*/
-/*--------------------------------------------------------------------------*/
-
-#endif /* OPTutils.h included */
-
-/*--------------------------------------------------------------------------*/
-/*---------------------- End File OPTutils.h -------------------------------*/
-/*--------------------------------------------------------------------------*/
diff --git a/impl/Operator.hpp b/impl/Operator.hpp
index 4a573b8..0b08866 100644
--- a/impl/Operator.hpp
+++ b/impl/Operator.hpp
@@ -153,8 +153,8 @@ struct MinCostFlow : public Operator<Domain> {
ends[i] = arcs[i].second;
}
- _solver.LoadNet(supplies.size(), supplies.size(), // max nodes/arcs
- supplies.size(), supplies.size(), // current nodes/arcs
+ _solver.LoadNet(supplies.size(), arcs.size(), // max nodes/arcs
+ supplies.size(), arcs.size(), // current nodes/arcs
NULL, NULL, // arcs have inf cap, zero cost (to begin)
deficits, // deficits for each node
starts, ends); // start/end of each edge