[BZOJ2876][Noi2012]骑行川藏（二分答案+拉格朗日乘数法）

Address

https://www.lydsy.com/JudgeOnline/problem.php?id=2876

Solution

为了方便某些萌新，讲下拉格朗日乘数法。
作用：求多元函数 $f(x_1,x_2,...,x_n)$ 在约束 $\phi(x_1,x_2,...,x_n)=0$ 下的极值。引入变量 $\lambda$ ，令拉格朗日函数：

$$F(x_1,x_2,...,x_n,\lambda)=f(x_1,x_2,...,x_n)+\lambda\phi(x_1,x_2,...x_n)$$
对 $F$ 的每个变量（包括 $\lambda$ ）分别求偏导。
（事实上对 $\lambda$ 求偏导后的结果就是 $\phi(x_1,x_2,...,x_n)$ ）
（多元函数 $F$ 对一个自变量 $x$ 的偏导数，定义为将 $F$ 除 $x$ 之外的自变量全部视为常数后对 $F$ 求导）
当 $F$ 对每个自变量的偏导数都等于 $0$ 时，
$f(x_1,x_2,...,x_n)$ 在约束 $\phi(x_1,x_2,...,x_n)=0$ 下取得极值。
回到原问题。
发现体力不是一个标准的约束。
但如果体力没用完，则一定可以在某一天速度更快。
所以， 存在一种最优方案使得体力恰好用完。
同样我们也能得到，存在一种最优方案使得对于每个 $i$ 有 $v_i\ge v_i'$ 。
（ $v_i$ 表示第 $i$ 天实际速度）
目标函数：
$$f(v_1,v_2,...,v_n)=\sum_{i=1}^n\frac{s_i}{v_i}$$
约束条件：
$$\phi(v_1,v_2,...,v_n)=-E_U\sum_{i=1}^nk_i(v_i-v_i')^2s_i=0$$
引入系数 $\lambda$ ，构造拉格朗日函数：
$$\sum_{i=1}^n(\frac{s_i}{v_i}+\lambda k_i(v_i-v_i')^2s_i)-E_U$$
对 $v_i$ 求偏导得到：
$$-\frac{s_i}{v_i^2}+\lambda k_i(2v_i-2v_i'+v_i'^2)s_i$$
要求上式为 $0$ 。
$$-\frac 1{v_i^2}+\lambda k_i(2v_i-2v_i')=0$$
$$2\lambda k_i(v_i-v_i')v_i^2=1$$
易得上式左边在 $v_i\ge v_i'$ 即上式左边的值大于等于 $0$ 时单调递增。
故在 $\lambda$ 确定时二分答案可以算出 $v_i$ 。
同时， $\lambda$ 越大， $v_i$ 越小。
所以 $\lambda$ 的值也可以二分，通过 $\phi(v_1,v_2,...,v_n)$ 的值来调整 $\lambda$ 。
最后得出了 $\lambda$ 的值之后就能计算取到极值时的 $v_i$ 。

Code

#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define For(i, a, b) for (i = a; i <= b; i++)
using namespace std;
const int N = 1e4 + 5;
const double eps = 1e-12;
int n;
double eu, k[N], s[N], vp[N], ans;
double vsol(double lambda, int i) {
    double l = max(vp[i], 0.0), r = 1e5;
    while (r - l > eps) {
        double mid = (l + r) / 2, val;
        val = 2.0 * lambda * k[i] * (mid - vp[i]) * mid * mid;
        if (val < 1) l = mid;
        else r = mid;
    }
    return (l + r) / 2;
}
double calc(double lambda) {
    int i; double res = 0, tmp;
    For (i, 1, n) tmp = vsol(lambda, i),
        res += k[i] * s[i] * (tmp - vp[i]) * (tmp - vp[i]);
    return res;
}
int main() {
    int i;
    double x, y, lambda, l = 0, r = 1e5;
    cin >> n >> eu;
    For (i, 1, n) scanf("%lf%lf%lf", &s[i], &k[i], &vp[i]);
    while (r - l > eps) {
        double mid = (l + r) / 2;
        if (calc(mid) > eu) l = mid;
        else r = mid;
    }
    For (i, 1, n) ans += s[i] / vsol((l + r) / 2, i);
    printf("%.10lf\n", ans);
    return 0;
}