矩阵压缩储存
矩阵是许多科学与工程计算问题中研究的数学对象。我们感兴趣的不是矩阵本身而是如何储存矩阵的元,从而使矩阵的各种运算能够有效地进行。
矩阵压缩需要矩阵具有某些相同的元素或者零元素,假若值相同的元素或者零元素在矩阵中的分布具有一定规律,则我们称此类矩阵为特殊矩阵;反之,成为稀疏矩阵。下面分别讨论他们的压缩储存。
特殊矩阵
对称矩阵若n阶矩阵A中的元满足
则称为n阶对称矩阵。
显然我们可以折半储存,每一对只存一次,对于主对角线对称的矩阵我们储存其下三角(包括对角线)中的元。我们用
s
a
[
n
(
n
+
1
)
/
2
]
sa[n(n+1)/2]
sa[n(n+1)/2]作为n阶对称矩阵A的储存结构,则
s
a
[
k
]
,
a
i
j
sa[k],a_{ij}
sa[k],aij之间存在着一一对应的关系。
a 11 a_{11} a11 | a 21 a_{21} a21 | a 22 a_{22} a22 | a 31 a_{31} a31 | ⋯ \cdots ⋯ | a n , 1 a_{n,1} an,1 | ⋯ \cdots ⋯ | a n , n a_{n,n} an,n |
---|
对角线矩阵所有的非零元素都集中在以主对角线为中心的带状区域中。对于这种的矩阵我们也可以按照某一原则将其压缩到一维数组上。
a 11 a_{11} a11 | a 12 a_{12} a12 | a 21 a_{21} a21 | a 22 a_{22} a22 | a 23 a_{23} a23 | ⋯ \cdots ⋯ | a n , n − 1 a_{n,n-1} an,n−1 | a n , n a_{n,n} an,n |
---|
这些特殊矩阵中非零元素的分布都具有一定明显的规律,所以我们可以将其压缩到一个一维数组中,并找到对应关系。
稀疏矩阵
在实际应用中我们还经常遇到另一类矩阵,其非零元素较零元素少,且分布没有一定规律,我们称之为稀疏矩阵。
人们无法给出准确的定义,它只是一个凭人们的直觉来了解的概念。假设在
m
∗
n
m*n
m∗n的矩阵中,有t个元素不为零。令
δ
=
t
m
∗
n
\delta=\frac{t}{m*n}
δ=m∗nt,称
δ
\delta
δ为矩阵的稀疏因子。通常认为
δ
≤
0.05
\delta\le0.05
δ≤0.05时称为稀疏矩阵。
三元组顺序表我们可以以顺序储存结构来表示三元组表,则可得稀疏矩阵的一种压缩方式。
typedef struct {
int i, j;
int e;
}Triple;//三元组的数据结构
typedef struct {
Triple data[MAX + 1];//非零三元组表
int mu, nu, tu;//矩阵的行数、列数和非零元的个数
}TSMatrix;
注意:为了还原数组需要记录原数组的行,列数。非零元个数是为了方便运算。
在此的非零元三元组是以行序为主序顺序排列的,我们会再后面看出这样做将有利于进行某种矩阵运算。
转置运算
朴素算法:(1)将矩阵的行列值相互交换;
(2)将每个元组中的i和j相互调换;
(3)重排三元组之间的次序便可实现矩阵的转置。
代码略。
搜索算法:根据转置后的三元组顺序寻找三元组转置。换句话说,按照矩阵M的列序来进行转置。为了找到M的每一列中所有的非零元素,需要对原三元组表从第一行起整个扫描一遍,由于原三元表是以M为行序为主序来存放每个非零元的,所以由此得到的恰好是转置后的应有的顺序。
void TransposeSMAtrix(TSMatrix M, TSMatrix &T) {
T.mu = M.nu, T.nu = M.mu, T.tu = M.tu;
if (T.tu) {
int q = 1;
for(int col=1;col<=M.nu;col++)
for(int p=1;p<=M.tu;p++)
if (M.data[p].j == col) {
T.data[q].i = M.data[p].j;
T.data[q].j = M.data[p].i;
T.data[q].e = M.data[p].e;
q++;
}
}
return;
}
我们知道一般矩阵的直接求转置矩阵时间复杂度是
O
(
m
u
∗
n
u
)
O(mu*nu)
O(mu∗nu)。
上面的算法主要是col与p的两层循环,当矩阵的非零元个数tu与
m
u
∗
n
u
mu*nu
mu∗nu同数量级时,算法的时间复杂度就成了
O
(
m
u
∗
n
u
2
)
O(mu*nu^2)
O(mu∗nu2),所以上面算法只适用于
t
u
≪
m
u
∗
n
u
tu\ll mu*nu
tu≪mu∗nu的情况。
转置高效算法
按照原矩阵的三元组的次序进行转置,并将转置后的三元组置入b中恰当的位置。如果能预先确定矩阵M中的每一列的第一个非零元素在
b
.
d
a
t
a
b.data
b.data中应有的位置,那么在对
a
.
d
a
t
a
a.data
a.data中的三元组依次转置时,便可直接放到
b
.
d
a
t
a
b.data
b.data中恰当的位置上去。为了确认这些位置,在转置前,应先求得M的每一列中非零元的个数,进而求得每一列中非零元的个数,进而求得每一列的第一个非零元在
b
.
d
a
t
a
b.data
b.data中应有的位置。
在此,需要附设num和cpot两个数组。
n
u
m
[
c
o
l
]
num[col]
num[col]表示矩阵M中第col列中非零元的个数,
c
p
o
t
[
c
o
l
]
cpot[col]
cpot[col]指示M中第col列的一个非零元在
b
.
d
a
t
a
b.data
b.data中的恰当位置。有
void FastTransposeSMatrix(TSMatrix M, TSMatrix&T) {
T.mu = M.nu, T.nu = M.mu, T.tu = M.tu;
if (T.tu) {
for (int col = 1; col <= M.nu; col++)num[col] = 0;
for (int t = 1; t <= M.tu; t++)num[M.data[t].j]++;
cpot[1] = 1;
for (int col = 2; col <= M.nu; col++)cpot[col] = cpot[col - 1] + num[col - 1];
for (int p = 1; p <= M.tu; p++) {
int col = M.data[p].j;
int q = cpot[col];
T.data[q].i = M.data[p].j;
T.data[q].j = M.data[p].j;
T.data[q].e = M.data[p].e;
cpot[col]++;
}
}
}
这个算法仅比前一个算法多用了两个辅助数组,从时间复杂度上看算法中有四个单循环,循环次数分别是nu,tu因而总的时间复杂度是
O
(
n
u
+
t
u
)
O(nu+tu)
O(nu+tu)在非零元素tu与
n
u
∗
m
u
nu*mu
nu∗mu同一数量级时,其时间复杂度为
O
(
m
u
∗
n
u
)
O(mu*nu)
O(mu∗nu),和经典算法相同。
矩阵乘法
经典算法是使用二维数组实现矩阵乘法。
for(i=1;i<=m1;++i)
for(j=1;j<=n2;++j){
Q[i][j]=0;
for(k=1;k<=n1;++k)Q[i][j]+=M[i][k]*N[k][j];
}
此算法的时间复杂度是
O
(
m
1
∗
n
1
∗
n
2
)
O(m_1*n_1*n_2)
O(m1∗n1∗n2)。
那如果使用顺序三元组如何实现矩阵乘法呢?
我们想要得到的数组也用顺序三元组表示,这就要利用原本三元组的顺序,每次取M数组的j,找到对应的N数组的i,相乘累加到暂存数组中。待一行的乘积和都完成后按顺序存储。最终就得到了顺序三元组。
void MultSMatrix(RLSMatrix M, RLSMatrix N, RLSMatrix &Q) {
int ctemp[MAX];
int tp, t, ccol;
if (M.nu != N.mu) { cout << "错误" << endl; return; }
Q.nu = M.nu; Q.mu = N.mu; Q.tu = 0;
if (M.tu*N.tu != 0) {
for (int arow = 1; arow <= M.mu; arow++) {
memset(ctemp, 0, sizeof(ctemp));
Q.rpos[arow] = Q.tu + 1;
if (arow < M.mu)tp = M.rpos[arow + 1];
else tp = M.tu + 1;
for (int p = M.rpos[arow]; p < tp; p++) {
int brow = M.data[p].j;
if (brow < N.mu)t = N.rpos[brow + 1];
else t = N.tu + 1;
for (int q = N.rpos[brow]; q < t; q++) {
ccol = N.data[q].j;
ctemp[ccol] += M.data[p].e*N.data[q].e;
}
}
for (ccol = 1; ccol <= Q.nu; ccol++) {
if (++Q.tu > MAX) { cout<<"错误"<<endl; return; }
Q.data[Q.tu] = { arow, ccol, ctemp[ccol] };
}
}
}
}
M数组m行n列,N数组是n行p列,则算法总复杂度是 O ( m ∗ p ) O(m*p) O(m∗p),这显然非常高效。
十字链表
当矩阵的非零元个数和的个数和位置在操作过程中变化太大时,就不宜采用顺序储存结构来表示三元组的线性表。例如将一个矩阵加到另一个矩阵上的操作时,由于非零元的删除或插入将会引起三元组表中的元素的移动。为此,对于这种类型的矩阵,我们可以使用链表储存。
如果只需要列遍历或行遍历可以使用邻接链表储存,但如果既要列遍历又要行遍历就需要使用十字链表,储存两个索引分别指向行与列。
void CreateSMatrix_OL(CrossList &M){
int m, n, t;
cin >> m >> n >> t;
M.mu = m; M.nu = n; M.tu = t;
M.rhead = new OLink;
M.chead = new OLink;
M.rhead = M.chead = NULL;
int i, j, e;
for (cin >> i >> j >> e; i != 0; cin >> i >> j >> e) {
OLink p = new OLNode;
OLink q = new OLNode;
p->i = i; p->j = j; p->e = e;
if (M.rhead[i] == NULL || M.rhead[i]->j > j) {
p->right = M.rhead[i]; M.rhead[i] = p;
}
else {
for (q = M.rhead[i]; (q->right) && q->right->j < j; q = q->right);
p->right = q->right; q->right = p;
}
if (M.chead[j] == NULL || M.chead[j]->i > i) {
p->down = M.chead[j]; M.chead[j] = p;
}
else {
for (q = M.chead[j]; (q->down) && q->down->i < i; q = q->down);
p->down = q->down; q->down = p;
}
}
}
对于m行n列且有t个非零元的稀疏矩阵,算法的时间复杂度是
O
(
t
∗
s
)
,
s
=
m
a
x
(
m
,
n
)
O(t*s),s=max(m,n)
O(t∗s),s=max(m,n),这是因为每建立一个非零元的结点时都要寻找他在行表和列表中的插入位置,此算法对于非零元输入的先后顺序没任何要求。反之,若按照行序为主序的次序依次输入三元组,则可将建立十字链表的算法改写成
O
(
t
)
O(t)
O(t)的数量级的(t为非零元的个数)。
将一个矩阵加到另一个矩阵与构建矩阵的算法类似。