实例来自于《CCSDS标准的LDPC编译码仿真》中代码(实际上有点点差别),代码优化从以下几个方面进行
稀疏
仿真中的第一个困难在于ccsdscheckmatrix函数在输入SIZE_M很大的时候,先不说运行时间,直接就爆内存了。(输入参数4096,2/3)
先分析分析内存的问题,实际上这个函数的最后输出结果就是一个矩阵,这个矩阵的大小是12288×28672,计算double型的内存占用也就2G左右。但是函数运行过程中产生了很多中间变量没有清除。当然最后的解决办法也没有去管这些东西,由于矩阵H是稀疏矩阵,所以之际采用sparse后,这个运行就没有任何问题了。
对于矩阵H和H_sparse = spares(H),占用内存如下(当然H要是稀疏的,不然得不偿失)
Name Size Bytes Class Attributes
H 12288x28672 2818572288 double
H_sparse 12288x28672 1736712 double sparse
也可以对比稀疏矩阵和原始矩阵的运行时间(和稀疏程度有关)
代码:tic;H*message\';toc;
结果:时间已过 0.288934 秒。
代码:tic;H_sparse*message\';toc;
结果:时间已过 0.001210 秒。
类型转换
MATLAB中的运算符支持多种类型,譬如矩阵乘法中多用double型变量,但如果一个矩阵是逻辑输入也没有关系。但运算速度差异较大,譬如
>> Gc_logic = Gc>0;
>> a=randi([0 1],1,16384);
>> tic;b = a*Gc;toc
时间已过 0.107618 秒。
>> tic;b = a*Gc_logic;toc
时间已过 0.503132 秒。
观测结果类型为double,我们可以大胆推测实际上逻辑型变量在运算过程中先转化为了double型(逻辑怎么乘呢?)另一个实验结果是
>> tic;Gc_logic=double(Gc_logic);b = a*Gc_logic;toc
时间已过 0.546412 秒。
这一定程度上证明了我们的假设。所以在运算过程中数据类型是重要的,如果上述乘法出现在循环内,那么实现转化矩阵类型是必要的。即使只运行一次,那么显式的转化矩阵类型(特制新建变量)也有好处。譬如
>> tic;Gt=double(Gc_logic);b = a*Gt;toc
时间已过 0.373506 秒。
通过创建新变量,运行速度些许。
向量化
向量化实际上是原代码修改中获益最大的方法,这实际上是因为原先的译码程序写了太多的循环。向量化后运行时间变成了原先的1/40 。当然,原先的代码通用性强,而向量化这个过程实际上是运用了H的一些结构的。译码函数太复杂,此处不做举例。
此处分析差分调制中的例子(实际上对这个程序没有什么影响)
原来的代码是这个样子的(更新值为其本身和前一个值的异或)
encodeData_extend = [1 encodeData];
for num = 2:length(encodeData_extend)
encodeData_extend(num) = xor(encodeData_extend(num),encodeData_extend(num-1));
end
向量化的结果为(累加模二代替异或)
encodeData1 = [1 encodeData];
encodeData1_sum = cumsum(encodeData1);
encodeData_2 = mod(encodeData1_sum,2);
运行时间分别为
时间已过 0.023424 秒。
时间已过 0.015003 秒。
虽然后者没有快很多,但这取决于向量的长度,长度大的话会有较大差距。
其他
MATLAB中提及的都能对代码运行速度带来细微的改进,包括
-
将长脚本拆开成小段,调用执行;
-
将大的代码块分开为独立的函数;
-
将过分复杂的函数或是表达式采用简单的来代替;
-
采用函数,而不是脚本;
上述测试脚本(和以上运行条件有差别)
%% 稀疏矩阵测试
M=4096;
theta=[1 1 2 3 1 1 2 3 1 2 3 0 2 3 0 2 3 0 1 3 0 1 3 0 1 2];
fai=[1787 1077 1753 697 1523 5 2035 331 1920 130 4 85 551 15 1780 1960 3 145 1019 691 132 42 393 502 201 1064
1502 602 749 1662 1371 9 131 1884 1268 1784 19 1839 81 2031 76 336 529 74 68 186 905 1751 1516 1285 1597 1712
1887 521 590 1775 1738 2032 2047 85 1572 78 26 298 1177 1950 1806 128 1855 129 269 1614 1467 1533 925 1886 2046 1167
1291 301 1353 1405 997 2032 11 1995 623 73 1839 2003 2019 1841 167 1087 2032 388 1385 885 707 1272 7 1534 1965 588];
A = zeros(M);
B = eye(M);
L = 0:M-1;
for matrixNum = 1:14
t_k = theta(matrixNum);
f_4i_M = floor(4*L/M);
f_k = fai(f_4i_M+1,matrixNum)\';
col_1 = M/4*(mod((t_k+f_4i_M),4)) + ...
mod((f_k+L),M/4);
row_col = col_1+1 + L*M;
C_temp = zeros(M);
C_temp(ind2sub([M,M],row_col)) = 1;
C{matrixNum} = sparse(C_temp)\';
end
H = [A A A B B+C{1};B+C{8} B+C{7}+C{6} A A B;A B B+C{5} A C{4}+C{3}+C{2}];
H_23 = [A A;B C{11}+C{10}+C{9};C{14}+C{13}+C{12} B];
H=[H_23 H];
H_full = full(H);
whos H H_full
%% 稀疏矩阵乘法测试
message = randi([0 1],1,28672);
tic;H*message\';toc;
tic;H_full*message\';toc;
%% 数据类型测试
Gc = randn(16384);
Gc_logic = Gc>0;
a=randi([0 1],1,16384);
tic;b = a*Gc;toc
tic;b = a*Gc_logic;toc %逻辑型运行花费时间
tic;Gc_logic=double(Gc_logic);b = a*Gc_logic;toc %类型转换
tic;Gt=double(Gc_logic);b = a*Gt;toc %建立新变量
%% 向量化测试
encodeData = randi([0 1],1,1000000);
tic;
encodeData_extend = [1 encodeData];
for num = 2:length(encodeData_extend)
encodeData_extend(num) = xor(encodeData_extend(num),encodeData_extend(num-1));
end
toc;
tic;
encodeData1 = [1 encodeData];
encodeData1_sum = cumsum(encodeData1);
encodeData_2 = mod(encodeData1_sum,2);
toc;
View Code
|