---

22.2稀疏矩阵上的线性代数¶

Octave包括一个稀疏矩阵的多态解算器,其中用于分解矩阵的精确解算器取决于稀疏矩阵本身的性质。通常,相对于矩阵本身的因子分解成本,确定矩阵类型的成本很小,但在任何情况下,矩阵类型在计算后都会被缓存,因此每次在线性方程中使用时都不会重新确定。

如何求解线性方程的选择树是

1. 如果矩阵是对角的,直接求解并转到8
2. 如果矩阵是一个排列对角线,则直接考虑排列进行求解。转到8
3. 如果矩阵是正方形的,带状的,并且带状密度小于spparms(“bandden”)继续,否则转到4。
   1. 如果矩阵是三对角的,并且右侧不是稀疏连续的,则转到3b。
      1. 如果矩阵是埃尔米特矩阵,具有正实对角线,则尝试使用Cholesky因子分解LAPACK xPTSV。
      2. 如果上述失败或矩阵不是具有正实对角线的埃尔米特矩阵,则使用带枢轴的高斯消去法LAPACK xGTSV和goto 8。
   2. 如果矩阵是具有正实对角线的埃尔米特矩阵,则尝试使用LAPACK xBTRF。
   3. 如果以上失败或矩阵不是具有正实对角线的埃尔米特矩阵,则使用带枢轴的高斯消去法LAPACK xGBTRF和goto 8。
4. 如果矩阵是上三角或下三角,则执行稀疏前向或后向替换,并且goto 8
5. 如果矩阵是具有列排列的上三角矩阵或具有行排列的下三角矩阵,则执行稀疏前向或后向替换,并且goto 8
6. 如果矩阵是正方形的,则具有实正对角线的埃尔米特矩阵,尝试稀疏Cholesky因子分解使用CHOLMOD.
7. 如果稀疏Cholesky因子分解失败,或者矩阵不是具有实正对角线的Hermitian矩阵,并且矩阵是平方的,则使用UMFPACK.
8. 如果矩阵不是正方形的,或者之前的任何解算器都标记了奇异或接近奇异矩阵,则使用CXsparse¹⁰.

频带密度定义为频带中非零值的数量除以整个频带中值的总数。可以使用完全禁用带域矩阵解算器spparms设置类登至1(即。,spparms(“bandden”,1)).

QR解算器使用Dulmage Mendelsohn分解将问题分解为多个块,这些块可以作为过度确定的块、多个确定良好的块和一个最终确定的块来处理。对于具有强连通节点块的矩阵,这是一个巨大的胜利,因为LU分解可以用于许多块。它还显著提高了为过度确定的问题找到解决方案的机会,而不仅仅是返回向量NaNs

上面的所有解算器都可以计算条件数的估计值。这可以用于检测解中的数值稳定性问题,并强制使用最小范数解。然而,对于窄带、三角形或对角矩阵,计算条件数的成本是巨大的,事实上可能超过分解矩阵的成本。因此,在这些情况下不计算条件数,Octave依赖于简单的技术来检测奇异矩阵或底层矩阵LAPACK 带状矩阵情况下的代码。

用户可以使用强制矩阵的类型matrix_type作用这克服了查找矩阵类型的成本。然而,需要注意的是,错误地识别矩阵类型会导致不可预测的结果,因此matrix_type应小心使用。

 

:巢= normest (A)¶

:巢= normest (A,tol)¶

:[巢,iter] = normest (…)¶

矩阵的2-范数的估计A使用幂级数分析。

这通常用于大型矩阵,其中计算成本标准A)是禁止的,并且2-范数的近似是可接受的。

tol是计算2-范数的公差。默认情况下tol是1e-6。

可选输出iter返回所需的迭代次数normest以收敛。

详见: 标准1,标准,cond,condest.

 

:巢= 标准1 (A)¶

:巢= 标准1 (A,t)¶

:巢= 标准1 (A,t,x0)¶

:巢= 标准1 (Afcn,t,x0,p1,第2页, …)¶

:[巢,v] = 标准1 (A, …)¶

:[巢,v,w] = 标准1 (A, …)¶

:[巢,v,w,iter] = 标准1 (A, …)¶

估计矩阵的1-范数A使用块算法。

标准1对于只需要估计范数的大型稀疏矩阵是最好的。对于中小型矩阵,考虑使用标准A1.此外标准1可以用于估计线性算子的1-范数A当矩阵向量乘积A * x和A' *x可以计算。在这种情况下,而不是矩阵A,一个函数Afcn(旗帜,x)使用;它必须返回:

• 尺寸n属于A如果旗帜是暗淡的
• true ifA是一个真正的运算符,如果旗帜是真实的
• 结果A * x如果旗帜是“nottransp”
• 结果A' *x如果旗帜是“transp”

一个典型的案例是A从定义b^m,其中有结果A * x可以计算,甚至不需要明确的形式b^m通过

y=x;对于我1.m
  y=b * y;外循环

参数p1,第2页,…是的参数Afcn(旗帜,x,p1,第2页, …).

的默认值t是2。该算法需要具有大小的矩阵矩阵乘积nxn和nxt.

初始矩阵x0应该具有单位1-范数的列。默认初始矩阵x0具有第一列个(n1.n并且,如果t>1、具有随机元素的剩余列1.n,1.n,除以n.

在输出时,巢是期望的估计,v和w是这样的向量w=A * v具有标准w1.=c标准v1.. iter包含在iter1.迭代次数(最大值ishardcoded为5)和iter2.产品总数A * x或A' *x从算法执行。

算法注释:标准1在评估过程中使用随机数。因此,如果需要一致的结果状态的应在调用前修复标准1.

参考文献:N.J.Higham和F.Tisseur,矩阵1-模估计的块算法及其在1-模伪谱中的应用,SIAM J.矩阵分析。Appl。,第1185-1201页,第21卷,第4期,2000年。

详见: normest,标准,cond,condest.

 

:cest= condest (A)¶

:cest= condest (A,t)¶

:cest= condest (A,Ainvfcn)¶

:cest= condest (A,Ainvfcn,t)¶

:cest= condest (A,Ainvfcn,t,p1,第2页, …)¶

:cest= condest (Afcn,Ainvfcn)¶

:cest= condest (Afcn,Ainvfcn,t)¶

:cest= condest (Afcn,Ainvfcn,t,p1,第2页, …)¶

:[cest,v] = condest (…)¶

平方矩阵的1-范数条件数的估计A使用t测试向量和随机1-范数估计器。

可选输入t指定测试向量的数量(默认为5)。

输入可以是矩阵A(该算法特别适用于大型稀疏矩阵)。或者,矩阵的行为可以从函数隐式定义。当使用隐含定义时,condest需要以下函数:

• Afcn(旗帜,x)必须返回
  • 尺寸n属于A如果旗帜是暗淡的
  • true ifA是一个真正的运算符,如果旗帜是真实的
  • 结果A * x如果旗帜是“nottransp”
  • 结果A' *x如果旗帜是“transp”
• Ainvfcn(旗帜,x)必须返回
  • 尺寸n属于inv(A)如果旗帜是暗淡的
  • true ifinv(A)是一个真正的运算符,如果旗帜是真实的
  • 结果inv(A) *x如果旗帜是“nottransp”
  • 结果inv(A)' *x如果旗帜是“transp”

任何参数p1,第2页,…是的附加参数Afcn(旗帜,x,p1,第2页, …)和Ainvfcn(旗帜,x,p1,第2页, …).

主要输出是1-范数条件数估计cest.

可选的第二个输出v是近似零向量;它满足方程标准A*v,1)==范数(A,1)*标准(v1.cest.

算法注释:condest使用随机化算法来近似1-范数。因此,如果需要一致的结果状态的应在调用之前修复condest.

参考文献:

• N.J.Higham和F.Tisseur,矩阵1-模估计的块算法及其在1-模伪谱中的应用SIMAX第21卷第4号,第1185–1201页。https://dx.doi.org/10.1137/S0895479899356080
• N.J.Higham和F.Tisseur,矩阵1-模估计的块算法及其在1-模伪谱中的应用. https://citeseer.ist.psu.edu/223007.html

详见: cond,rcond,标准,标准1,normest.

 

:spparms ()¶

:vals= spparms ()¶

:[钥匙,vals] = spparms ()¶

:val= spparms (钥匙)¶

:spparms (vals)¶

:spparms 默认¶

:spparms 牢固的¶

:spparms (钥匙,val)¶

查询或设置稀疏解算器和因子分解函数使用的参数。

上面的前四个调用获取有关当前设置的信息,而其他调用则更改当前设置。参数存储为成对的键和值,其中值都是浮点值,键是以下字符串之一:

‘斯普莫尼’

打印解算器的调试信息级别(默认为0)

‘ths_rel’

包括兼容性。未使用。(默认值1)

‘ths_abs’

包括兼容性。未使用。(默认值1)

‘精确_d’

包括兼容性。未使用。(默认为0)

‘supernd’

包括兼容性。未使用。(默认值3)

‘减少’

包括兼容性。未使用。(默认值3)

‘wh_frac’

包括兼容性。未使用。(默认0.5)

‘自动mmd’

符号LU/QR以及“\”和“/”运算符是否将自动使用稀疏性保留mmd函数(默认值1)

‘autoamd’

符号LU以及“\”和“/”运算符是否将自动使用保留稀疏性的amd函数(默认值1)

‘piv_tol’

的枢轴公差UMFPACK解算器(默认值为0.1)

‘sym_tol’

的枢轴公差UMFPACK对称解算器(默认值0.001)

‘类登’

带状矩阵中非零元素在被处理之前的密度LAPACK 带状解算器(默认值为0.5)

‘umfpack’

符号是否UMFPACK或mmd解算器用于LU、“\”和“/”操作(默认值1)

可以使用设置单个关键点的值spparms(钥匙,val)。默认值可以使用特殊关键字恢复默认.特殊关键字tight可用于设置mmd解算器,以尝试以较长运行时间为潜在代价的稀疏解算。

详见: chol,colamd,lu,qr,symamd.

 

:</p>= sprank (S)¶

计算稀疏矩阵的结构体秩S.

注意,在这个计算中,只有矩阵的结构体是基于Dulmage-Mentelsohn排列来块三角形的S以为界sprank(S等级S)。忽略浮点错误sprank(S等级S).

详见: dmperm.

 

:[计数,h,父级对象亲,邮递,R] = 符号事实 (S)¶

:[…] = 符号事实 (S,type)¶

:[…] = 符号事实 (S,type,mode)¶

对稀疏矩阵执行符号分解分析S.

输入变量为

S

S是实数或复数稀疏矩阵。

type

是因子分解的类型,可以是

“sym”默认

因式分解S.假设S是对称的并且使用矩阵的上三角部分。

“col”

因式分解S' *S.

一行

因式分解S * S'.

“哦”

因式分解S'.假设S是对称的,并使用矩阵的下三角部分。

mode

什么时候mode是未指定的,返回的Cholesky因子分解R如果mode是lower或L然后返回共轭转置R'这是一个较低的三角因子。共轭转置版本更快,使用更少的内存,但仍然为所有其他输出返回相同的值:计数,h,父级对象亲和邮递.

输出变量为:

计数

Cholesky因子分解的行数从type.使用执行真因子分解的计算难度chol是总和计数2..

h

消除树的高度。

父级对象亲

消除树本身。

邮递

一个稀疏布尔矩阵,其结构体是从下式确定的Holesky因子分解的结构体type.

详见: chol,etree,树状布局.

对于非平方矩阵,用户还可以使用spaugment函数来找到线性方程的最小二乘解。

 

:s= spaugment (A,c)¶

创建的增广矩阵A.

这是从给出的

[ceyem,m),A;A',零(n,n)]

这与的最小二乘解有关A\b通过

s* [r/c; x] =[b,零(n,列(b)) ]

这里的r是残差

r=b-A * x

作为矩阵s是对称不定的,可以用分解lu,因此可以在不需要的情况下找到最小范数解qr因子分解。从于残余误差零(m,m)对于不确定的问题,例如

m=11;n=10;mn=最大值(m,n);A=spdiags([ones(mn,1),10*ones(mn1),-ones(mn,1])],[-1,0,1],m,n);x0=A\ones(m,1);s=间距(A);[L,U,P,Q]=lu(s);x1=Q*(U\(L\(P*[一(m,1);零(n,1)]));x1=x1(结束-n+1:结束);

要找到一个超定问题的解,需要估计剩余误差r因此,使用spaugment作用

通常,左除法运算符比使用spaugment作用

详见: mldivide.

最后,函数eigs可以用于基于选择标准来计算有限数量的特征值和特征向量,对于svds其计算有限数量的奇异值和向量。

 

:d= eigs (A)¶

:d= eigs (A,k)¶

:d= eigs (A,k,西格玛)¶

:d= eigs (A,k,西格玛,opts)¶

:d= eigs (A,B)¶

:d= eigs (A,B,k)¶

:d= eigs (A,B,k,西格玛)¶

:d= eigs (A,B,k,西格玛,opts)¶

:d= eigs (Af,n)¶

:d= eigs (Af,n,k)¶

:d= eigs (Af,n,k,西格玛)¶

:d= eigs (Af,n,k,西格玛,opts)¶

:d= eigs (Af,n,B)¶

:d= eigs (Af,n,B,k)¶

:d= eigs (Af,n,B,k,西格玛)¶

:d= eigs (Af,n,B,k,西格玛,opts)¶

:[五、,D] = eigs (…)¶

:[五、,D,旗帜] = eigs (…)¶

根据选择标准计算有限数量的特征值和特征向量。

默认情况下,eigs求解相应的特征向量所在的方程。如果给定正定矩阵B然后eigs求解一般特征值方程

输入A是维度的平方矩阵n通过n典型的A也是大而稀疏的。

输入B对于广义特征值问题,是一个大小与A(n通过n). 典型的B也是大而稀疏的。

要计算的特征值和特征向量的数量从下式给出k并且默认为6。

参数西格玛确定返回哪些特征值。西格玛可以是标量或字符串。当西格玛是向量k最接近的特征值西格玛返回。如果西格玛是字符串,它必须是以下值之一。

“lm”

最大幅值(默认值)。

sm

最小幅值。

“la”

最大代数(仅适用于实对称问题)。

“sa”

最小代数(仅适用于实对称问题)。

是

两端,如果k是奇数(仅对真正的对称问题有效)。

“lr”

最大实部(仅适用于复杂或不对称问题)。

“sr”

最小实部(仅适用于复杂或不对称问题)。

“李”

最大的想象部分(仅适用于复杂或不对称问题)。

“si”

最小的想象部分(仅适用于复杂或不对称问题)。

如果opts是一个定义可能参数的结构体eigs应该使用。的字段opts结构体为:

issym

如果Af则该标志(真/假)确定函数是否Af定义了一个对称问题。如果矩阵A给出了。默认值为false。

以色列

如果Af则该标志(真/假)确定函数是否Af定义了一个真正的问题。如果矩阵A给出了。默认值为true。

tol

定义所需的收敛公差,计算为tol*标准(A)。默认为eps.

maxit

最大迭代次数。默认值为300。

</p>

要使用的Lanczos基向量的数量。更多的向量将导致更快的收敛,但会更多地使用内存。的最佳值</p>取决于问题,应在范围内k1.到n。默认值为2.k.

v0

算法的起始向量。接近最终解算器的初始向量将加快收敛速度。默认值为ARPACK随机生成一个起始向量。如果指定,v0必豆n-by-1向量,其中n=行(A).

disp

diag打印输出的级别(0|1|2)。如果disp为0,则禁用diag。默认值为0。

cholB

如果正在计算广义特征值问题,则此标志(真/假)指定B输入表示chol(B)或者只是矩阵B。默认值为false。

permB

的Cholesky因子分解的置换向量B如果cholB是真的。它是通过[R,~,permB]=chol(B,“向量”)。默认为1.n.

也可以表示A通过表示为Af. Af必须后跟标量参数n定义接受的向量自变量的长度Af. Af可以是functionhandle、内联函数或字符串。当Af是一个没有要使用的函数名称的字符串。

Af是形式的函数y=Af(x)其中的要求返回值Af从的值决定西格玛。四种可能的形式是

A*x

如果西格玛不是给定的,或者是除“sm”之外的字符串。

A\x

如果西格玛是0或“sm”。

(A-西格玛*I)\x

如果西格玛是不等于0的标量;我身份矩阵的大小是否与A.

(A-西格玛*B)\x

对于一般的特征值问题。

返回参数及其形式取决于指定的返回参数的数量。对于单个返回参数,列向量d的长度k返回时包含k已经查找的特征值。对于两个返回自变量,五、是一个n通过kmatrix其列为k与返回的特征值相对应的特征向量。特征值本身返回D以的形式k通过k矩阵,其中对角线上的元素是本征值。

第三个返回参数旗帜返回收敛的状态。如果旗帜如果为0,则所有特征值都已收敛。任何其他值都表示无法收敛。

编程注意事项:对于小问题,n<500,考虑使用eig(完整(A)).

如果ARPACK无法收敛,请考虑增加Lanczos向量的数量(opt-p),增加迭代次数(opt.maxiter),或减小公差(opt.tol)。

参考:此函数基于ARPACK包,从R编写。Lehoucq、K.Maschhoff、D.Sorensen和C.Yang。有关更多信息,详见http://www.caam.rice.edu/software/ARPACK/.

详见: eig,svds.

 

:s= svds (A)¶

:s= svds (A,k)¶

:s= svds (A,k,西格玛)¶

:s= svds (A,k,西格玛,opts)¶

:[u,s,v] = svds (…)¶

:[u,s,v,旗帜] = svds (…)¶

求矩阵的几个奇异值A.

奇异值的计算使用

[m,n大小A);s=eigs([稀疏(m,m),A;A稀疏的n,n)])

返回的特征值eigs对应于的奇异值A。要计算的奇异值的数量从k并且默认为6。

参数西格玛指定要查找的奇异值。什么时候西格玛是字符串L,默认情况下,的最大奇异值A找到。否则西格玛必须是实数标量并且奇异值最接近西格玛找到。作为推论,西格玛= 0找到最小的奇异值。请注意,对于相对较小的值西格玛,则有可能找不到所要求的奇异值的数目。那样的话西格玛应该增加。

opts是一个定义参数的结构体svds将传递给eigs。该结构体的可能字段记录在eigs默认情况下,svds设置以下三个字段:

tol

奇异值所需的收敛容差。默认值为1e-10。eigs已通过tol/sqrt(2).

maxit

最大迭代次数。默认值为300。

disp

diag打印输出的级别(0|1|2)。如果disp为0,则禁用diag。默认值为0。

如果指定多个输出,则svds将返回的奇异值分解的近似A

A_大约=u*s*v'

这里的A_近似是一个大小矩阵A但只有等级k.

旗帜如果算法已成功收敛,则返回0,否则返回1。收敛性的测试是

标准A*v-u*s1.tol标准A1.

svds最好是从一个大的稀疏矩阵中只找到几个奇异值。否则,svd(完整(A))可能会更有效率。

详见: svd,eigs.