Talk with spookywooky

最后更新于：2022年6月25日 下午

Problem: mobilization

This problem aims to compute \(\displaystyle \max (\sum h_i t_i)(\sum p_i t_i)\) with constraint \(\displaystyle \sum c_i t_i \leq d, \; t_i \geq 0\), is a classical nonlinear programming.(note that " You can acquire any combination of the troop types, even fractional, such that the total cost is no more than the amount you have to spend" in this problem, so \(x_i\) can be any non-negtive "real" number).

It is impossible to have a algorithm dealing with all nonlinear problem, we can only deal with it one (type) by one (type).

For this particular problem, Convex Hull Method or Convex Hull trick ("斜率优化"(Slope optimization) in Chinese) may be a good choose.

First we simply the above problem as \(\displaystyle \max (\sum a_i x_i)(\sum b_i x_i)\) with constraint \(\displaystyle \sum x_i \leq 1\).

where \(a_i = \frac{h_i d}{c_i}\), \(b_i = \frac{p_i d}{c_i}\), \(x_i = \frac{t_i c_i}{d}\).

Note that all \((\sum a_i x_i, \sum b_i x_i)\) form a Convex Hull of \((a_i, b_i)\), so the max point must belong to the side of Convex Hull. So just get the Convex Hull and computer on each side and we get the answer.

The convex function defined on the bounded closed convex set must obtain the maximum value at the boundary

Solution Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#include <bits/stdc++.h>
#define watch(x) std::cout << (#x) << " is " << (x) << std::endl
using LL = long long;
using pii = std::pair<int, int>;
using pll = std::pair<LL, LL>;
using pdd = std::pair<double, double>;

const double eps = 1e-12;
int dcmp(double x) {
  return fabs(x) < eps ? 0 : (x > 0 ? 1 : -1);
}
bool crossLeft(const pii &op, const pii &sp, const pii &ep) {
  return (sp.first - op.first) * (ep.second - op.second)
  < (sp.second - op.second) * (ep.first - op.first);
}
std::vector<pdd> convexHull(std::vector<pdd> p) {
  std::sort(p.begin(), p.end());
  p.erase(std::unique(p.begin(), p.end()), p.end());
  int n = p.size();
  std::vector<pdd> q(n + 1);
  int top = 0;
  for (int i = 0; i < n; ++i) {
    while (top > 1 && crossLeft(q[top - 1], p[i], q[top - 2])) --top;
    q[top++] = p[i];
  }
  int len = top;
  for (int i = n - 2; i >= 0; --i) {
    while (top > len && crossLeft(q[top - 1], p[i], q[top - 2])) --top;
    q[top++] = p[i];
  }
  top -= n > 1;
  q.resize(top);
  return q;
}
int main() {
  // freopen("in","r",stdin);
  std::ios::sync_with_stdio(false);
  std::cin.tie(nullptr);
  int n, d, c;
  std::cin >> n >> d;
  std::vector<pdd> p(n);
  for (auto &[x, y] : p) {
    std::cin >> c >> x >> y;
    x = x * d / c;
    y = y * d / c;
  }
  auto q = convexHull(p);
  auto cal = [](pdd p) {
    return p.first * p.second;
  };
  auto deal = [](pdd a, pdd b) -> double {
    a.first -= b.first; a.second -= b.second;
    if (a.first * a.second < 0) {
      double t = -(b.first / a.first + b.second / a.second) / 2;
      if (t > 0 && t < 1) return (a.first * t + b.first) * (a.second * t + b.second);
    }
    return 0;
  };
  double r = -1;
  for (auto &x : q) r = std::max(r, cal(x));
  for (int i = 0; i != q.size(); ++i) {
    r = std::max(r, deal(q[i], q[(i + 1) % q.size()]));
  }
  std::cout.precision(12);
  std::cout << std::fixed << r << std::endl;
  return 0;
}

Problem: Cheese, If You Please

This problem aims to compute \(\displaystyle \max \sum t_i x_i\) with constraint \(x_i \geq 0\) and \[ \begin{pmatrix} p_{11} & \cdots & p_{1n} \\ \vdots & \ddots & \\ p_{m1} & \cdots & p_{mn} \end{pmatrix} \begin{pmatrix} x_1 \\ \vdots \\ x_n \end{pmatrix} \leq \begin{pmatrix} w_1 \\ \vdots \\ w_n \end{pmatrix} \] , is a classical linear programmingby Simplex algorithm. You can find it in wikipedia or Code template here (sorry I only know this place which writed in chinese, you can find the code nevertheless)

Solution code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
#include <bits/stdc++.h>
#define watch(x) std::cout << (#x) << " is " << (x) << std::endl
using LL = long long;
using pii = std::pair<int, int>;
using pll = std::pair<LL, LL>;

using VD = std::vector<double>;
const double eps = 1e-10;
const double inf = 1e10;
// make sure that A = (I, A') and b >= 0, compute max cx
VD simplexCore(VD c, std::vector<VD> A, VD b) {
  int n = A.size(), m = c.size();
  std::vector<int> p(m);
  std::iota(p.begin(), p.end(), 0);
  for (int i = 0; i < n; ++i) A[i].emplace_back(b[i]);
  c.emplace_back(0);
  A.emplace_back(c);
  for (int j = n; j <= m; ++j) {
    for (int i = 0; i < n; ++i) {
      A[n][j] -= A[n][i] * A[i][j];
    }
  }
  auto check = [&]() -> bool {
    for (int j = n; j < m; ++j) if (A[n][j] > eps) {
      bool flag = false;
      for (int i = 0; i < n; ++i) if (A[i][j] > eps) {
        flag = true;
        break;
      }
      if (!flag) return false;
    }
    return true;
  };
  while (1) {
    int ch = std::max_element(A[n].begin() + n, A[n].begin() + m) - A[n].begin(), hc;
    if (A[n][ch] < eps) break;
    if (!check()) {
      std::cout << "\nNo max solution!\n";
      return VD();
    }
    double theta = DBL_MAX;
    for (int i = 0; i < n; ++i) if (A[i][ch] > eps && A[i].back() / A[i][ch] < theta) {
      theta = A[i].back() / A[i][ch];
      hc = i;
    }
    std::swap(p[ch], p[hc]);
    double tmp = 1 / A[hc][ch];
    for (int j = n; j <= m; ++j) A[hc][j] *= tmp;
    for (int i = 0; i <= n; ++i) if (i != hc) {
      for (int j = n; j <= m; ++j) if (j != ch) {
        A[i][j] -= A[i][ch] * A[hc][j];
      }
    }
    for (int i = 0; i <= n; ++i) A[i][ch] *= -tmp;
    A[hc][ch] = tmp;
  }
  VD x(m);
  for (int i = 0; i < n; ++i) x[p[i]] = A[i].back();
  // watch(-A.back().back()); // max_val
  return x; // point Corresponds to max_val
}
// compute max cx, with Aqx = bq and Alq x <= blq, end of 0 can be ommit in A and Aq
VD simplex(VD c, std::vector<VD> Aq, VD bq, std::vector<VD> Alq, VD blq) {
  assert(Aq.size() == bq.size());
  assert(Alq.size() == blq.size());
  int n = Aq.size() + Alq.size();
  int m = c.size();
  for (int i = 0; i < bq.size(); ++i) if (bq[i] < -eps) {
    for (auto &x : Aq[i]) x = -x;
    bq[i] = -bq[i];
  }
  for (int i = 0; i < blq.size(); ++i) if (blq[i] < -eps) {
    for (auto &x : Alq[i]) x = -x;
    ++m;
  }
  std::vector<VD> A(n, VD(n + m));
  VD f(n + m), b(n);
  int now = n + c.size();
  for (int i = 0; i < n; ++i) A[i][i] = 1;
  for (int i = 0; i < Aq.size(); ++i) {
    for (int j = 0; j < Aq[i].size(); ++j) A[i][n + j] = Aq[i][j];
    b[i] = bq[i];
    f[i] = -inf;
  }
  for (int i = 0; i < Alq.size(); ++i) {
    for (int j = 0; j < Alq[i].size(); ++j) A[i + Aq.size()][n + j] = Alq[i][j];
    if (blq[i] < -eps) {
      A[i + Aq.size()][now++] = -1;
      f[i + Aq.size()] = -inf;
    }
    b[i + Aq.size()] = fabs(blq[i]);
  }
  for (int i = 0; i < c.size(); ++i) f[n + i] = c[i];
  auto x = simplexCore(f, A, b);
  return VD(x.begin() + n, x.begin() + n + c.size());
}
int main() {
  // freopen("in","r",stdin);
  std::ios::sync_with_stdio(false);
  std::cin.tie(nullptr);
  int n, m;
  std::cin >> n >> m;
  std::vector<std::vector<double>> A(n, std::vector<double>(m, 0));
  std::vector<double> b(n), c(m);
  for (auto &x : b) std::cin >> x;
  for (int j = 0; j < m; ++j) {
    for (int i = 0; i < n; ++i) std::cin >> A[i][j];
    std::cin >> c[j];
  }
  auto x = simplex(c, std::vector<VD>(), VD(), A, b);
  double ans = 0;
  for (int i = 0; i < m; ++i) ans += c[i] * x[i];
  std::cout.precision(2);
  std::cout << std::fixed << ans * 100 << std::endl;
  return 0;
}