/usr/share/julia/test/linalg/bidiag.jl is in julia-common 0.4.7-6.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 | # This file is a part of Julia. License is MIT: http://julialang.org/license
using Base.Test
import Base.LinAlg: BlasReal, BlasFloat
debug = false
n = 10 #Size of test matrix
srand(1)
debug && println("Bidiagonal matrices")
for relty in (Float32, Float64, BigFloat), elty in (relty, Complex{relty})
debug && println("elty is $(elty), relty is $(relty)")
dv = convert(Vector{elty}, randn(n))
ev = convert(Vector{elty}, randn(n-1))
b = convert(Matrix{elty}, randn(n, 2))
c = convert(Matrix{elty}, randn(n, n))
if (elty <: Complex)
dv += im*convert(Vector{elty}, randn(n))
ev += im*convert(Vector{elty}, randn(n-1))
b += im*convert(Matrix{elty}, randn(n, 2))
end
debug && println("Test constructors")
@test Bidiagonal(dv,ev,'U') == Bidiagonal(dv,ev,true)
@test_throws ArgumentError Bidiagonal(dv,ev,'R')
@test_throws DimensionMismatch Bidiagonal(dv,ones(elty,n),true)
@test_throws ArgumentError Bidiagonal(dv,ev)
#getindex and size
BD = Bidiagonal(dv,ev,true)
@test_throws BoundsError BD[n+1,1]
@test_throws ArgumentError size(BD,0)
@test size(BD,3) == 1
debug && println("show")
dstring = sprint(Base.print_matrix,BD.dv')
estring = sprint(Base.print_matrix,BD.ev')
@test sprint(show,BD) == "$(summary(BD)):\n diag:$dstring\n super:$estring"
BD = Bidiagonal(dv,ev,false)
@test sprint(show,BD) == "$(summary(BD)):\n diag:$dstring\n sub:$estring"
debug && println("Test upper and lower bidiagonal matrices")
for isupper in (true, false)
debug && println("isupper is: $(isupper)")
T = Bidiagonal(dv, ev, isupper)
@test size(T, 1) == size(T, 2) == n
@test size(T) == (n, n)
@test full(T) == diagm(dv) + diagm(ev, isupper?1:-1)
@test Bidiagonal(full(T), isupper) == T
@test big(T) == T
@test full(abs(T)) == abs(diagm(dv)) + abs(diagm(ev, isupper?1:-1))
@test full(real(T)) == real(diagm(dv)) + real(diagm(ev, isupper?1:-1))
@test full(imag(T)) == imag(diagm(dv)) + imag(diagm(ev, isupper?1:-1))
z = zeros(elty, n)
debug && println("Idempotent tests")
for func in (conj, transpose, ctranspose)
@test func(func(T)) == T
end
debug && println("triu and tril")
@test istril(Bidiagonal(dv,ev,'L'))
@test !istril(Bidiagonal(dv,ev,'U'))
@test tril!(Bidiagonal(dv,ev,'U'),-1) == Bidiagonal(zeros(dv),zeros(ev),'U')
@test tril!(Bidiagonal(dv,ev,'L'),-1) == Bidiagonal(zeros(dv),ev,'L')
@test tril!(Bidiagonal(dv,ev,'U'),-2) == Bidiagonal(zeros(dv),zeros(ev),'U')
@test tril!(Bidiagonal(dv,ev,'L'),-2) == Bidiagonal(zeros(dv),zeros(ev),'L')
@test tril!(Bidiagonal(dv,ev,'U'),1) == Bidiagonal(dv,ev,'U')
@test tril!(Bidiagonal(dv,ev,'L'),1) == Bidiagonal(dv,ev,'L')
@test tril!(Bidiagonal(dv,ev,'U')) == Bidiagonal(dv,zeros(ev),'U')
@test tril!(Bidiagonal(dv,ev,'L')) == Bidiagonal(dv,ev,'L')
@test_throws ArgumentError tril!(Bidiagonal(dv,ev,'U'),n+1)
@test istriu(Bidiagonal(dv,ev,'U'))
@test !istriu(Bidiagonal(dv,ev,'L'))
@test triu!(Bidiagonal(dv,ev,'L'),1) == Bidiagonal(zeros(dv),zeros(ev),'L')
@test triu!(Bidiagonal(dv,ev,'U'),1) == Bidiagonal(zeros(dv),ev,'U')
@test triu!(Bidiagonal(dv,ev,'U'),2) == Bidiagonal(zeros(dv),zeros(ev),'U')
@test triu!(Bidiagonal(dv,ev,'L'),2) == Bidiagonal(zeros(dv),zeros(ev),'L')
@test triu!(Bidiagonal(dv,ev,'U'),-1) == Bidiagonal(dv,ev,'U')
@test triu!(Bidiagonal(dv,ev,'L'),-1) == Bidiagonal(dv,ev,'L')
@test triu!(Bidiagonal(dv,ev,'L')) == Bidiagonal(dv,zeros(ev),'L')
@test triu!(Bidiagonal(dv,ev,'U')) == Bidiagonal(dv,ev,'U')
@test_throws ArgumentError triu!(Bidiagonal(dv,ev,'U'),n+1)
debug && println("Linear solver")
Tfull = full(T)
condT = cond(map(Complex128,Tfull))
x = T \ b
tx = Tfull \ b
@test_throws DimensionMismatch Base.LinAlg.naivesub!(T,ones(elty,n+1))
@test norm(x-tx,Inf) <= 4*condT*max(eps()*norm(tx,Inf), eps(relty)*norm(x,Inf))
@test_throws DimensionMismatch T \ ones(elty,n+1,2)
@test_throws DimensionMismatch T.' \ ones(elty,n+1,2)
@test_throws DimensionMismatch T' \ ones(elty,n+1,2)
if relty != BigFloat
x = T.'\c.'
tx = Tfull.' \ c.'
@test norm(x-tx,Inf) <= 4*condT*max(eps()*norm(tx,Inf), eps(relty)*norm(x,Inf))
@test_throws DimensionMismatch T.'\b.'
x = T'\c.'
tx = Tfull' \ c.'
@test norm(x-tx,Inf) <= 4*condT*max(eps()*norm(tx,Inf), eps(relty)*norm(x,Inf))
@test_throws DimensionMismatch T'\b.'
x = T\c.'
tx = Tfull\c.'
@test norm(x-tx,Inf) <= 4*condT*max(eps()*norm(tx,Inf), eps(relty)*norm(x,Inf))
@test_throws DimensionMismatch T\b.'
end
debug && println("Round,float,trunc,ceil")
if elty <: BlasReal
@test floor(Int,T) == Bidiagonal(floor(Int,T.dv),floor(Int,T.ev),T.isupper)
@test trunc(Int,T) == Bidiagonal(trunc(Int,T.dv),trunc(Int,T.ev),T.isupper)
@test round(Int,T) == Bidiagonal(round(Int,T.dv),round(Int,T.ev),T.isupper)
@test ceil(Int,T) == Bidiagonal(ceil(Int,T.dv),ceil(Int,T.ev),T.isupper)
end
debug && println("Generic Mat-vec ops")
@test_approx_eq T*dv Tfull*dv
@test_approx_eq T'*dv Tfull'*dv
if relty != BigFloat # not supported by pivoted QR
@test_approx_eq T/dv' Tfull/dv'
end
debug && println("Diagonals")
@test diag(T,2) == zeros(elty, n-2)
@test_throws ArgumentError diag(T,n+1)
debug && println("Eigensystems")
d1, v1 = eig(T)
d2, v2 = eig(map(elty<:Complex ? Complex128 : Float64,Tfull))
@test_approx_eq isupper?d1:reverse(d1) d2
if elty <: Real
Test.test_approx_eq_modphase(v1, isupper?v2:v2[:,n:-1:1])
end
debug && println("Singular systems")
if (elty <: BlasReal)
@test_approx_eq full(svdfact(T)) full(svdfact!(copy(Tfull)))
@test_approx_eq svdvals(Tfull) svdvals(T)
u1, d1, v1 = svd(Tfull)
u2, d2, v2 = svd(T)
@test_approx_eq d1 d2
if elty <: Real
Test.test_approx_eq_modphase(u1, u2)
Test.test_approx_eq_modphase(v1, v2)
end
@test_approx_eq_eps 0 vecnorm(u2*diagm(d2)*v2'-Tfull) n*max(n^2*eps(relty), vecnorm(u1*diagm(d1)*v1' - Tfull))
@inferred svdvals(T)
@inferred svd(T)
end
debug && println("Binary operations")
@test -T == Bidiagonal(-T.dv,-T.ev,T.isupper)
@test convert(elty,-1.0) * T == Bidiagonal(-T.dv,-T.ev,T.isupper)
@test T * convert(elty,-1.0) == Bidiagonal(-T.dv,-T.ev,T.isupper)
for isupper2 in (true, false)
dv = convert(Vector{elty}, randn(n))
ev = convert(Vector{elty}, randn(n-1))
T2 = Bidiagonal(dv, ev, isupper2)
Tfull2 = full(T2)
for op in (+, -, *)
@test_approx_eq full(op(T, T2)) op(Tfull, Tfull2)
end
end
debug && println("Inverse")
@test_approx_eq inv(T)*Tfull eye(elty,n)
end
@test Matrix{Complex{Float64}}(BD) == BD
end
# Issue 10742 and similar
let A = Bidiagonal([1,2,3], [0,0], true)
@test istril(A)
@test isdiag(A)
end
#test promote_rule
A = Bidiagonal(ones(Float32,10),ones(Float32,9),true)
B = rand(Float64,10,10)
C = Tridiagonal(rand(Float64,9),rand(Float64,10),rand(Float64,9))
@test promote(B,A) == (B,convert(Matrix{Float64},full(A)))
@test promote(C,A) == (C,Tridiagonal(zeros(Float64,9),convert(Vector{Float64},A.dv),convert(Vector{Float64},A.ev)))
|