34 Julia中的矩阵计算功能(*)

您所在的位置:网站首页 3个矩阵相乘的逆矩阵 34 Julia中的矩阵计算功能(*)

34 Julia中的矩阵计算功能(*)

2023-09-11 17:07| 来源: 网络整理| 查看: 265

34 Julia中的矩阵计算功能(*) 34.1 基本矩阵运算

Julia是一个向Matlab软件一样强调科学计算, 尤其是线性代数运算的编程语言, 所以其矩阵功能比R还要强大。

Julia的基本矩阵构造、加减法、乘法、求逆、求解线性方程组等运算, 见第@(jintro)章。

size(A)返回包含行数和列数的二元组(tuple)。

eye(n)生成n阶的元素类型为Float64的单位阵。 ones(n, m)生成\(n \times m\)的元素类型为Float64的都是1的矩阵, zeros(n,m)生成\(n \times m\)的元素类型为Float64的都是0的矩阵。

A'表示A的共轭转置。 A + B和A - B表示矩阵加减, A * B表示矩阵乘法。 inv(A)返回逆矩阵。 A \ b求解\(A x = b\)中的\(x\)。

dot(x, y)计算两个向量的内积, 与sum(x .* y)相同。 norm(x)求向量的欧式模, 等于sqrt(sum(x .^ 2))。

factorize(A)对A计算适合的分解, 自动判断A的特殊结构。

det(A)求行列式。 trace(A)计算迹。

eigvals(A)求特征值。 eig(A)返回由特征值、特征向量组成的二元组, 特征值储存在一个向量中, 特征向量储存在矩阵的各列中。

eig(A, B)求解广义特征值问题\(A v = \lambda B v\)。

rank(A)求秩。

Julia的线性代数模块提供了丰富的矩阵类型和矩阵计算功能, 详见Base.LinAlg模块的文档。

34.2 用Givens变换作QR分解示例 using Distributions # 随机数函数 srand(101) # 指定起始种子 A=round.(rand(4,4), 2) ## 4×4 Array{Float64,2}: ## 0.09 1.0 0.97 0.42 ## 0.83 0.86 0.62 0.8 ## 0.63 0.8 0.75 0.23 ## 0.16 0.95 0.73 0.07 ## Givens变换矩阵,在i1,i2两行和两列作用,将(a, b)长度集中到a所在元素 function givensmat(n, a, b, i1, i2) submat = [a b; -b a]/sqrt(a*a+b*b) resmat = eye(n) resmat[[i1,i2], [i1,i2]] = submat resmat end n=4 Q1 = givensmat(n, A[1,1], A[2,1], 1, 2) ## 4×4 Array{Float64,2}: ## 0.107802 0.994172 0.0 0.0 ## -0.994172 0.107802 0.0 0.0 ## 0.0 0.0 1.0 0.0 ## 0.0 0.0 0.0 1.0 A12 = Q1 * A ## 4×4 Array{Float64,2}: ## 0.834865 0.96279 0.720955 0.840615 ## 1.38778e-17 -0.901463 -0.89751 -0.331311 ## 0.63 0.8 0.75 0.23 ## 0.16 0.95 0.73 0.07 Q2 = givensmat(n, A12[1,1], A12[3,1], 1, 3) ## 4×4 Array{Float64,2}: ## 0.798229 0.0 0.602354 0.0 ## 0.0 1.0 0.0 0.0 ## -0.602354 0.0 0.798229 0.0 ## 0.0 0.0 0.0 1.0 A13 = Q2 * A12 ## 4×4 Array{Float64,2}: ## 1.0459 1.25041 1.02725 0.809545 ## 1.38778e-17 -0.901463 -0.89751 -0.331311 ## 0.0 0.0586429 0.164402 -0.322755 ## 0.16 0.95 0.73 0.07 Q3 = givensmat(n, A13[1,1], A[4,1], 1, 4) ## 4×4 Array{Float64,2}: ## 0.9885 0.0 0.0 0.15122 ## 0.0 1.0 0.0 0.0 ## 0.0 0.0 1.0 0.0 ## -0.15122 0.0 0.0 0.9885 A14 = Q3 * A13 ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## 1.38778e-17 -0.901463 -0.89751 -0.331311 ## 0.0 0.0586429 0.164402 -0.322755 ## -2.77556e-17 0.749989 0.566264 -0.0532239 Q4 = givensmat(n, A14[2,2], A14[3,2], 2, 3) ## 4×4 Array{Float64,2}: ## 1.0 0.0 0.0 0.0 ## 0.0 -0.997891 0.0649159 0.0 ## 0.0 -0.0649159 -0.997891 0.0 ## 0.0 0.0 0.0 1.0 A23 = Q4 * A14 ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## -1.38485e-17 0.903368 0.906289 0.30966 ## -9.00889e-19 -6.93889e-18 -0.105793 0.343581 ## -2.77556e-17 0.749989 0.566264 -0.0532239 Q5 = givensmat(n, A23[2,2], A23[4,2], 2, 4); A24 = Q5 * A23 ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## -2.83844e-17 1.17412 1.05901 0.204255 ## -9.00889e-19 -6.93889e-18 -0.105793 0.343581 ## -1.25092e-17 0.0 -0.143223 -0.238751 Q6 = givensmat(n, A24[3,3], A24[4,3], 3, 4); A34 = Q6 * A24 ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## -2.83844e-17 1.17412 1.05901 0.204255 ## 1.05971e-17 4.1227e-18 0.178059 -0.012095 ## 6.70761e-18 -5.58136e-18 0.0 0.418215

到A34已经变成了上三角矩阵。各个Givens旋转合并起来也是一个正交阵:

Q = Q6*Q5*Q4*Q3*Q2*Q1 ## 4×4 Array{Float64,2}: ## 0.085061 0.784451 0.595427 0.15122 ## 0.751748 -0.189333 -0.0183152 0.631421 ## 0.43877 -0.351868 0.556258 -0.611757 ## 0.484892 0.474318 -0.579403 -0.451877 Q'*Q ## 4×4 Array{Float64,2}: ## 1.0 -5.55112e-17 -5.55112e-17 0.0 ## -5.55112e-17 1.0 -2.22045e-16 0.0 ## -5.55112e-17 -2.22045e-16 1.0 0.0 ## 0.0 0.0 0.0 1.0 B = Q * A ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## 0.0 1.17412 1.05901 0.204255 ## 4.16334e-17 0.0 0.178059 -0.012095 ## -4.16334e-17 5.55112e-17 1.11022e-16 0.418215 maximum(B - A34) ## 2.220446049250313e-16

上面相当于\(Q A = R\),所以\(A = Q' R\)。 Julia的qr()函数可以计算QR分解:

Q1, R1 = qr(A); R1 ## 4×4 Array{Float64,2}: ## -1.05806 -1.37969 -1.12583 -0.81082 ## 0.0 1.17412 1.05901 0.204255 ## 0.0 0.0 -0.178059 0.012095 ## 0.0 0.0 0.0 -0.418215

qr()函数给出的上三角矩阵不满足对角线元素为正的要求。需要时可以将负的主对角元素所在行反号并将对应的Q的列反号。如:

R2 = R1; Q2 = Q1 R2[1,:] = -R2[1,:]; Q2[:,1] = -Q2[:,1] R2[3,:] = -R2[3,:]; Q2[:,3] = -Q2[:,3] R2[4,:] = -R2[4,:]; Q2[:,4] = -Q2[:,4] R2 ## 4×4 Array{Float64,2}: ## 1.05806 1.37969 1.12583 0.81082 ## 0.0 1.17412 1.05901 0.204255 ## -0.0 -0.0 0.178059 -0.012095 ## -0.0 -0.0 -0.0 0.418215 Q2 ## 4×4 Array{Float64,2}: ## 0.085061 0.751748 0.43877 0.484892 ## 0.784451 -0.189333 -0.351868 0.474318 ## 0.595427 -0.0183152 0.556258 -0.579403 ## 0.15122 0.631421 -0.611757 -0.451877 Q' ## 4×4 Array{Float64,2}: ## 0.085061 0.751748 0.43877 0.484892 ## 0.784451 -0.189333 -0.351868 0.474318 ## 0.595427 -0.0183152 0.556258 -0.579403 ## 0.15122 0.631421 -0.611757 -0.451877


【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3