【拉格朗日乘数法】[NOI2012]骑行川藏

题目描述

分析

首先，我们来考虑一下部分分。

$N=1$

直接算即可，$v_1=\sqrt{\frac{E}{k_1\times s_1}}+v_1'\\T=\frac{s_1}{v_1}$

$N=2$

$v_1=\sqrt{\frac{E_1}{k_1\times s_1}}+v_1'\\v_2=\sqrt{\frac{E-E_1}{k_2\times s_2}}+v_2' \\T=\frac{s_1}{v_1}+\frac{s_2}{v_2}$
我们猜测$T$关于$E_1$的函数是一个单谷函数，对$E_1$三分一下就可以了

$N\le 100$&$N\le 1000$

我不会QAQ

$N\le 10000$

首先，我们将题目要求转化为数学式子

$$E\ge\sum_{i=1}^N k_is_i(v_i-v_i')^2\\ T=\sum_{i=1}^N \frac{s_i}{v_i}$$
求 $T_{min}$
基于贪心的思想，显然 $E=\sum_{i=1}^N k_is_i(v_i-v_i')^2$时，答案最优。

拉格朗日乘数法

具体证明我不会。。。。
例题

在$g(x1, x2, x3, \dots, xn) = c$的约束条件下，求$f(x1, x2, x3, \dots, xn)$的最值。

我们先来讲一下结论。

结论

当$f(x1, x2, x3, \dots, xn)$取到最值时，$f(x1, x2, x3, \dots, xn)$和$g(x1, x2, x3,\dots, xn)=c$的法向量（也就是梯度向量）平行。

我们将例题具体化，比如$f(x,y)=x^2+y^2\\g(x,y)=xy=3$，求$f(x,y)_{min}$
先来肉眼观察，很容易发现$f(x,y)_{min}=6$
然后我们用拉格朗日乘数法试试。
设

$$\nabla f(x,y)=\lambda\nabla g(x,y)$$
我们就可以列出方程 $$\begin{cases}\begin{align} 2x&=\lambda y \\ 2y&=\lambda x \\xy&=3\end{align}\end{cases}$$
解得 $$\lambda=\pm 2$$
所以 $\begin{align}\begin{cases}x&=\sqrt 3\\y&=\sqrt 3\end{cases} or \begin{cases}x&=-\sqrt3 \\y&=-\sqrt 3\end{cases}\end{align}$
此时 $f(x,y)=6$，和我们观察出的答案一致。

怎么求多元函数梯度
我也不太会，但是如果这个多元函数可以分解为多个一元函数相加的话，我们求每个一元函数的导数，然后组成的向量就是多元函数的梯度了。比如$\nabla f(x,y)=(2x,2y)$，这对于这道题来说已经够用了。

使用拉格朗日乘数法解决这道题

我们再来看那两个式子

$$E=\sum_{i=1}^N k_is_i(v_i-v_i')^2\\ T=\sum_{i=1}^N \frac{s_i}{v_i}$$
就是要求在 $E=\sum_{i=1}^N k_is_i(v_i-v_i')^2$的约数条件下， $T_{min}$
套用拉格朗日乘数法，列出方程。
$$\begin{align}\begin{cases}\frac{-s_1}{v_1^2}&=2\lambda k_1s_1(v_1-v_1')\\ \frac{-s_2}{v_2^2}&=2\lambda k_2s_2(v_2-v_2')\\\cdots\\\frac{-s_n}{v_n^2}&=2\lambda k_ns_n(v_n-v_n')\end{cases}\end{align}$$
移一下项，可以得到 $2\lambda k_is_i(v_i-v_i')v_i^2=-1$
由于 $\begin{cases}v_i\gt0\\v_i\gt v_i'\end{cases}$，显然可以得到 $\lambda \lt 0$。
$2 k_is_i(v_i-v_i')v_i^2=\frac{-1}{\lambda}$， $v_i$随 $\lambda$的增大而增大。
令 $h(v_i)=2\lambda k_is_i(v_i-v_i')v_i^2$， $h(v_i)$随 $v_i$的增大而减小。

我们可以先二分$\lambda$的大小，再根据$h(v_i)=-1$二分求出$v_i$的值，再根据此时的$E$来调整$\lambda$的大小。

代码

#include<cstdio>
#include<cmath>
#include<algorithm>
#define MAXN 10000
#define EPS 1e-7
using namespace std;
int n;
double s[MAXN+10],k[MAXN+10],v[MAXN+10],vv[MAXN+10],E,ans;
inline double sqr(double x){
    return x*x;
}
inline double cal(int i,double vi,double lambda){
    return vi*vi*2*lambda*k[i]*(vi-v[i]);
}
double Divide_Conqure(int i,double l,double r,double lambda){
    double mid;
    int t=0;
    while(++t<60){
        mid=(l+r)/2;
        if(cal(i,mid,lambda)<-1)
            r=mid;
        else
            l=mid;
    }
    return l;
}
double check(double lambda){
    int i;
    double En=0;
    for(i=1;i<=n;i++){
        vv[i]=Divide_Conqure(i,max(v[i],0.0),1e3,lambda);
        En+=k[i]*s[i]*sqr(vv[i]-v[i]);
    }
    return En;
}
double Divide_Conqure(double l,double r){
    double mid;
    int t=0;
    while(++t<70){
        mid=(l+r)/2;
        if(check(mid)<E)
            l=mid;
        else
            r=mid;
    }
    return l;
}
//　f(v1, v2, v3, ..., vn) = ∑ si / vi 
//　　g(v1, v2, v3, ..., vn) = ∑ kisi(vi-vi')^2
// (si'vi-sivi')/vi^2=-si/vi^2=2*lambda*kisi(vi-vi')
//vi^2*2*lambda*ki(vi-vi')和vi成反比&&=-1 
//vi^3*2*lambda*ki-vi^2*lambda*vi'*2*ki=-1/lambda vi和lambda成正比
int main()
{
    scanf("%d%lf",&n,&E);
    int i;
    for(i=1;i<=n;i++)
        scanf("%lf%lf%lf",s+i,k+i,v+i);
    Divide_Conqure(-1e3,0);
    for(i=1;i<=n;i++)
        ans+=s[i]/vv[i];
    printf("%.8lf\n",ans);
}