模拟退火入门 / P10228 [UVA] A Star not a Tree?

给出 \(N\) 个点,让你求其的费马点 \(Ans\),即 \(Ans\) 到所有其它点的距离总和最短。

很显然位于 \(Ans\) 旁边的点肯定都是距离要大的,且再旁边一点肯定是距离更大的,那么这就是一个单峰的函数,可以用三分法来做。不过我还是把这道题当退火的入门题来做了。

\[\texttt{模拟退火} \]

它的作用是,给你一个集合 \((x_1,x_2,x_3...x_p)\) 所形成的函数 \(y\),其中 \(y\) 是一个没有规律 (或有一些规律) 的函数,模拟退火可以在一定时间内求出 \(y\) 的最值。

简单来说,模拟退火基于一个参数 \(Temp\),当 \(Temp\) 为一个非常小的值的时候,答案就会变得越来越精确且稳定。

定义:

• \(Temp\) : 则现在的温度
• \(Delta\) : 降温系数
• \(Final\_Temp\) : 结束的温度

我们一开始温度可以为一个 \(10^6\),使得温度足够高。每次我们依照原来的 \(x,y\) (现在的最优解) 随机生成一个答案 \(l,r\),如: (\(y\) 同理)

\[l=x+(2 \times (random(ransum))-ransum)*Temp \]

中间的一坨随机数是生成一个 \(-ransum\)~\(ransum\) 的数,\(ransum\) 一般设置为最大的坐标,这里就是 \(10^4\)。

一开始 \(l,r\) 可能会很乱,但是当 \(Temp\) 降到一定程度的时候,\(l,r\) 就可能会和答案非常接近。其次我们需要知道我们是否要这个答案 \(l,r\)。

我们定义答案为 \(Ans\) , \(dix=Ans(l,r)-Ans(x,y)\)。当 \(dix<0\) (也就是最优) 或者
\(dix>=0\),且 \(\exp(-dix/Temp) \times ransum>ransum\) 的时候,也就是有一定的几率会走这一步。其中我们可以知道,当 \(Temp\) 越来越小的时候,这个条件越来越难达成。如图:

最后要降温一下,即 \(Temp=Temp \times Delta\),\(Delta=0.999\) (我认为这样子答案的正确性比较好),当 \(Temp<Final\_Temp\) 的时候停止。
我们说说模拟退火的过程,顺便讲一讲就算答案稍劣也要选择的必要性。

一开始我们在一个任意的点,当然你可以把 \(x,y\) 搞得很接近答案。

其次它会随机一个更低的点跳。

一直跳以后我们会发现它进入了一个并不是很友好的洞。

我们就有一定可能走出一个蓝色的劣点,但这样子可以帮助我们走出劣势,回到应该有的状态。

其次我们会发现模拟退火不一定会到达真正的答案,我们可以多做几遍,或者把 \(Delta\) 调小等各种方法。

在一维的时候,我们可以直接扫过去一遍就可以找出答案。但是在二维甚至 \(p\) 维的时候,时间复杂度爆炸,我们需要用到模拟退火。可以发现的事,\(p\) 越大,答案的准确度越低。

Uses math;

Const
    total=110;

var
    point:array[-1..total,1..2] of longint;
    n,i,test,ransum:longint;
    x,y,ans:extended;

function Contribution(x,y:extended):extended;
var i:longint;
begin
    Contribution:=0;
    for i:=1 to n do Contribution:=Contribution+sqrt(sqr(point[i,1]-x)+sqr(point[i,2]-y));
end;

procedure Simulate_Anneal(var x,y,ans:extended);
var
    l,r,dix,tmp,temp,delta,final_temperature:extended;
begin
    delta:=0.999; temp:=1 << 25;
    final_temperature:=0.00000000001;
    while (temp>=final_temperature) do
    begin
        l:=x+((random(ransum) << 1)-ransum)*temp;
        r:=y+((random(ransum) << 1)-ransum)*temp;
        tmp:=Contribution(l,r); dix:=tmp-ans;
        if dix<0 then begin x:=l; y:=r; ans:=tmp; end else
        if (exp(-dix/temp)*ransum>random(ransum)) then begin x:=l; y:=r; end;
        temp:=temp*delta;
    end;
end;

begin
    read(test);
    while test>0 do 
    begin
        read(n); ans:=maxlongint; randomize;
        for i:=1 to n do
        begin
            read(point[i,1],point[i,2]);
            x:=x+point[i,1]; y:=y+point[i,2];
        end;
        x:=x/n; y:=y/n; ransum:=10000; Simulate_Anneal(x,y,ans);
        writeln(Contribution(x,y):0:0); dec(test);
    end;
    writeln;
end.