MATLAB模板

这些模板即开即用，或作为通用化的起点， 例如处理多变量且带有正交约束的问题，或是代码优化。 在最简单的模式下，用户只需提供一个函数来最小化F(Y)（以及可选的第一和第二导数）在F.m（在dF.m和ddF.m中）和一个初始猜测Y0。调用序列则是一个单一的调用sg_min（以Stiefel和Grassmann命名）：

  [fopt, yopt] =  sg_min(Y0).

例如，如果函数F.m的形式为

  function f=F(Y)
  f=trace( Y' * diag(1:10) * Y * diag(1:3) );

我们可以调用sg_min(rand(10,3))，它指定了一个随机的起始点。

我们强烈建议也提供第一导数信息：

  function df=dF(Y)
  df =  2 * diag(1:10) * Y * diag(1:3);

代码可以通过有限差分进行计算，但这非常慢且存在问题。用户也可以提供第二导数信息（这对于速度的提升远不如提供第一导数信息重要，但可能提高精度）：

  function ddf=ddF(Y,H)
  ddf =  2 * diag(1:10) * H * diag(1:3);

一个示例测试调用，其中F(Y)已知其最优值为10：

>> rand('state',0);	% 初始化随机数生成器
>> fmin = sg_min(rand(10,3))
iter    grad            F(Y)              flops         step type
0       1.877773e+01    3.132748e+01         2470       none
  invdgrad: Hessian非正定，CG提前终止
1       1.342879e+01    2.011465e+01       163639       Newton step
  invdgrad: Hessian非正定，CG提前终止
2       1.130915e+01    1.368044e+01       344198       Newton step
  invdgrad: Hessian非正定，CG提前终止
3       5.974819e+00    1.063045e+01       515919       Newton step
  invdgrad: 通过CG求逆Hessian达到最大迭代次数
4       1.135417e+00    1.006835e+01       789520       Newton step
5       5.526359e-02    1.000009e+01      1049594       Newton step
6       5.072118e-05    1.000000e+01      1337540       Newton step
7       1.276025e-06    1.000000e+01      1762455       Newton step

fmin =
   10.0000

sg_min的完整调用序列是

  [fopt, yopt]=sg_min(Y0,rc,mode,metric,motion,verbose,...
                      gradtol,ftol,partition),

其中Y0是必需的，可选参数包括（见表9.1）：

表9.1： sg_min可选参数的简短列表

参数	描述
rc	{'real', 'complex' }
mode	{'newton', 'dog', 'prcg', 'frcg' }
metric	{'flat', 'euclidean', 'canonical' }
motion	{'approximate', 'exact' }
verbose	{'verbose', 'quiet'}
ftol	首次收敛容差
gradtol	第二次收敛容差
partition	表示f对称性的p分区

rc

指定矩阵是否具有复数项。尽管大部分代码对此不敏感，但rc对于正确计算问题的维度至关重要。省略时，sg_min会根据Y0是否具有非零虚部进行猜测。

mode

选择将使用的优化方法。'newton'选择带有共轭梯度法的海森堡矩阵逆牛顿法。'dog'选择插值最速下降和牛顿法步长的狗腿算法。'frcg'选择Fletcher-Reeves共轭梯度法。最后，'prcg'选择Polak-Ribiere共轭梯度法，其优点是不需要非常精确的海森堡矩阵（因此，如果使用有限差分近似实现ddF.m，这是这些方法中最安全的一种）。虽然'newton'是默认选项，但这并不意味着我们推荐它优于其他方法，其他方法在某些问题上可能更准确且成本更低。

metric

选择赋予约束表面的几何类型。这最终会影响协变海森堡矩阵的定义。'flat'将无约束海森堡矩阵的结果投影到切空间上，而'euclidean'和'canonical'则添加了与其几何相关的连接项。'euclidean'是默认选项。

motion

选择沿流形移动时是使用所选度量的测地线运动方程的解析解，还是使用计算成本较低的近似解（默认）。（对于'flat'度量，没有测地线方程，因此在这种情况下此参数无效。）

verbose

确定函数在执行时是否显示每次迭代的报告。当此参数为'verbose'（默认）时，将显示数据并记录在全局SGdata中。当此参数为'quiet'时，则不显示或记录收敛数据。

gradtol和ftol

如果满足以下两个条件之一，我们声明收敛：
grad/gradinit < gradtol（默认1e-7）或(f-fold)/f < ftol（默认1e-10）， 其中gradinit是梯度的初始大小，fold是上一次迭代中F(Y)的值。

partition

是一个单元数组，其元素是表示1:p不相交分区的索引向量。 如果F没有对称性，则分区为num2cell(1:p)。 如果F(Y)=F(YQ)对所有正交Q成立，则分区为{1:p}。 如果F的对称性为F(Y)=F(YQ)，其中正交Q具有稀疏结构

\begin{bmatrix} \times & \times & \\ \times & \times & \\ & & \times & \times & \times \\ & & \times & \times & \times \\ & & \times & \times & \times \end{bmatrix}

用户同样可以指定{3:5,1:2}或{[5 3 4],[2 1]}作为分区；即，分区是一组集合。 分区的目的是消除F(Y)的任何块对角正交对称性导致的零空间。 如果省略此参数，我们的代码将通过确定F的对称性（使用其partition.m函数）来选择分区。 然而，用户应注意，错误分区可能会破坏收敛。