gc_pwl_func_c++.cpp
/ * 2019年版权,Gurobi优狗万app足彩化,LLC这个例子考虑以下非凸非线性问题最大化2 x + y受制于exp (x) + 4倍根号(y) < = 9 x, y > = 0,我们给你两种方法解决:1)使用分段线性方法来处理一般函数约束(如exp和sqrt)。a)将两个变量u = exp(x) v = sqrt(y)计算出u = exp(x)的一些步长(例如x = 0,1e - 3,2e -3，…)对于步长(例如，y = 0,1e - 3,2e -3，…)的点(y, v)ymax)。我们需要计算xmax和ymax(这在本例中很容易，但在一般情况下并不适用)。c)利用这些点添加两个分段线性类型的一般约束。2)直接使用Gurobis内置的通用函数约束(EXP和POW)。在这里，我们不需要计算点和最大可能值，这将由Gurobi在内部完成。在这种方法中，我们将展示如何“放大”最佳解决方案并收紧公差以提高解决方案的质量。*/ #if define (WIN32) || define (WIN64) #include #endif #include "gurobi_c++.h" #include using namespace std;Static double f(double u) {return exp(u);} static double g(double u) {return sqrt(u); } static void printsol(GRBModel& m, GRBVar& x, GRBVar& y, GRBVar& u, GRBVar& v) { cout << "x = " << x.get(GRB_DoubleAttr_X) << ", u = " << u.get(GRB_DoubleAttr_X) << endl; cout << "y = " << y.get(GRB_DoubleAttr_X) << ", v = " << v.get(GRB_DoubleAttr_X) << endl; cout << "Obj = " << m.get(GRB_DoubleAttr_ObjVal) << endl; // Calculate violation of exp(x) + 4 sqrt(y) <= 9 double vio = f(x.get(GRB_DoubleAttr_X)) + 4 * g(y.get(GRB_DoubleAttr_X)) - 9; if (vio < 0.0) vio = 0.0; cout << "Vio = " << vio << endl; } int main(int argc, char* argv[]) { double* xpts = NULL; double* ypts = NULL; double* vpts = NULL; double* upts = NULL; try { // Create environment GRBEnv env = GRBEnv(); // Create a new model GRBModel m = GRBModel(env); // Create variables double lb = 0.0, ub = GRB_INFINITY; GRBVar x = m.addVar(lb, ub, 0.0, GRB_CONTINUOUS, "x"); GRBVar y = m.addVar(lb, ub, 0.0, GRB_CONTINUOUS, "y"); GRBVar u = m.addVar(lb, ub, 0.0, GRB_CONTINUOUS, "u"); GRBVar v = m.addVar(lb, ub, 0.0, GRB_CONTINUOUS, "v"); // Set objective m.setObjective(2*x + y, GRB_MAXIMIZE); // Add linear constraint m.addConstr(u + 4*v <= 9, "l1"); // Approach 1) PWL constraint approach double intv = 1e-3; double xmax = log(9.0); int len = (int) ceil(xmax/intv) + 1; xpts = new double[len]; upts = new double[len]; for (int i = 0; i < len; i++) { xpts[i] = i*intv; upts[i] = f(i*intv); } GRBGenConstr gc1 = m.addGenConstrPWL(x, u, len, xpts, upts, "gc1"); double ymax = (9.0/4.0)*(9.0/4.0); len = (int) ceil(ymax/intv) + 1; ypts = new double[len]; vpts = new double[len]; for (int i = 0; i < len; i++) { ypts[i] = i*intv; vpts[i] = g(i*intv); } GRBGenConstr gc2 = m.addGenConstrPWL(y, v, len, ypts, vpts, "gc2"); // Optimize the model and print solution m.optimize(); printsol(m, x, y, u, v); // Approach 2) General function constraint approach with auto PWL // translation by Gurobi // restore unsolved state and get rid of PWL constraints m.reset(); m.remove(gc1); m.remove(gc2); m.update(); GRBGenConstr gcf1 = m.addGenConstrExp(x, u, "gcf1"); GRBGenConstr gcf2 = m.addGenConstrPow(y, v, 0.5, "gcf2"); m.set(GRB_DoubleParam_FuncPieceLength, 1e-3); // Optimize the model and print solution m.optimize(); printsol(m, x, y, u, v); // Zoom in, use optimal solution to reduce the ranges and use a smaller // pclen=1e-5 to solve it double xval = x.get(GRB_DoubleAttr_X); double yval = y.get(GRB_DoubleAttr_X); x.set(GRB_DoubleAttr_LB, max(x.get(GRB_DoubleAttr_LB), xval-0.01)); x.set(GRB_DoubleAttr_UB, min(x.get(GRB_DoubleAttr_UB), xval+0.01)); y.set(GRB_DoubleAttr_LB, max(y.get(GRB_DoubleAttr_LB), yval-0.01)); y.set(GRB_DoubleAttr_UB, min(y.get(GRB_DoubleAttr_UB), yval+0.01)); m.update(); m.reset(); m.set(GRB_DoubleParam_FuncPieceLength, 1e-5); // Optimize the model and print solution m.optimize(); printsol(m, x, y, u, v); } catch(GRBException e) { cout << "Error code = " << e.getErrorCode() << endl; cout << e.getMessage() << endl; } catch(...) { cout << "Exception during optimization" << endl; } if (xpts) delete[] xpts; if (ypts) delete[] ypts; if (upts) delete[] upts; if (vpts) delete[] vpts; return 0; }
