洛谷P1337 【[JSOI2004]平衡点 / 吊打XXX】(模拟退火)
很可惜,充满Mo力的Mo拟退火并不是正解。不过这是一道最适合开始入手Mo拟退火的好题。
对模拟退火还不是很清楚的可以看一下
这道题还真和能量有点关系。达到平衡稳态的时候,物体的总能量应该是最小的。而总的能量来源于每个物体的重力势能之和。要想让某个物体势能减小,那就让拉着它的绳子在桌面下方的长度尽可能的长,也就是桌面上的要尽可能短。由此看来,某个物体的势能与桌面上的绳子的长度、物体重量都成正比。
于是,为了找到平衡点,我们要找一个点使得\(\sum_{i=1}^n d_i*w_i\)最小(\(d_i\)为\(i\)点到该点的距离)。
函数已经求出来了,接下来就是正常的模拟退火过程。注意一些细节就好了。
首先,初始解可以设成\(({\sum_{i=1}^n x_i\over n},{\sum_{i=1}^n y_i\over n})\),可以更接近正解
解变动值最好随机两个值,\(\Delta x\)和\(\Delta y\)。随机变动距离和角度可能常数有点大。
然后是写法问题。蒟蒻又学习了一招,知道有一个常数叫RAND_MAX,用于生成\([0,1)\)的随机值还是很方便的,在不同机子下可移植性也很强。(难怪Windows上rand出来的都是短整形数)
另外,因为退火的时候有可能本来找到了最优解,但后来接受了较劣解并逐步降温稳定了,那么也就是说本来找到的最优解到最后却给弄丢了。
为了避免这种情况,蒟蒻参考YL的做法,设一个全局最优解best,它只接受比自己更优的解。
而且这样还有一个好处,就是进行多次退火的过程中,也能保证best是单调递减的,无论如何也不会变劣。
最后就是调参数的问题了。
发现自己写了个假的模拟退火以后,蒟蒻马上一改,却发现参数又需要重新调。。。。。。(YL据老也背了锅)
拼命的二分,二分。。。。。。参数大得都要TLE了还是不能保证正确率?!
结果小号提交记录第二次刷屏了(第一次还是打表过洛谷愚人节T3。。。。。。)
后来突然猛地一想:\(n,x,y,w\)都不小,乘来乘去不会掉精度了吧?
于是改用long double继续,情况好多了
update:
可能发现在哪里咕咕了,应该不是\(n,x,y,w\)的问题。
随机变动距离写的是T*(rand()*2-RAND_MAX)
,也就意味着生成了一个\([-T*RAND\_MAX,T*RAND\_MAX)\)的数。
这没有错,因为RD的值域确实与T成正比,但不对劲的是它太大了,难道是因为这里爆double精度?
可是看到其它巨佬用double没出锅诶。。。
其实还有一个优化就是,退火结束后从best处向小范围内rand一些点,rand个千把次,找最优的一个点作为答案。
这样比多次退火效率高了不少,就更容易弥补之前没找到最优解的遗憾了。
听XZY大佬的讲课学的,stO XZY Orz
代码蒟蒻就懒得重写了。。。。。。(连个理由都懒得找)
下面的参数可供参考,因为D算是很小的了,所以跑得快多了,只要28ms
#include<cstdio>
#include<cmath>
#include<ctime>
#include<cstdlib>
#define RG register
#define R RG long double
#define RD T*(rand()*2-RAND_MAX)//生成一个[-T*RAND_MAX,T*RAND_MAX)的随机变动距离
//以前写的是T*((double)rand()/RAND_MAX*2-1),总是WA(可能是先除后乘掉精度了),然后就参考了YL的写法
const int N=1009;
const long double D=0.97,EPS=1e-14;//更新后二分出了保证极高正确率和较高效率的参数
double x[N],y[N],w[N];
int n;
inline long double calc(R x0,R y0){
R res=0,dx,dy;
for(RG int i=1;i<=n;++i){//函数求值
dx=x[i]-x0;dy=y[i]-y0;
res+=sqrt(dx*dx+dy*dy)*w[i];
}
return res;
}
int main(){
R T,x0,y0,x1,y1,res,ans,best,bx=0,by=0;
RG int i,times=1;//一次的正确率很高的,不放心也可以调大
scanf("%d",&n);
for(i=1;i<=n;++i){
scanf("%lf%lf%lf",&x[i],&y[i],&w[i]);
bx+=x[i];by+=y[i];
}//初始横纵坐标均选择平均值
best=ans=calc(bx/=n,by/=n);
srand(time(NULL));
while(times--/*clock()<CLOCKS_PER_SEC*0.9*/){//比赛的时候试试注释里写的,卡时还是靠谱些
ans=best;x0=bx;y0=by;
for(T=100000;T>EPS;T*=D){
x1=x0+RD;y1=y0+RD;
res=calc(x1,y1);
if(best>res)
best=res,bx=x1,by=y1;//更新最优答案
if(ans>res||exp((ans-res)/T)>(long double)rand()/RAND_MAX)
ans=res,x0=x1,y0=y1;//接受新解
}
};
printf("%.3Lf %.3Lf\n",bx,by);
return 0;
}