#include using namespace std; #define ALL(c) c.begin(), c.end() #define FOR(i,c) for(typeof(c.begin())i=c.begin();i!=c.end();++i) #define REP(i,n) for(int i=0;i= -EPS) return -T[m][n+m]; // optimal if (t < T[q][n+m]) { // tight on c -> primal update REP(j, m) if (T[j][p] >= EPS) if (T[j][p]*(T[q][n+m]-t) >= T[q][p]*(T[j][n+m]-t)) q = j; if (T[q][p] <= EPS) return INF; // primal infeasible } else { // tight on b -> dual update REP(i, n+m+1) T[q][i] *= -1; REP(i, n+m) if (T[q][i] >= EPS) if (T[q][i]*(T[m][p]-t) >= T[q][p]*(T[m][i]-t)) p = i; if (T[q][p] <= EPS) return -INF; // dual infeasible } REP(i, m+n+1) if (i != p) T[q][i] /= T[q][p]; T[q][p] = 1; // pivot(q,p) REP(j, m+1) if (j != q) { double alpha = T[j][p]; REP(i, n+m+1) T[j][i] -= T[q][i] * alpha; } } } #include int doit(int seed) { srand(seed); int n = 2, m = 2; double c[n]; double A[m*n]; double b[m]; c[0] = -1000; c[1] = -2000; REP(i, m) cin >> b[i]; A[0] = 3/4.; A[2] = 1/4.; A[1] = 2/7.; A[3] = 5/7.; /* cout << "c = ["; REP(i,n) cout << c[i] << (i == n-1 ? "];" : ","); cout << endl; cout << "b = ["; REP(i,m) cout << b[i] << (i == m-1 ? "];" : ";"); cout << endl; cout << "A = ["; REP(j, m) { REP(i, n) cout << A[j*n+i] << (i == n-1 ? j == m-1 ? "];" : ";" : ","); } cout << endl; cout << "[sol,opt] = linprog(c,[A;-eye(size(c,2))],[b;zeros(size(c,2),1)])" << endl; cout << endl; */ cout << fixed << setprecision(20) << -simplexMethodPD(c, n, b, m, A) << endl; } int main() { doit( time(0) ); }