访问稀疏矩阵
非零元素
以下多条命令可以提供有关稀疏矩阵的非零元素的概要信息:
要尝试上述中的一些命令,请加载提供的稀疏矩阵 west0479
,该矩阵是 Harwell-Boeing 集合之一。
load west0479
whos
Name Size Bytes Class Attributes west0479 479x479 34032 double sparse
该矩阵为八个阶段的化工精馏塔建模。
尝试以下命令。
nnz(west0479)
ans = 1887
format short e west0479
west0479 = (25,1) 1.0000e+00 (31,1) -3.7648e-02 (87,1) -3.4424e-01 (26,2) 1.0000e+00 (31,2) -2.4523e-02 (88,2) -3.7371e-01 (27,3) 1.0000e+00 (31,3) -3.6613e-02 (89,3) -8.3694e-01 (28,4) 1.3000e+02 . . .
nonzeros(west0479)
ans = 1.0000e+00 -3.7648e-02 -3.4424e-01 1.0000e+00 -2.4523e-02 -3.7371e-01 1.0000e+00 -3.6613e-02 -8.3694e-01 1.3000e+02 . . .
注意
使用 Ctrl+C 随时停止列出 nonzeros
。
请注意,最初在默认情况下,nnz
与 nzmax
的值相同。即,非零元素数等于为非零元素分配的存储位置数。但是,如果将其他的数组元素置零,MATLAB® 不会动态释放内存。将某些矩阵元素的值更改为零时会更改 nnz
的值,但不会更改 nzmax
的值。
但是,可以根据需要将尽可能多的非零元素添加到矩阵中。不受 nzmax
原始值的限制。
索引和值
对于任何矩阵,无论是满矩阵还是稀疏矩阵,find
函数都会返回非零元素的索引和值。其语法是
[i,j,s] = find(S);
find
返回向量 i
中的非零值的行索引、向量 j
中的列索引以及向量 s
中的自身非零值。下面的示例使用 find
查找稀疏矩阵中的非零索引和值。sparse
函数同时使用 find
输出和矩阵大小重新创建矩阵。
S1 = west0479; [i,j,s] = find(S1); [m,n] = size(S1); S2 = sparse(i,j,s,m,n);
稀疏矩阵运算中的索引
由于稀疏矩阵是以压缩稀疏列格式存储的,因此为稀疏矩阵进行索引的相关成本与为满矩阵进行索引的相关成本不同。在只需更改稀疏矩阵中的若干元素时,这类成本可忽略不计,因此,在这类情况下,通常使用常规数组索引来重新分配值:
B = speye(4); [i,j,s] = find(B); [i,j,s]
ans = 1 1 1 2 2 1 3 3 1 4 4 1
B(3,1) = 42; [i,j,s] = find(B); [i,j,s]
ans = 1 1 1 3 1 42 2 2 1 3 3 1 4 4 1
42
位于 (3,1)
位置,MATLAB 会在非零值向量和下标向量中插入额外的一行,然后移动 (3,1)
后面的所有矩阵值。如果线性索引超过 2^48-1
(即当前矩阵中允许的元素数上限),使用线性索引在大型稀疏矩阵中访问或指定元素将失败。
S = spalloc(2^30,2^30,2); S(end) = 1
Maximum variable size allowed by the program is exceeded.
要访问其线性索引大于 intmax
的元素,请使用数组索引:
S(2^30,2^30) = 1
S = (1073741824,1073741824) 1
尽管在稀疏矩阵中进行索引以更改单个元素的成本可忽略不计,但该成本在循环环境下会增加,而且在大型矩阵中该操作可能会使执行速度变得很慢。因此,在需要更改大量稀疏矩阵元素的情况下,最好使用向量化方法而不要使用循环方法来执行该操作。例如,考虑稀疏单位矩阵:
n = 10000; A = 4*speye(n);
A
的元素慢于类似的向量化运算:tic A(1:n-1,n) = -1; A(n,1:n-1) = -1; toc
Elapsed time is 0.003344 seconds.
tic for k = 1:n-1 C(k,n) = -1; C(n,k) = -1; end toc
Elapsed time is 0.448069 seconds.
A
中的多个条目。如果为稀疏矩阵预分配内存,然后以类似的逐个元素的方式填充,会使对稀疏数组进行索引产生大量开销:
S1 = spalloc(1000,1000,100000); tic; for n = 1:100000 i = ceil(1000*rand(1,1)); j = ceil(1000*rand(1,1)); S1(i,j) = rand(1,1); end toc
Elapsed time is 2.577527 seconds.
构建索引和值向量则无需为稀疏数组进行索引,因此这种方法的速度快得多:
i = ceil(1000*rand(100000,1)); j = ceil(1000*rand(100000,1)); v = zeros(size(i)); for n = 1:100000 v(n) = rand(1,1); end tic; S2 = sparse(i,j,v,1000,1000); toc
Elapsed time is 0.017676 seconds.
因此,最好使用构造函数(例如 sparse
或 spdiags
函数)一次构造所有稀疏矩阵。
例如,假定需要稀疏形式的坐标矩阵 C
:
使用 sparse
函数,以及行下标、列下标和值组成的三联对组,直接构造该五列矩阵:
i = [1 5 2 5 3 5 4 5 1 2 3 4 5]'; j = [1 1 2 2 3 3 4 4 5 5 5 5 5]'; s = [4 1 4 1 4 1 4 1 -1 -1 -1 -1 4]'; C = sparse(i,j,s)
C = (1,1) 4 (5,1) 1 (2,2) 4 (5,2) 1 (3,3) 4 (5,3) 1 (4,4) 4 (5,4) 1 (1,5) -1 (2,5) -1 (3,5) -1 (4,5) -1 (5,5) 4
可视化稀疏矩阵
以图的形式查看非零元素在稀疏矩阵内的分布通常很有用。MATLAB spy
函数生成稀疏结构的模板视图,其中图表中的每点代表一个非零数组元素的位置。
例如:
加载提供的稀疏矩阵 west0479
,该矩阵是 Harwell-Boeing 集合之一。
load west0479
查看稀疏结构体。
spy(west0479)