diff --git a/inc/solver.hpp b/inc/solver.hpp index e3383b2..0fdd6cb 100644 --- a/inc/solver.hpp +++ b/inc/solver.hpp @@ -5,6 +5,9 @@ using std::vector; class Expr; class Var; +#define MDL_MINIMIZE 1 +#define MDL_MAXIMIZE -1 + Expr operator+(const Expr& x, const Expr& y); Expr operator-(const Expr& x, const Expr& y); Expr operator+(const Expr& x); @@ -25,7 +28,6 @@ Expr operator/(const Expr& x, double a); typedef vector> Matrix; - enum class Rtn { LOADED, OPTIMAL, @@ -47,10 +49,212 @@ enum class Rtn { PIVOT, }; +enum class IntAttr { + NumConstrs, + NumVars, + NumSOS, + NumQConstrs, + NumGenConstrs, + NumNZs, + NumQNZs, + NumQCNZs, + NumIntVars, + NumBinVars, + NumPWLObjVars, + ModelSense, + IsMIP, + IsQP, + IsQCP, + IsMultiObj, + Status, + ConcurrentWinMethod, + SolCount, + BarIterCount, + VBasis, + CBasis, + PWLObjCvx, + BranchPriority, + VarPreStat, + BoundVioIndex, + BoundSVioIndex, + ConstrVioIndex, + ConstrSVioIndex, + ConstrResidualIndex, + ConstrSResidualIndex, + DualVioIndex, + DualSVioIndex, + DualResidualIndex, + DualSResidualIndex, + ComplVioIndex, + IntVioIndex, + IISMinimal, + IISLB, + IISUB, + IISConstr, + IISSOS, + IISQConstr, + IISGenConstr, + IISLBForce, + IISUBForce, + IISConstrForce, + IISSOSForce, + IISQConstrForce, + IISGenConstrForce, + TuneResultCount, + Lazy, + VarHintPri, + ObjNPriority, + NumObj, + GenConstrType, + NumStart, + Partition, + LicenseExpiration, + NumScenarios, + FuncPieces, + BatchErrorCode, + BatchStatus, + Fingerprint, + PoolIgnore, + FuncNonlinear, + ModelColumn +}; + +enum class DoubleAttr{ + Runtime, + Work, + ObjCon, + LB, + UB, + Obj, + Start, + PreFixVal, + RHS, + QCRHS, + Coeff, + MaxCoeff, + MinCoeff, + MaxBound, + MinBound, + MaxObjCoeff, + MinObjCoeff, + MaxRHS, + MinRHS, + MaxQCRHS, + MinQCRHS, + MaxQCCoeff, + MinQCCoeff, + MaxQCLCoeff, + MinQCLCoeff, + MaxQObjCoeff, + MinQObjCoeff, + ObjVal, + ObjBound, + ObjBoundC, + MIPGap, + IterCount, + NodeCount, + X, + RC, + Pi, + QCPi, + Slack, + QCSlack, + MaxVio, + BoundVio, + BoundSVio, + BoundVioSum, + BoundSVioSum, + ConstrVio, + ConstrSVio, + ConstrVioSum, + ConstrSVioSum, + ConstrResidual, + ConstrSResidual, + ConstrResidualSum, + ConstrSResidualSum, + DualVio, + DualSVio, + DualVioSum, + DualSVioSum, + DualResidual, + DualSResidual, + DualResidualSum, + DualSResidualSum, + ComplVio, + ComplVioSum, + IntVio, + IntVioSum, + Kappa, + KappaExact, + SAObjLow, + SAObjUp, + SALBLow, + SALBUp, + SARHSLow, + SAUBLow, + SAUBUp, + SARHSUp, + Xn, + FarkasProof, + FarkasDual, + UnbdRay, + PStart, + DStart, + BarX, + VarHintVal, + ObjN, + ObjNCon, + ObjNWeight, + ObjNRelTol, + ObjNAbsTol, + ObjNVal, + OpenNodeCount, + PoolObjBound, + PoolObjVal, + ScenNLB, + ScenNUB, + ScenNObj, + ScenNRHS, + ScenNX, + ScenNObjBound, + ScenNObjVal, + FuncPieceError, + FuncPieceLength, + FuncPieceRatio, + DNumNZs +}; + +enum class StringAttr { + ModelName, + VarName, + ConstrName, + QCName, + GenConstrName, + ObjNName, + ScenNName, + BatchID, + VTag, + CTag, + QCTag, + BatchErrorMessage +}; + +enum class ConstrOper { + LESS_EQUAL, + GREATER_EQUAL, + EQUAL +}; + class Var { public: - Var(double coef = 1) :col(0), coeffs(coef) {}; + Var(double coef = 1) :col(0), val(0), coeffs(coef) {}; + double get(DoubleAttr attr); + int get(IntAttr attr); + friend class Model; + friend class Expr; +private: double coeffs; + double val; int col; }; @@ -101,18 +305,17 @@ class Model public: Model(); Var* addVars(int col); - void addConstr(const Expr& expr, char sense, double rhs); - void setObjective(Expr obje, int sense = 0); + void addConstr(const Expr& expr, ConstrOper sense, double rhs); + void setObjective(Expr obje, int sense = MDL_MAXIMIZE); void print(); - Rtn optimize(bool isMax = true); - double get(); - + Rtn optimize(); + double get(DoubleAttr attr); + int get(IntAttr attr); private: double _simplex(); Rtn _pivot(std::pair& p); void _gaussian(std::pair p); - Var* vars; Matrix matrix; @@ -120,4 +323,5 @@ private: vector basic; Rtn rtn; double res; + int sense; }; \ No newline at end of file diff --git a/src/main.cpp b/src/main.cpp index 1e5889c..5430942 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -9,24 +9,26 @@ int main(int argc, char* argv[]) { Model mdl; - Var* vars = mdl.addVars(3); + Var* vars = mdl.addVars(2); - mdl.addConstr(vars[0] + vars[1] + vars[2], '<', 4); - mdl.addConstr(vars[0], '<', 2); - mdl.addConstr(vars[2], '<', 3); - mdl.addConstr(3 * vars[1] + vars[2], '<', 6); + mdl.addConstr(vars[0] + vars[1], ConstrOper::LESS_EQUAL, 2); + mdl.addConstr(vars[0] + vars[1], ConstrOper::GREATER_EQUAL, 1); + mdl.addConstr(5 * vars[0] - 2 * vars[1], ConstrOper::GREATER_EQUAL, -2); - mdl.setObjective(vars[0] + 14 * vars[1] + 6 * vars[2]); + mdl.setObjective(vars[0] + 2 * vars[1], MDL_MAXIMIZE); - switch(mdl.optimize(false)) { + switch(mdl.optimize()) { case Rtn::OPTIMAL: - cout << "find solution: " << mdl.get() << endl; + cout << "OPTIMAL SOLUTION: " << mdl.get(DoubleAttr::Obj) << endl; + for (int i = 0; i < mdl.get(IntAttr::NumVars); i++) { + cout << "var[" << i << "] : " << vars[i].get(DoubleAttr::X) << endl; + } break; case Rtn::INFEASIBLE: - cout << "no solution" << endl; + cout << "INFEASIBLE MODEL" << endl; break; case Rtn::UNBOUNDED: - cout << "unbounded solution" << endl; + cout << "UNBOUNDED SOLUTION" << endl; break; default: assert(false); diff --git a/src/solver.cpp b/src/solver.cpp index 200495d..b437bc9 100644 --- a/src/solver.cpp +++ b/src/solver.cpp @@ -32,16 +32,16 @@ Expr operator+(const Expr& x) Expr operator+(Var x, Var y) { Expr exp; - exp.coeffs.resize(std::max(x.col + 1, y.col + 1), 0); - exp.coeffs[x.col] = x.coeffs; - exp.coeffs[y.col] = y.coeffs; + exp.coeffs.resize(std::max(x.get(IntAttr::ModelColumn) + 1, y.get(IntAttr::ModelColumn) + 1), 0); + exp.coeffs[x.get(IntAttr::ModelColumn)] = x.get(DoubleAttr::Coeff); + exp.coeffs[y.get(IntAttr::ModelColumn)] = y.get(DoubleAttr::Coeff); return exp; } Expr operator+(Var x, double a) { Expr exp; - exp.coeffs.resize(x.col + 1); + exp.coeffs.resize(x.get(IntAttr::ModelColumn) + 1); exp.constant = a; return exp; } @@ -85,9 +85,9 @@ Expr operator-(Var x) Expr operator-(Var x, Var y) { Expr exp; - exp.coeffs.resize(std::max(x.col + 1, y.col + 1), 0); - exp.coeffs[x.col] = x.coeffs; - exp.coeffs[y.col] = -y.coeffs; + exp.coeffs.resize(std::max(x.get(IntAttr::ModelColumn) + 1, y.get(IntAttr::ModelColumn) + 1), 0); + exp.coeffs[x.get(IntAttr::ModelColumn)] = x.get(DoubleAttr::Coeff); + exp.coeffs[y.get(IntAttr::ModelColumn)] = -y.get(DoubleAttr::Coeff); return exp; } @@ -104,8 +104,8 @@ Expr operator-(double a, Var x) Expr operator*(double a, Var x) { Expr exp; - exp.coeffs.resize(x.col + 1, 0); - exp.coeffs[x.col] = a * x.coeffs; + exp.coeffs.resize(x.get(IntAttr::ModelColumn) + 1, 0); + exp.coeffs[x.get(IntAttr::ModelColumn)] = a * x.get(DoubleAttr::Coeff); return exp; } @@ -147,6 +147,7 @@ Expr operator/(const Expr& x, double a) Model::Model(): cn(0), bn(0), + sense(0), vars(nullptr) { @@ -165,22 +166,50 @@ Var* Model::addVars(int num) return vars; } -void Model::addConstr(const Expr& expr, char sense, double rhs) +void Model::addConstr(const Expr& expr, ConstrOper sense, double rhs) { - matrix.push_back(vector(1, rhs)); // TODO：没比大小，没考虑常数 - matrix.back().insert(matrix.back().end(), expr.coeffs.begin(), expr.coeffs.end()); + + if (sense == ConstrOper::LESS_EQUAL) { + matrix.push_back(vector(1, rhs - expr.constant)); + matrix.back().insert(matrix.back().end(), expr.coeffs.begin(), expr.coeffs.end()); + } + else if (sense == ConstrOper::GREATER_EQUAL) { + matrix.push_back(vector(1, expr.constant - rhs)); + for (int coeff : expr.coeffs) { + matrix.back().push_back(-coeff); + } + } + else { + addConstr(expr, ConstrOper::LESS_EQUAL, rhs); + addConstr(expr, ConstrOper::GREATER_EQUAL, rhs); + } + for (int c = matrix.back().size(); c <= cn; c++) { matrix.back().push_back(0); } } -void Model::setObjective(Expr obje, int sense) +void Model::setObjective(Expr obje, int _sense) { - matrix.insert(matrix.begin(), obje.coeffs); - matrix.front().insert(matrix.front().begin(), obje.constant); // TODO：写的不优雅 - for (int c = matrix.front().size(); c <= cn; c++) { - matrix.front().push_back(0); + assert(_sense == 1 || _sense == -1); + if (sense == 0) { + matrix.insert(matrix.begin(), obje.coeffs); + matrix.front().insert(matrix.front().begin(), -obje.constant); + for (int c = obje.coeffs.size() + 1; c <= cn; c++) { + matrix.front().push_back(0); + } } + else { + matrix.front().front() = -obje.constant; + for (int c = 0; c < obje.coeffs.size(); c++) { + matrix.front()[c + 1] = obje.coeffs[c]; + } + for (int c = obje.coeffs.size() + 1; c <= cn; c++) { + matrix.front()[c] = 0; + } + + } + sense = _sense; } void Model::print() @@ -193,7 +222,7 @@ void Model::print() } } -Rtn Model::optimize(bool isMax) +Rtn Model::optimize() { bn = matrix.size(); for (int row = 1; row < bn; row++) { @@ -203,12 +232,10 @@ Rtn Model::optimize(bool isMax) } } cn = matrix.front().size(); - if (!isMax) { - for (int i = 0; i < cn; i++) { - matrix.front()[i] = -matrix.front()[i]; - } + for (int i = 0; i < cn; i++) { + matrix.front()[i] = sense * matrix.front()[i]; } - + for (size_t i = 1; i < bn; i++) { basic.push_back(cn - bn + i); } @@ -246,9 +273,9 @@ Rtn Model::optimize(bool isMax) auto iter = find(basic.begin(), basic.end(), cn); if (iter != basic.end()) { for (int col = 1; col < matrix.front().size(); col++) { - if (matrix.front()[col] != 1) { - t = make_pair(iter - basic.begin(), col); - _pivot(t); + if (fabs(matrix.front()[col]) > 1e-10) { + t = make_pair(iter - basic.begin() + 1, col); + _gaussian(t); break; } } @@ -270,15 +297,25 @@ Rtn Model::optimize(bool isMax) } } - res = isMax ? _simplex() : _simplex(); + res = -sense * _simplex(); + for (int row = 1; row < bn; row++) { + if (basic[row - 1] - 1 < cn - bn) { + vars[basic[row - 1] - 1].val = matrix[row].front(); + } + } return rtn; } -double Model::get() +double Model::get(DoubleAttr attr) { return res; } +int Model::get(IntAttr attr) +{ + return cn - bn; +} + double Model::_simplex() { pair t; @@ -289,7 +326,6 @@ double Model::_simplex() } _gaussian(t); } - return matrix.front().front(); } @@ -381,6 +417,7 @@ void Expr::operator+=(const Expr& expr) for (int c = 0; c < expr.coeffs.size(); c++) { coeffs[c] += expr.coeffs[c]; } + constant += expr.constant; } void Expr::operator-=(const Expr& expr) @@ -389,6 +426,7 @@ void Expr::operator-=(const Expr& expr) for (int c = 0; c < expr.coeffs.size(); c++) { coeffs[c] -= expr.coeffs[c]; } + constant -= expr.constant; } void Expr::operator*=(double mult) @@ -396,6 +434,7 @@ void Expr::operator*=(double mult) for (int c = 0; c < coeffs.size(); c++) { coeffs[c] *= mult; } + constant *= mult; } void Expr::operator/=(double a) @@ -403,6 +442,7 @@ void Expr::operator/=(double a) for (int c = 0; c < coeffs.size(); c++) { coeffs[c] /= a; } + constant /= a; } Expr Expr::operator+(const Expr& rhs) @@ -411,6 +451,7 @@ Expr Expr::operator+(const Expr& rhs) for (int c = 0; c < rhs.coeffs.size(); c++) { coeffs[c] += rhs.coeffs[c]; } + constant += rhs.constant; return *this; } @@ -420,5 +461,22 @@ Expr Expr::operator-(const Expr& rhs) for (int c = 0; c < rhs.coeffs.size(); c++) { coeffs[c] -= rhs.coeffs[c]; } + constant -= rhs.constant; return *this; } + +double Var::get(DoubleAttr attr) +{ + switch (attr) { + case DoubleAttr::Coeff: + return coeffs; + case DoubleAttr::X: + return val; + } + return val; +} + +int Var::get(IntAttr attr) +{ + return col; +}