CS205A HW6

距离上次更新已经 2133 天了，文章内容可能已经过时。

课程主页：https://graphics.stanford.edu/courses/cs205a-13-fall/schedule.html

这次回顾作业6。

Problem 1

（备注，这题没有讲明，但是从后面的讨论中可以推出这里为无向图）

假设$|V_0|=m$，并且

$${\overrightarrow{v}}_i={\overrightarrow{v}}_i^0,i=n-m+1,\dots ,n$$

记$B$为邻接矩阵，即

$$B_{ij}=\{\begin{array}{ll}1& (i,j)\in E\\ 0& \mathrm{其}\mathrm{他}\end{array}$$

将其记为如下分块形式：

$$B=\left[\begin{array}{cc}{B}_1& B_2\\ B_3& B_4\end{array}\right]\in {\mathbb{R}}^{n\times n}$$

其中

$$B_1\in {\mathbb{R}}^{(n-m)\times (n-m)},B_2\in {\mathbb{R}}^{(n-m)\times m},B_3\in {\mathbb{R}}^{m\times (n-m)},B_4\in {\mathbb{R}}^{m\times m}$$

记$1_k\in {\mathbb{R}}^k$表示全$1$的$k$维列向量，假设

$${\overrightarrow{v}}_i=[x_i,y_i{]}^T$$

那么记

$$\begin{array}{rl}\overrightarrow{x}& =\left[\begin{array}{c}{x}_1\\ ⋮\\ x_{n-m}\end{array}\right],\overrightarrow{y}=\left[\begin{array}{c}{y}_1\\ ⋮\\ y_{n-m}\end{array}\right]\\ {\overrightarrow{x}}^{\prime }& =\left[\begin{array}{c}{x}_{n-m+1}\\ ⋮\\ x_n\end{array}\right],{\overrightarrow{y}}^{\prime }=\left[\begin{array}{c}{y}_{n-m+1}\\ ⋮\\ y_n\end{array}\right]\\ V_1& =\left[\begin{array}{cc}\overrightarrow{x}& \overrightarrow{y}\end{array}\right],V_2=\left[\begin{array}{cc}{\overrightarrow{x}}^{\prime }& {\overrightarrow{y}}^{\prime }\end{array}\right]\end{array}$$

(a)对能量式进行化简

$$\begin{array}{rl}E({\overrightarrow{v}}_1,\dots ,{\overrightarrow{v}}_n)& =\sum _{(i,j)\in E}{\Vert {\overrightarrow{v}}_i-{\overrightarrow{v}}_j\Vert }_2^2\\ & =\sum _{(i,j)\in E}({\overrightarrow{v}}_i^T{\overrightarrow{v}}_i-2{\overrightarrow{v}}_i^T{\overrightarrow{v}}_j+{\overrightarrow{v}}_j^T{\overrightarrow{v}}_j)\end{array}$$

$\mathrm{\forall }1\le k\le n-m$，我们有

$$\begin{gathered} \begin{array}{rl}{\mathrm{\nabla }}_{{\overrightarrow{v}}_k}\mathrm{\Lambda }& =\sum _{(k,j)\in E}(2{\overrightarrow{v}}_k-2{\overrightarrow{v}}_j)+\sum _{(i,k)\in E}(-2{\overrightarrow{v}}_i+2{\overrightarrow{v}}_k)\\ & =2(\sum _{(k,j)\in E}{\overrightarrow{v}}_k+\sum _{(i,k)\in E}{\overrightarrow{v}}_k-\sum _{(k,j)\in E \\ 1\le j\le n-m}{\overrightarrow{v}}_j-\sum _{(k,j)\in E \\ n-m+1\le j\le n}{\overrightarrow{v}}_j-\sum _{(i,k)\in E \\ 1\le i\le n-m}{\overrightarrow{v}}_i-\sum _{(i,k)\in E \\ n-m+1\le i\le n}{\overrightarrow{v}}_i)\end{array} \end{gathered}$$

令上式为$\overrightarrow{0}$得到

$$\begin{gathered} \sum _{(k,j)\in E}{\overrightarrow{v}}_k+\sum _{(i,k)\in E}{\overrightarrow{v}}_k-\sum _{(k,j)\in E \\ 1\le j\le n-m}{\overrightarrow{v}}_j-\sum _{(i,k)\in E \\ 1\le i\le n-m}{\overrightarrow{v}}_i-\sum _{(k,j)\in E \\ n-m+1\le j\le n}{\overrightarrow{v}}_j-\sum _{(i,k)\in E \\ n-m+1\le i\le n}{\overrightarrow{v}}_i=\overrightarrow{0} \end{gathered}$$

写成矩阵形式为

$$\begin{array}{rl}(\left[\begin{array}{cc}{B}_1& B_2\end{array}\right]+{\left[\begin{array}{c}{B}_1\\ B_3\end{array}\right]}^T)1_n{V}_1-(B_1+B_1^T)V_1& =(B_2+B_3^T)V_2\\ (\text{diag}\left((\left[\begin{array}{cc}{B}_1& B_2\end{array}\right]+{\left[\begin{array}{c}{B}_1\\ B_3\end{array}\right]}^T)1_n\right)-(B_1+B_1^T))V_1& =(B_2+B_3^T)V_2\end{array}$$

记

$$\begin{array}{rl}A& =\text{diag}\left((\left[\begin{array}{cc}{B}_1& B_2\end{array}\right]+{\left[\begin{array}{c}{B}_1\\ B_3\end{array}\right]}^T)1_n\right)-(B_1+B_1^T)\\ {\overrightarrow{b}}_x& =(B_2+B_3^T){\overrightarrow{x}}^{\prime }\\ {\overrightarrow{b}}_y& =(B_2+B_3^T){\overrightarrow{y}}^{\prime }\end{array}$$

那么将上式写为分量形式得到

$$\begin{array}{rl}A\overrightarrow{x}& ={\overrightarrow{b}}_x\\ A\overrightarrow{y}& ={\overrightarrow{b}}_y\end{array}$$

最后验证$A$为对称正定矩阵，对称性显然，下证正定性，首先记

$$C_1=\left[\begin{array}{cc}{B}_1& B_2\end{array}\right],C_2=\left[\begin{array}{c}{B}_1\\ B_3\end{array}\right]$$

$\mathrm{\forall }\overrightarrow{y}=(y_1,\dots ,y_{n-m}{)}^T$，计算${\overrightarrow{y}}^{T}A\overrightarrow{y}$：

$$\begin{array}{rl}{\overrightarrow{y}}^{T}A\overrightarrow{y}& ={\overrightarrow{y}}^T(\text{diag}\left((C_1+C_2^T)1_n\right)-(B_1+B_1^T))\overrightarrow{y}\\ & ={\overrightarrow{y}}^T(C_1+C_2^T)1_n\overrightarrow{y}-{\overrightarrow{y}}^T(B_1+B_1^T)\overrightarrow{y}\\ & =\sum _{i=1}^{n-m}\sum _{j=1}^n(B_{ij}+B_{ji})y_i^2-\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})y_i{y}_j\end{array}$$

因为$B_{ij}\in \{0,1\}$，所以

$$\begin{array}{rl}\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})y_i{y}_j& \le \frac{1}{2}\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})(y_i^2+y_j^2)\\ & =\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})y_i^2\end{array}$$

因此

$$\begin{array}{rl}{\overrightarrow{y}}^{T}A\overrightarrow{y}& =\sum _{i=1}^{n-m}\sum _{j=1}^n(B_{ij}+B_{ji})y_i^2-\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})y_i{y}_j\\ & \ge \sum _{i=1}^{n-m}\sum _{j=1}^n(B_{ij}+B_{ji})y_i^2-\sum _{i=1}^{n-m}\sum _{j=1}^{n-m}(B_{ij}+B_{ji})y_i^2\\ & =\sum _{i=1}^{n-m}\sum _{j=n-m+1}^n(B_{ij}+B_{ji})y_i^2\\ & \ge 0\end{array}$$

当且仅当$y_i=0$时等号成立，所以$A$对称正定。

(b)(i)将之前讨论的部分实现即可，代码如下

B = zeros(totalVertices);
[m, n] = size(edges);
for i = 1:m
    x = edges(i, 1);
    y = edges(i, 2);
    B(x, y) = 1;
    B(y, x) = 1;
end

B1 = B(unconstrainedVertices, unconstrainedVertices);
B2 = B(unconstrainedVertices, constrainedVertices);
B3 = B2';
B4 = B(constrainedVertices, constrainedVertices);

A = sparse(diag(([B1, B2] + [B1; B3]') * ones(totalVertices, 1)) - B1 - B1');
rhs = (B2 + B3') * constraints;

(ii)算法如下：

$$\begin{array}{rl}{\overrightarrow{d}}_k& =\overrightarrow{b}-A{\overrightarrow{x}}_{k-1}\\ {\alpha }_k& =\frac{{\overrightarrow{d}}_k^T{\overrightarrow{d}}_k}{{\overrightarrow{d}}_k^{T}A{\overrightarrow{d}}_k}\\ {\overrightarrow{x}}_k& ={\overrightarrow{x}}_{k-1}+{\alpha }_k{\overrightarrow{d}}_k\end{array}$$

所以对应代码如下：

for i=1:maxIters
    % TODO:  Update curX
    d = rhs - A * curX;
    a1 = sum(d .* d, 1);
    a2 = diag(d' * A * d)';
    
    alpha = a1 / a2;
    curX = curX + alpha * d;
    
    
    % Display the current iterate
    curResult(unconstrainedVertices,:) = curX;
    plotGraph(curResult, edges, f);
    title(sprintf('Gradient descent iteration %d',i));
    drawnow;
    pause(.1);
end

(iii)算法如下：

$$\begin{array}{rl}\text{Update search direction: }& {\overrightarrow{v}}_k={\overrightarrow{r}}_{k-1}-\frac{{⟨{\overrightarrow{r}}_{k-1},{\overrightarrow{v}}_{k-1}⟩}_A}{{⟨{\overrightarrow{v}}_{k-1},{\overrightarrow{v}}_{k-1}⟩}_A}{\overrightarrow{v}}_{k-1}\\ \text{Line search: }& {\alpha }_k=\frac{{\overrightarrow{v}}_k^{\mathrm{\top }}{\overrightarrow{r}}_{k-1}}{{\overrightarrow{v}}_k^{\mathrm{\top }}A{\overrightarrow{v}}_k}\\ \text{Update estimate: }& {\overrightarrow{x}}_k={\overrightarrow{x}}_{k-1}+{\alpha }_k{\overrightarrow{v}}_k\\ \text{Update residual: }& {\overrightarrow{r}}_k={\overrightarrow{r}}_{k-1}-{\alpha }_{k}A{\overrightarrow{v}}_k\end{array}$$

代码如下：

%初始化
r = rhs - A * curX;
v = zeros(size(r)) + 1e-3;

for i=1:maxIters
    % TODO:  Update curX
    r1 = diag(r' * A * v)';
    v1 = diag(v' * A * v)';
    v = r - r1 ./ v1 .* v;
    alpha = diag(v' * r) ./ diag(v' * A * v);
    curX = curX + alpha' .* v;
    r = r - alpha' .* (A * v);

    % Display the current iterate
    curResult(unconstrainedVertices,:) = curX;
    plotGraph(curResult, edges, f);
    title(sprintf('Conjugate gradients iteration %d',i));
    drawnow;
    pause(.1);
end

(iv)共轭梯度法很快就收敛了。

Problem 2

(a)定义

$$\begin{array}{}\text{(1)}& {\overrightarrow{x}}^{\prime }& & =A^{-1}\overrightarrow{b}\text{(2)}& {\overrightarrow{x}}_k& & =M^{-1}(N{\overrightarrow{x}}_{k-1}+\overrightarrow{b})\text{(3)}& {\overrightarrow{e}}_k& & ={\overrightarrow{x}}_k-{\overrightarrow{x}}^{\prime }\end{array}$$

下面考虑${\overrightarrow{e}}_k$和${\overrightarrow{e}}_{k-1}$的递推关系：

$$\begin{array}{rl}{\overrightarrow{e}}_k& ={\overrightarrow{x}}_k-{\overrightarrow{x}}^{\prime }\\ & =M^{-1}(N{\overrightarrow{x}}_{k-1}+\overrightarrow{b})-{\overrightarrow{x}}^{\prime }\\ & =M^{-1}N({\overrightarrow{x}}_{k-1}+N^{-1}(\overrightarrow{b}-M{\overrightarrow{x}}^{\prime }))\end{array}$$

注意

$$A=M-N$$

所以由(1)可得

$$\begin{array}{rl}A{\overrightarrow{x}}^{\prime }& =(M-N){\overrightarrow{x}}^{\prime }=\overrightarrow{b}\\ \overrightarrow{b}-M{\overrightarrow{x}}^{\prime }& =-N{\overrightarrow{x}}^{\prime }\\ N^{-1}(\overrightarrow{b}-M{\overrightarrow{x}}^{\prime })& =-{\overrightarrow{x}}^{\prime }\end{array}$$

因此

$$\begin{array}{rl}{\overrightarrow{e}}_k& =M^{-1}N({\overrightarrow{x}}_{k-1}+N^{-1}(\overrightarrow{b}-M{\overrightarrow{x}}^{\prime }))\\ & =G({\overrightarrow{x}}_{k-1}-{\overrightarrow{x}}^{\prime })\\ & =G{\overrightarrow{e}}_{k-1}\end{array}$$

递推可得

$${\overrightarrow{e}}_k=G^k{\overrightarrow{e}}_0$$

假设$G$的特征值为${\lambda }_1,\dots ,{\lambda }_n$，对应的特征向量为${\overrightarrow{v}}_1,\dots ,{\overrightarrow{v}}_n$，记

$$V=\left(\begin{array}{c}{\overrightarrow{v}}_1,\dots ,{\overrightarrow{v}}_n\end{array}\right),\mathrm{\Lambda }=\text{diag}\{{\lambda }_1,\dots ,{\lambda }_n\}$$

题目假设${\overrightarrow{v}}_1,\dots ,{\overrightarrow{v}}_n$张成${\mathbb{R}}^n$，那么

$$\begin{array}{rl}G& =V\mathrm{\Lambda }{V}^{-1}\\ G^k& =V{\mathrm{\Lambda }}^k{V}^{-1}\end{array}$$

因为${\lambda }_i<1$，所以

$$\begin{array}{rl}{\mathrm{\Lambda }}^k& \to 0\\ G^k=V{\mathrm{\Lambda }}^k{V}^{-1}& \to 0\\ {\overrightarrow{e}}_k& \to 0\end{array}$$

因此

$${\overrightarrow{x}}_k\to {\overrightarrow{x}}^{\prime }$$

(b)利用圆盘定理即可：

圆盘定理

假设$A\in {\mathbb{R}}^{n\times n}$为$n$阶矩阵，令

$$R_i=\sum _{j\ne i}^n\left|a_{ij}\right|=\left|a_{i1}\right|+\cdots +\left|a_{i,i-1}\right|+\left|a_{i,i+1}\right|+\cdots +\left|a_{in}\right|$$

那么$A$的特征值$z$在如下圆盘中

$$|z-a_{ii}|⩽R_i,i=1,2,\cdots ,n$$

证明：任取$A$的特征值${\lambda }_0$，对应的特征向量为$\overrightarrow{x}$，那么

$$A\overrightarrow{x}={\lambda }_0\overrightarrow{x}$$

写成方程形式为

$$\{\begin{array}{l}{a}_{11}{x}_1+a_{12}{x}_2+\cdots +a_{1n}{x}_n={\lambda }_0{x}_1\\ a_{21}{x}_1+a_{2z}{x}_2+\cdots +a_{2n}{x}_n={\lambda }_0{x}_2\\ \dots \dots \dots \\ a_{n1}{x}_1+a_{n2}{x}_2+\cdots +a_{nn}{x}_n={\lambda }_0{x}_n\end{array}$$

设

$$|x_r|=max\{\left|x_1\right|,\left|x_2\right|,\cdots ,\left|x_n\right|\}$$

那么

$$({\lambda }_0-a_{rr})x_r=a_{r1}{x}_1+\cdots +a_{r,r-1}{x}_{r-1}+a_{r,r+1}{x}_{r+1}+\cdots +a_{rn}{x}_n$$

于是

$$\begin{array}{rl}|{\lambda }_0-a_{rr}|\left|x_r\right|& \le \left|a_{r1}\right|\left|x_{\mathrm{l}}\right|+\cdots +\left|a_{r,r-1}\right|\left|x_{r-1}\right|+\left|a_{r,r-1}\right|\left|x_{r+1}\right|+\cdots +\left|a_{rn}\right|\left|x_n\right|\\ & =(\left|a_{r1}\right|+\cdots +\left|a_{r,r-1}\right|+\left|a_{r,r\to 1}\right|+\cdots +\left|a_{rn}\right|)\left|x_r\right|\end{array}$$

即

$$|{\lambda }_0-a_{rr}|\left|x_r\right|\le R_r\left|x_r\right|$$

但是显然$|x_r|>0$，所以

$$|{\lambda }_0-a_{rr}|\le R_r$$

回到原题，我们有

$$\begin{array}{rl}{R}_i& =\sum _{j\ne 1}^n\left|g_{ij}\right|\\ & =\left|g_{i1}\right|+\cdots +\left|g_{i,i-1}\right|+\left|g_{i,i+1}\right|+\cdots +\left|g_{in}\right|\\ & =\frac{1}{|a_{ii}|}(\left|a_{i1}\right|+\cdots +\left|a_{i,i-1}\right|+\left|a_{i,i+1}\right|+\cdots +\left|a_{in}\right|)\\ & <1\end{array}$$

而$g_{ii}=0$，所以$G$的特征值满足

$$|\lambda |<R_i<1$$

因此收敛。

Problem 3

如果

$$\sum _{i=1}^n{a}_i{\overrightarrow{x}}_i=\overrightarrow{0}$$

左乘${\overrightarrow{x}}_k^{T}A,k=1,\dots ,n$可得

$$\begin{array}{}\text{(1)}& \sum _{i=1}^n{a}_i{\overrightarrow{x}}_k^{T}A{\overrightarrow{x}}_i=a_k{\overrightarrow{x}}_k^{T}A{\overrightarrow{x}}_k+\sum _{i\ne k}{a}_i{\overrightarrow{x}}_k^{T}A{\overrightarrow{x}}_i=0\end{array}$$

由$A$正交的定义可得

$${\overrightarrow{x}}_i^{T}A{\overrightarrow{x}}_j=0,i\ne j$$

所以(1)即为

$$a_k({\overrightarrow{x}}_k^{T}A{\overrightarrow{x}}_k)=0$$

如果$A$正定，${\overrightarrow{x}}_k$非零，那么

$$a_k=0,k=1,\dots ,n$$

此时${\overrightarrow{x}}_i$线性无关。

如果$A$半正定，那么因为${\overrightarrow{x}}_k^{T}A{\overrightarrow{x}}_k$可能为$0$，所以无法判断$a_k$，即此时无法推出${\overrightarrow{x}}_i$线性无关。