[LT记录]UVA10228 A Star not a Tree?
gcwixsxr
·
·
题解
题目传送门：UVA10228 A Star not a Tree?
题目大意：给定一个 n 边形所有顶点坐标 (x_i,y_i) ，求一个点使得到所有顶点距离和最小。
数据范围：1\leq n\leq 100
似乎题解区里并没有明确给出是单峰函数的证明。那么接下来给出一种证明方法（用到一些偏导数的知识）。
对单峰函数的证明
由题意可知，目标函数是一个二元函数，可以表示为：
f(x,y)=\sum_{i=1}^n\sqrt{(x-x_i)^2+(y-y_i)^2}
那么对 x,y 分别求偏导数，有：
\begin{aligned}
f^\prime_x(x,y)=\sum_{i=1}^n\dfrac {x-x_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}\\
f^\prime_y(x,y)=\sum_{i=1}^n\dfrac {y-y_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}\\
\end{aligned}
然后再求一次一阶偏导数（即最初式子的二阶偏导数）：
\begin{aligned}
f^{\prime\prime}_ {xx}(x,y)=\sum_{i=1}^n\dfrac {(y-y_i)^2}{\left(\sqrt{(x-x_i)^2+(y-y_i)^2}\right)^3}\\
f^{\prime\prime}_ {yy}(x,y)=\sum_{i=1}^n\dfrac {(x-x_i)^2}{\left(\sqrt{(x-x_i)^2+(y-y_i)^2}\right)^3}\\
\end{aligned}
注意到其值恒大于 0 ，则原函数是个单峰函数。
（图源网络，侵删）
若有无穷多个费马点当且仅当 $n$ 为偶数且 $n$ 个点全部共线。其余情况，当且仅当一阶偏导数为 $0$ 时取得最低点。
## 不动点迭代法
大部分人都给出了模拟退火的做法，但实际上由于是单峰函数，爬山做法也可以保证正确性。另外也有提出基于单峰函数的三分套三分的做法。那么这里给出一种迭代的方法。
### 不动点与吸性不动点
有一类定义在实数域上的函数 $\varphi$。定义 $\varphi(x)=x$ 的解叫做函数 $\varphi$ 的不动点。
那么有一类特殊的不动点，叫做吸性不动点。我们记该不动点为 $x^* $，那么对于一个区间 $[a,b]$ 中的一个任何 $x$ 值，迭代序列 $\{x,\varphi(x),\varphi(\varphi(x)),\varphi(\varphi(\varphi(x))),\dots\}$ 收敛于 $x^* $。
举个例子，我们要求 $\cos(x)$ 的不动点，该不动点是一个在 $\mathbb R$ 内的吸性不动点。也就是说从任何一个 $x$ 值开始迭代，都会收敛至该不动点。
我们不妨从 $0$ 或者 $1145141919810$ 开始。那么其迭代序列为：
可以看出其值不断向不动点 $x^* = 0.7390851$ 收敛。
在图像上的体现如下（以$0$ 为例，$y=\cos x$ 为绿线，$y=x$ 为橙线，黑线为迭代过程）：
然而并非所有的函数都有吸性不动点，如 $f(x)=2x^3-1$。
我们分别从 $0.99$ 和 $1.01$ 开始迭代（$1$ 为其不动点），迭代序列如下：
注意到,虽然初始值离不动点很近，但是并没有向不动点收敛，反而发散。
在图像上的体现如下（以$1.01$ 为例，$y=2x^3-1$ 为绿线，$y=x$ 为橙线，黑线为迭代过程）：
那么问题来了：什么时候一个函数的多次迭代会向它的不动点收敛，什么时候发散呢？
### 不动点迭代法的收敛定理
注意到上图，$y=\cos x$ 在不动点附近较为平缓，而 $y=2x^3-1$ 则在不动点附近较为陡峭。那么可以猜测，迭代能否收敛是否与函数在不动点附近的斜率有关？
先给出不动点迭代的收敛定理：
设迭代函数 $\varphi(x)$ 在 $[a,b]$ 上连续且可导。满足以下两个条件：
1. 对于任意的 $x\in[a,b]$ ，有 $\varphi(x)\in[a,b]$。
2. 存在一个整数 $L$ ，满足 $0<L<1$， 且对于任意的 $x\in[a,b]$ ，均满足 $|\varphi^{\prime}(x)|\leq L$。
那么这样的迭代函数有以下性质：
1. 方程 $\varphi(x)=x$ 在 $[a,b]$ 内有唯一解 $x^* $（即不动点）。
2. 对于任意初值 $x_0\in[a,b]$ ，迭代法 $x_{k+1}=\varphi(x_k)$ 均收敛于 $x^* $。
3. $|x_k-x^* |\leq \dfrac{L}{1-L}|x_k-x_{k-1}|
-
|x_k-x^* |\leq \dfrac{L^k}{1-L}|x_1-x_0|
对第一条性质的证明
设 f(x)=x-\varphi(x)，则 f(x) 在 [a,b] 上连续可导。由条件1可知 a\leq\varphi(a),\varphi(b)\leq b，那么 f(a)=a-\varphi(a)\leq 0，同理也有 f(b)=b-\varphi(b)\geq 0。由根的存在性定理，方程 f(x)=0 在 [a,b] 上至少有一根。
接下来由于 |\varphi^{\prime}(x)|\leq L<1，那么 f^\prime(x)=1-\varphi^{\prime}(x)>0。则 f(x) 在 [a,b] 上单增。
综上所述，f(x) 在 [a,b] 上有且仅有一根。
证毕。
对第二条性质的证明
对于函数 \varphi(x_k)=x_{k+1}，由拉格朗日中值定理可知，存在一个实数 \xi\in(x_{k},x^* )（我们这里假设 x_{k}<x^* ，实际上若不满足两者交换即可）使得满足
\begin{aligned}
x_{k+1}-x^* =\varphi(x_k)-\varphi(x^* )=\varphi^\prime(\xi)(x_k-x^* )
\end{aligned}
那么又由于 |\varphi(x)|\leq L ，故
\begin{aligned}
|x_{k+1}-x^* |=|\varphi^\prime(\xi)(x_k-x^* )|\leq L|x_k-x^* |
\end{aligned}
接下来递推可得
\begin{aligned}
|x_{k}-x^* |&\leq L|x_{k-1}-x^* |\\
&\leq L^2|x_{k-2}-x^* |\\
&\dots\\
&\leq L^{k}|x_0-x^* |
\end{aligned}
而 L<1 ，那么 \lim_{k\rightarrow\infty}(|x_k-x^* |)=0。
证毕。
对第三、四条性质的证明
同样可以使用拉格朗日中值定理，存在一个实数 \xi\in(x_{k-1},x_{k} )（假设 x_{k-1}<x_k ）使得满足
\begin{aligned}
x_{k+1}-x_k =\varphi(x_k)-\varphi(x_{k-1} )=\varphi^\prime(\xi)(x_k-x_{k-1} )
\end{aligned}
由 |\varphi(x)|\leq L ，可得
\begin{aligned}
|x_{k+1}-x_k |=|\varphi^\prime(\xi)(x_k-x_{k-1} )|\leq L|x_k-x_{k-1} |
\end{aligned}
由性质二，
\begin{aligned}
|x_{k+1}-x^* |&\leq L|x_k-x^* |\\
&=L|x_{k+1}-x^* -(x_{k+1}-x_k)|\\
&\leq L|x_{k+1}-x^* |+L|x_{k+1}-x_k|
\end{aligned}
移项可以得到
\begin{aligned}
|x_{k+1}-x^* |\leq\dfrac{L}{1-L}|x_{k+1}-x_{k}|
\end{aligned}
第三条性质证毕。
接下来递推可得
\begin{aligned}
|x_{k}-x^* |&\leq \dfrac{L}{1-L}|x_k-x_{k-1} |\\
&\leq \dfrac{L^2}{1-L} |x_{k-1}-x_{k-2} |\\
&\dots\\
&\leq \dfrac{L^k}{1-L}|x_1-x_0 |
\end{aligned}
第四条性质证毕。
算法流程（求 f(x) 的不动点）
- 选取一个可以收敛至不动点的区间（满足 \varphi^\prime(x)\leq L<1）。
- 在选好的区间内任意选取一个初始值 x_0。
- 将其带入函数并计算 x_k=f(x_{k-1})。
- 重复第三步直到接近不动点。
那什么时候算是接近了不动点呢？
由第三条性质我们可以知道 |x_{k}-x^* |\leq{\large\frac{L}{1-L}}|x_{k}-x_{k-1}| ，那么我们可以给定一个误差极限 \varepsilon，当 |x_{k}-x_{k-1}|<\varepsilon 时，|x_{k}-x^* |<{\large\frac{L}{1-L}}\varepsilon ，此时我们可以认为 x_k\approx x^* ，将 x_k 作为方程的解即可。
算法效率（收敛速度）
设 e_k=|x_k-x^* |，若存在 p\geq 1 和 c>0 满足
\lim_{k\rightarrow\infty}\dfrac{e_{k+1}}{e_k^p}=c
那么我们称该迭代法为 p 阶收敛，显然 p 越大，收敛速度越快。
将 \varphi(x) 在 x^* 处泰勒展开，有：
\begin{aligned}
\varphi(x)=\varphi(x^* )+\varphi^\prime(x^* )(x-x^* )+\dfrac{\varphi^{\prime\prime}(x^* )}{2!}(x-x^* )^2+\dots
\end{aligned}
若 \varphi^\prime(x^* )=\varphi^{\prime\prime}(x^* )=\dots=\varphi^{(p-1)}(x^* )=0，而 \varphi^{(p)}(x^* )\neq0。
那么有
\begin{aligned}
\varphi(x)=\varphi(x^* )+\dfrac{\varphi^{(p)}(x^* )}{p!}(x-x^* )^p+\dfrac{\varphi^{(p+1)}(x^* )}{(p+1)!}(x-x^* )^{p+1}+\dots
\end{aligned}
则
\begin{aligned}
x_{k+1}-x^* &=
\varphi(x_k)-\varphi(x^* )\\&
=\dfrac{\varphi^{(p)}(x^* )}{p!}(x_k-x^* )^p+\dfrac{\varphi^{(p+1)}(x^* )}{(p+1)!}(x_k-x^* )^{p+1}+\dots
\end{aligned}
那么
\begin{aligned}
\dfrac{e_{k+1}}{e_k^p}=\dfrac{|x_{k+1}-x^* |}{|x_{k}-x^* |^p}
=\left|\dfrac{\varphi^{(p)}(x^* )}{p!}+\dfrac{\varphi^{(p+1)}(x^* )}{(p+1)!}(x_k-x^* )+\dots\right|
\end{aligned}
当 k\rightarrow\infty 时，|x_k-x^* | \rightarrow 0 ，此时 \dfrac{e_{k+1}}{e_k^p}\rightarrow\left|\dfrac{\varphi^{(p)}(x^* )}{p!}\right| 。
那么迭代法 x_{k+1}=\varphi(x_k) 的收敛阶为 p。
综上所述，若函数 \varphi(x) 在 x^* 附近满足 p 阶可导，且 \varphi^\prime(x^* )=\varphi^{\prime\prime}(x^* )=\dots=\varphi^{(p-1)}(x^* )=0，而 \varphi^{(p)}(x^* )\neq0，则迭代法 x_{k+1}=\varphi(x_k) 的收敛阶为 p。
应用
解方程，将其形式转化成如 f(x)=x 的形式。
由于不动点迭代法只能解决函数单调的部分，可以将一个函数划为若干段单调区间分别求每段的解。
值得注意的是，同一个方程不同的格式有不同的结果。
举个例子，求 2x^3-1=x 的解。
若构造函数 f(x)=2x^3-1 并求它的不动点。明显是不行的。用定理来解释，f^\prime(x)=6x^2-1 ，不动点 1 附近的斜率为 6\times1^2-1=5 ，是无法满足收敛的。上文的图也有提到。
若构造函数 f(x)=\sqrt[3]{\dfrac {x+1}{2}} ，求其不动点，此时 f^\prime(x)=\dfrac{\sqrt[3]{4}}{6}(x+1)^{-\frac23} ，此时 1 附近的斜率为 \dfrac16，是可以收敛的。
回到这道题，显然我们要求的是其一阶偏导数为 0 时的解，我们可以构造一个函数 g:\mathbb {R}^2\rightarrow \mathbb{R}^2 ，即由一个点到另外一个点的一一对应关系。
要满足如下式子
\begin{aligned}
f^\prime_x(x,y)=\sum_{i=1}^n\dfrac {x-x_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}=0\\
f^\prime_y(x,y)=\sum_{i=1}^n\dfrac {y-y_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}=0\\
\end{aligned}
那么将 x,y 单独提出来，
\begin{aligned}
x=\dfrac{\sum_{i=1}^n\dfrac {x_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}{\sum_{i=1}^n\dfrac {1}{\sqrt{(x-x_i)^2+(y-y_i)^2}}},y=\dfrac{\sum_{i=1}^n\dfrac {y_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}{\sum_{i=1}^n\dfrac {1}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}
\end{aligned}
那么 g 函数则可以表示为：
\begin{aligned}
g(x,y)=\left(\dfrac{\sum_{i=1}^n\dfrac {x_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}{\sum_{i=1}^n\dfrac {1}{\sqrt{(x-x_i)^2+(y-y_i)^2}}},\dfrac{\sum_{i=1}^n\dfrac {y_i}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}{\sum_{i=1}^n\dfrac {1}{\sqrt{(x-x_i)^2+(y-y_i)^2}}}\right)
\end{aligned}
我们要求的，就是方程 g(x,y)=(x,y) 的解。
注意到所有的坐标 x_i,y_i 都大于 0 ，那么这个函数的曲线也就是从 y=x 的左边穿到了右边，又由于二阶偏导数恒为正，最多只有一个解，那么与 y=x 相交时的曲线斜率就不会大于 1。那么这样的二元方程同样可以使用迭代法逼近其不动点。代码如下：
#include<cstdio>
#include<iostream>
#include<cmath>
using namespace std;
const int N=105;
int n,x[N],y[N];
double dis(double X1,double X2,double Y1,double Y2){
return sqrt((X1-X2)*(X1-X2)+(Y1-Y2)*(Y1-Y2));
}
int main(){
int T;
scanf("%d",&T);
while(T--){
scanf("%d",&n);
double ansx=0,ansy=0;
for(int i=1;i<=n;i++)
scanf("%d%d",x+i,y+i),ansx+=x[i],ansy+=y[i];
ansx/=n;ansy/=n;
while(1){
double xz=0,yz=0,fm=0;
for(int i=1;i<=n;i++){
double now=dis(ansx,x[i],ansy,y[i]);
xz+=x[i]/now;yz+=y[i]/now;fm+=1.0/now;
}
double nowansx=xz/fm,nowansy=yz/fm;
if(fabs(nowansx-ansx)<1e-8&&fabs(nowansy-ansy)<1e-8)break;
ansx=nowansx;ansy=nowansy;
}
double sum=0;
for(int i=1;i<=n;i++)sum+=dis(ansx,x[i],ansy,y[i]);
printf("%.0lf\n",sum);
T&&puts("");
}
return 0;
}