Writing Fast Matlab code 6-7

6 内联简单函数

“内联一个函数”指用一个调用取代函数代码本身。注意你定义的M函数不要与MATLAB本身自带的函数混淆。如果你需要修改函数，在操作台键入：

edit [函数名]

以下函数值得内联：

conv:卷积

cross:叉乘

fft2:二维快速傅里叶变换

fliplr:将左边列交换到右边

flipud:将上边行交换到下边

ifft / ifft2 / ifftn

ind2sub:得到线性索引的下标

setdiff:集合差集

setxor:集合异或（不在交集中的元素）

unique:去除重复行列

7 引用操作

7.1 下标 vs. 索引

Matlab中存矩阵数据的时候是按一维数组重点内容存储的，故索引指的是元素在这个一维数组中的位置；而用下标能直接指到元素在矩阵中所在的位置。

向下扫描列比行更快！！

向下扫描列能提高缓存效率。

下标与索引间的转换可以用 sub2ind 函数 和 ind2sub 函数。

当然直接计算比调用M函数更快。下面给出转换规则：



7.2 向量化下标

利用下标比直接调用元素本身更好。例如，如果A是一个二维矩阵，那它的向量化下标：

A(rowv, colv)

此处rowv和colv都是向量，且长度任意。如果其中有一个是矩阵，会被Reshape成一个向量。在这里用行向量或者列向量都是没差别的。

7.3 索引向量

索引向量可以引用多个元素：

A(indexv)

7.4 通配符

利用通配符“ ：”可以指代整列或整行。

7.5 逻辑索引

一个常见用法是基于针对每个元素的条件，例如：

A(abs(A) < 1e-3) = 0

就将所有小于10.^-3的元素置为0

逻辑索引也可以用于选择非方阵的元素：

A(logical(eye(size(A))))

注意！这个比用diag(A)去引用了A对角元素快多了！！！

7.6 利用[]删除子矩阵

A(2,:) = [];

8 解线性方程Ax = b

一般来说x = A \ b，但若A奇异

x = pinv(A) * b（求近似解）

8.1 迭代方法

matlab中有很多迭代解决方法：

[bicg]双共轭梯度法求解线性方程组

共轭梯度算法是一种迭代方法，它从一组初值向量出发，定义代价函数，进而转化为最优化问题。

它是寻找一个向量，这个向量使得计算的代价函数为极值，这里为最小值。

clear;clc%利用是双共轭梯度算法求解线性方程组
A=[23.6 -6.8 4.5 6.3 7.2 8.6
4.2 7.62 0.8 -3.1 -4.5 5.3
9.8 -6.8 45.2 3.07 -4.37 7.82;
8.4 -9.05 24.4 2.65 6.53 -3.64;
-9.5 4.54 15.8 4.89 4.8 -7.48;
26.32 3.08 4.6 -3.5 -7.8 25.6];
b=[-9.50 8.12 42.30 -21.40 4.70 3.47]';
[x,flag,relres,iter,resverc]=bicg(A,b,1e-11)

[gmres]广义最小残差法求解线性方程组

当迭代值x所引起的相对残差小于或等于给定的收敛值时，停止迭代。

clear;clc%利用gmres函数求解线性方程组
A=[25.4 -6.5 5.5 6.3 6.2 4.6
4.5 7.62 0.7 -3.1 -4.5 5.1
8.9 -6.8 5.2 3.07 -4.37 7.8
8.4 -5.05 4.4 2.65 6.53 -4.4
-9.5 4.54 6.3 4.89 4.8 -7.4
6.3 4.08 1.6 -3.5 -7.8 5.6];
b=[-4.5 7.2 4.3 -2.4 24.7 5.7]';
[x,flag,relres,iter,resvec]=gmres(A,b)

还有很多很多的方法但语法都差不多，具体情况下用哪个参照下图：



除了可以将它存成矩阵形式，更为便捷的方法是将A表示为一个线性方程。比如



可表达为A * x = cumsum(x)\

另，有些特殊矩阵可以很快算出：

1.循环行列式矩阵；在傅里叶域中可对角化，x = ifft(fft(b) ./ fft(h))

2.三角阵和带状矩阵：可用LU分解：[L, U] = lu(sparse(A)); x = U\ (L\b);

3.泊松问题。可查找matlan multigrid