「NOI2018模拟7」Yja - 拉格朗日乘数法+二分 | Bill Yang's Blog

「NOI2018模拟7」Yja - 拉格朗日乘数法+二分

题目大意

    在平面上找$n$个点,要求这$n$个点离原点的距离分别为$r_1,r_2,\ldots,r_n$。最大化这$n$个点构成的凸包面积,凸包上的点的顺序任意。
    注意:不要求点全部在凸包上。


题目分析

根本不会做!
疯狂%AChen。

转化模型

首先,用$O(\sum\limits_{i=3}^8i!)$的时间枚举一下构成凸包的点。
注意,由贪心,组成凸包的点的$r$均比较大,因此不需要枚举所有的子排列。

考虑,$n$个点的凸包的面积可以表示为:

其中$\theta_i$是凸包上三角形$\triangle Op_{i}p_{i+1}$的$\angle O$。($O$是原点)
需要满足的限制是:

然后就要套用拉格朗日乘数法了。

拉格朗日乘数法

引用这里面所讲的最简单的形式:

我们先一般化一个二元最优化问题为:

将目标函数$f(x,y)$和等式约束条件$g(x,y)=c$画出来就如下图所示:
20171126104225520.png
其中的$f(x,y)$虚线为等高线,而红线为$g(x,y)=c$这个约束函数曲线与$f(x,y)$的交点的连线在$x−y$平面的映射。其中,假设有$d_3\gt d_2\gt d_1$,$d_1$点为最小值点(最优值点)。从直观上可以发现,在$g(x,y)=c$与$f(x,y)$的非最优化交点$A,B,C,D$上,其$f(x,y)$和$g(x,y)$的法线方向并不是共线的,注意,这个相当关键,因为如果不是共线的,说明$g(x,y)=c$与$f(x,y)$的交点中,还存在可以取得更小值的点存在。对于$A$点来说,$B$点就是更为小的存在。因此,我们从直觉上推论出只有当$g(x,y)=c$与$f(x,y)$的法线共线时,才是最小值点的候选点(鞍点)。推论到多元变量的问题的时候,法线便用梯度表示。于是,我们有原问题取得最优值的必要条件:

其中的$\lambda$表示两个梯度共线。
可以简单的变形为:

让我们去掉梯度算子,得出:

这个时候$\lambda$取个负号也是不影响的,所以上式通常写作:

看!我们得出了我们高数中经常见到的等式约束下的拉格朗日乘数函数的表示方法

此题我们可以设:

然后对$L(\theta_1,\theta_2,\ldots,\theta_n,\lambda)$求偏导:
这里以三维情况举个例:

偏导等于$0$。
故$\lambda+\frac 12r_1r_2\cos(\theta_1)=\lambda+\frac 12r_2r_3\cos(\theta_2)=\lambda+\frac 12r_3r_1\cos(\theta_3)=0$。
改变一下$\lambda$的值无影响,故得到:

二分

显然$0\le\theta_i\le\pi$,故$\cos$函数单调递减,$\lambda$具有单调性。
不妨二分$\lambda$,然后即可算出$\theta_i$,也就可以得到凸包面积了。


代码

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
#include<bits/stdc++.h>
using namespace std;
inline int Get_Int() {
int num=0,bj=1;
char x=getchar();
while(!isdigit(x)) {if(x=='-')bj=-1;x=getchar();}
while(isdigit(x)) {num=num*10+x-'0';x=getchar();}
return num*bj;
}
const double eps=1e-8,pi=acos(-1);
int dcmp(double x) {return fabs(x)<=eps?0:(x>eps?1:-1);}
int n;
double r[15],theta[15],ans=0;
double Cal(double lambda) {
double sum=0;
r[n+1]=r[1];
for(int i=1; i<=n; i++)theta[i]=acos(lambda/(r[i]*r[i+1])),sum+=theta[i];
return sum;
}
int main() {
n=Get_Int();
for(int i=1; i<=n; i++)r[i]=Get_Int();
sort(r+1,r+n+1);
for(; n>=3; n--) {
double L=-r[1]*r[2],R=r[1]*r[2];
do {
double Left=L,Right=R;
while(Right-Left>=eps) {
double mid=(Left+Right)/2;
if(Cal(mid)>=2*pi)Left=mid;
else Right=mid;
}
double mid=(Left+Right)/2;
if(dcmp(Cal(mid)-2*pi)==0) {
double sum=0;
theta[n+1]=theta[1];
for(int i=1; i<=n; i++)sum+=sin(theta[i])*r[i]*r[i+1];
ans=max(ans,sum);
}
} while(next_permutation(r+1,r+n+1));
sort(r+1,r+n+1);
for(int i=1; i<n; i++)r[i]=r[i+1];
}
printf("%0.8lf\n",ans/2);
return 0;
}
本篇文章有用吗?有用还不快点击!
0%