纯模拟退火题解
jia_shengyuan · · 题解
前言
考场上T1脑抽花了1h,导致这题没时间观察性质,没发现交换方差,只能写了个朴素退火,
这是一篇纯模拟退火的题解,相较于其他题解应该是很简洁的。
40\ pts 做法
按题目要求每次随机一个数进行操作,模拟退火即可。
80\ pts 做法
首先观察性质,
观察到交换差分的性质后,无论是推样例还是数感,都不难感觉到方差最小时,差分应该先减后增(你可以认为是让尽可能多的数集中,从而减小方差)。
因此,我们可以先构造出一种差分先减后增的初始解,之后每次随机两个位置并交换差分,在此基础上再模拟退火,即可获得
100\ pts 做法
我们可以把交换差分改为随机把一个差分插入到另一个地方。观察到许多插入操作的最优方案是一定的。例如对于差分序列 3 1 0 2 4
,若要移动最左测的
因此,我们每次可以随机一个值,并把它放到差分拐点的另一侧,使得差分序列仍然满足先减后增的性质,在此基础上模拟退火,并调出较好的参数,即可通过本题。
由于涉及到插入与删除,且有吸氧,我们不妨使用STL的 deque
简化代码。
我们可以用样例
#include <algorithm>
#include <climits>
#include <cstdio>
#include <random>
#include <deque>
#include <ctime>
#include <cmath>
#include <cassert>
#define ll long long
// #define err printf
#define err(...)
const int maxn = 1e4+9;
using namespace std;
typedef deque<int>::iterator iter;
int a[maxn],b[maxn],mini=INT_MAX,n;
deque<int> cf;//差分
iter lpt,rpt;
ll ans;
mt19937 rnd(time(0));
void Read(){
scanf("%d",&n);
for(int i=0; i<n; i++){
scanf("%d",&a[i]);
a[i]*=n; // 为防止出现小数,开始时把每一项乘以n即可
}
for(int i=0; i<n-1; i++){
b[i]=a[i+1]-a[i];
}
}
void Init(){
sort(b,b+n-1);
bool f=0;
// 构造初始解
for(int i=0; i<n-1; i++){
mini = min(mini,b[i]);
if(f) cf.push_front(b[i]);
else cf.push_back(b[i]);
f = !f;
}
// 计算最小值有几个,用于退火时判断当前数在拐点左侧还是右侧
f=0;
for(int i=0; i<n-1; i++){
if(cf[i]==mini){
if(!f){
lpt = cf.begin()+i;
rpt = cf.begin()+i;
f=1;
} else {
rpt = cf.begin()+i;
}
}
}
}
ll pf(ll x) { return x*x; }
// 计算方差。由于整体加减不影响方差,我们默认第一项是0
ll calc(){
ll sum=0,cur=0,pj,fc;
for(int solo : cf){
cur += solo;
sum += cur;
}
pj = sum/n;
fc = pf(pj);//a[0]的贡献
cur = 0;
for(int solo : cf){
cur += solo;
fc += pf(cur-pj);
}
return fc/n;
}
void SA(){
double T = 2e5;//实测这个初温有很好的出解率
ll now=ans,newans,de;
int rndsize = n-1-(rpt-lpt+1),val,len=rpt-lpt+1;
iter p1;
while(T>1e-5){
// 处理拐点区域
for(int i=0; i<n-1; i++){
if(cf[i] == mini){
lpt = cf.begin()+i;
rpt = lpt+len-1;
break;
}
}
deque<int> tmp = cf;
// 随机选一个数插入到别的地方
p1 = cf.begin()+(rnd()%rndsize);
if(p1 >= lpt){ //从右侧插入到左侧
p1 += len;
val = *p1;
cf.erase(p1);
iter now = cf.begin();
while(*now > val) now++;
cf.insert(now,val);
} else { // 从左侧插入到右侧
val = *p1;
cf.erase(p1);
iter now = cf.end();
while(*(now-1) > val) {
now--;
}
cf.insert(now,val);
}
// 更新答案
newans = calc(); de=newans-now;
if(de<0 || exp(-de/T)*UINT_MAX>rnd()) {
ans = min(ans,newans);
now = newans;
} else cf = tmp;
T*=0.97;
}
}
int main(){
Read();
Init();
ans = calc();
if(rpt-lpt+1 == n-1){ //如果差分全一样,直接输出防止诡异错误
printf("%lld",ans);
} else {
// 退火时尽可能把时间用满,增加出解率
clock_t s = clock();
SA();
clock_t len = clock()-s, limit = CLOCKS_PER_SEC*0.975;
while(clock()+len < limit) SA();
printf("%lld",ans);
}
return 0;
}