2016-08-19 89 views
2

我有列数等于行数。对角线等于零。我怎样才能建立这个矩阵?茱莉亚的下三角矩阵

#mat 
#  [,1] [,2] [,3] [,4] 
#[1,] 0 NA NA NA 
#[2,] 1 0 NA NA 
#[3,] 2 4 0 NA 
#[4,] 3 5 6 0 

我想这

x=rand(4,4) 
4x4 Array{Float64,2}: 
0.60064 0.917443 0.561744 0.135717 
0.106728 0.72391 0.0894174 0.0656103 
0.410262 0.953857 0.844697 0.0375045 
0.476771 0.778106 0.469514 0.398846 

c=LowerTriangular(x) 

4x4 LowerTriangular{Float64,Array{Float64,2}}: 
0.60064 0.0  0.0  0.0  
0.106728 0.72391 0.0  0.0  
0.410262 0.953857 0.844697 0.0  
0.476771 0.778106 0.469514 0.398846 

但我米寻找这样的事情

c=LowerTriangular(x) 
4x4 LowerTriangular{Float64,Array{Float64,2}}: 
0.0 NA  NA  NA  
0.106728 0.0 NA  NA  
0.410262 0.953857 0.0 NA  
0.476771 0.778106 0.469514 0 

对角线应该等于零。

+1

只是FYI茱莉亚不使用NA值。我们有NaN –

+0

x-ref:https://groups.google.com/forum/#!topic/julia-users/VGbdlfbRfbc – StefanKarpinski

+1

此外,请考虑仍使用零而不是NA(或NaN)。为什么? 1.线性代数的意义大部分是使用NA丢失的,这将阻止任何现有的线性代数代码的使用。已知新生感染有传染性和产生病毒。如果你更多地描述用例,这个选择可能会更清楚。 –

回答

1

这里是关于朱莉娅用户list的东西,其灵感来自代码由斯特凡斯基:

function vec2ltri_alt{T}(v::AbstractVector{T}, z::T=zero(T)) 
    n = length(v) 
    v1 = vcat(0,v) 
    s = round(Int,(sqrt(8n+1)-1)/2) 
    s*(s+1)/2 == n || error("vec2utri: length of vector is not triangular") 
    s+=1 
    [ i>j ? v1[round(Int, j*(j-1)/2+i)] : (i == j ? z : NaN) for i=1:s, j=1:s ] 
end 

julia> vec2ltri_alt(collect(1:6)) 
4x4 Array{Any,2}: 
0 NaN NaN NaN 
1 0 NaN NaN 
2 3 0 NaN 
3 4 6 0 

注:如果需要,请查看official documentation在三元运算符上更清楚地说明了这里的? ... :语法是怎么回事。

对于那些寻找一个更“标准”对角矩阵解决方案:

下面是创建一个更标准的解决方案的一个版本:如果你想上三角

function vec2ltri{T}(v::AbstractVector{T}, z::T=zero(T)) 
    n = length(v) 
    s = round(Int,(sqrt(8n+1)-1)/2) 
    s*(s+1)/2 == n || error("vec2utri: length of vector is not triangular") 
    [ i>=j ? v[round(Int, j*(j-1)/2+i)] : z for i=1:s, j=1:s ] 
end 

a = vec2ltri(collect(1:6)) 

julia> a = vec2ltri(collect(1:6)) 
3x3 Array{Int64,2}: 
1 0 0 
2 3 0 
3 4 6 

julia> istril(a) ## verify matrix is lower triangular 
true 

代替更低,只需将i<=j更改为i>=j即可。

其他随机工具还要注意函数tril!(a)将把给定的矩阵转换为下三角形,用零替换主对角线上的所有东西。有关此功能的更多信息,请参阅Julia documentation以及其他各种相关工具。

+0

thn = ank你这就是我要找的东西 – vincet

3

您可能想要使用列表理解。但如果你能在关于你想要做什么的问题中提供更多信息,那将是一件好事。

numrows =4 
numcols = 4 
[ x>y ? 1 : (x == y ? 0 : NaN) for x in 1:numrows, y in 1:numcols] 

,这将给:

0 NaN NaN NaN 
1 0 NaN NaN 
1 1 0 NaN 
1 1 1 0 

对于任何数量的行和列。然后你可以从那里工作。

见文档的列表内涵&条件语句:

http://docs.julialang.org/en/release-0.4/manual/arrays/#comprehensions

http://docs.julialang.org/en/release-0.4/manual/control-flow/#man-conditional-evaluation

+0

谢谢大家的澄清 – vincet

+0

现在,我想探索一个foor循环只是下三角矩阵。我怎样才能有效地形成快速优化的观点,因为我有一个5000乘5000的矩阵。谢谢 – vincet