Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

More special, eigen, and bitarray tests #12121

Merged
merged 4 commits into from
Jul 14, 2015
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
11 changes: 9 additions & 2 deletions base/linalg/eigen.jl
Original file line number Diff line number Diff line change
Expand Up @@ -87,14 +87,21 @@ function eigvals{T<:Number}(x::T; kwargs...)
return imag(val) == 0 ? [real(val)] : [val]
end

# TO DO: Put message about not being able to sort complex numbers back in!
#Computes maximum and minimum eigenvalue
function eigmax(A::Union{Number, StridedMatrix}; permute::Bool=true, scale::Bool=true)
v = eigvals(A, permute = permute, scale = scale)
iseltype(v,Complex) ? error("DomainError: complex eigenvalues cannot be ordered") : maximum(v)
if iseltype(v,Complex)
throw(DomainError())
end
maximum(v)
end
function eigmin(A::Union{Number, StridedMatrix}; permute::Bool=true, scale::Bool=true)
v = eigvals(A, permute = permute, scale = scale)
iseltype(v,Complex) ? error("DomainError: complex eigenvalues cannot be ordered") : minimum(v)
if iseltype(v,Complex)
throw(DomainError())
end
minimum(v)
end

inv(A::Eigen) = A.vectors * inv(Diagonal(A.values)) / A.vectors
Expand Down
7 changes: 7 additions & 0 deletions test/bitarray.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1168,6 +1168,7 @@ b1 = triu(bitrand(n2, n1))
b1 = bitrand(n1,n1)
b1 |= b1.'
@check_bit_operation issym(b1) Bool
@check_bit_operation ishermitian(b1) Bool

b1 = bitrand(n1)
b2 = bitrand(n2)
Expand Down Expand Up @@ -1236,6 +1237,12 @@ uA = bitunpack(A)
@test svd(A) == svd(uA)
@test qr(A) == qr(uA)

#gradient
A = bitrand(10)
fA = bitunpack(A)
@test gradient(A) == gradient(fA)
@test gradient(A,1.0) == gradient(fA,1.0)

#diag and diagm

v = bitrand(10)
Expand Down
5 changes: 5 additions & 0 deletions test/linalg1.jl
Original file line number Diff line number Diff line change
Expand Up @@ -149,6 +149,11 @@ debug && println("non-symmetric eigen decomposition")

num_fact = eigfact(one(eltya))
@test num_fact.values[1] == one(eltya)
h = a + a'
@test_approx_eq minimum(eigvals(h)) eigmin(h)
@test_approx_eq maximum(eigvals(h)) eigmax(h)
@test_throws DomainError eigmin(a - a')
@test_throws DomainError eigmax(a - a')
end

debug && println("symmetric generalized eigenproblem")
Expand Down
41 changes: 32 additions & 9 deletions test/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -106,27 +106,41 @@ end
@test_approx_eq erfi(1+2im) -0.011259006028815025076+1.0036063427256517509im
@test_approx_eq dawson(1+2im) -13.388927316482919244-11.828715103889593303im

for x in logspace(-200, -0.01)
@test_approx_eq_eps erf(erfinv(x)) x 1e-12*x
@test_approx_eq_eps erf(erfinv(-x)) -x 1e-12*x
@test_approx_eq_eps erfc(erfcinv(2*x)) 2*x 1e-12*x
if x > 1e-20
xf = Float32(x)
@test_approx_eq_eps erf(erfinv(xf)) xf 1e-5*xf
@test_approx_eq_eps erf(erfinv(-xf)) -xf 1e-5*xf
@test_approx_eq_eps erfc(erfcinv(2xf)) 2xf 1e-5*xf
for elty in [Float32,Float64]
for x in logspace(-200, -0.01)
@test_approx_eq_eps erf(erfinv(x)) x 1e-12*x
@test_approx_eq_eps erf(erfinv(-x)) -x 1e-12*x
@test_approx_eq_eps erfc(erfcinv(2*x)) 2*x 1e-12*x
if x > 1e-20
xf = Float32(x)
@test_approx_eq_eps erf(erfinv(xf)) xf 1e-5*xf
@test_approx_eq_eps erf(erfinv(-xf)) -xf 1e-5*xf
@test_approx_eq_eps erfc(erfcinv(2xf)) 2xf 1e-5*xf
end
end
@test erfinv(one(elty)) == Inf
@test erfinv(-one(elty)) == -Inf
@test_throws DomainError erfinv(2.0*one(elty))

@test erfcinv(zero(elty)) == Inf
@test_throws DomainError erfcinv(-one(elty))
end

@test erfinv(one(Int)) == erfinv(1.0)
@test erfcinv(one(Int)) == erfcinv(1.0)

# airy
@test_approx_eq airy(1.8) 0.0470362168668458052247
@test_approx_eq airy(1.8) airyai(1.8)
@test_approx_eq airyprime(1.8) -0.0685247801186109345638
@test_approx_eq airyaiprime(1.8) airyprime(1.8)
@test_approx_eq airybi(1.8) 2.595869356743906290060
@test_approx_eq airybiprime(1.8) 2.98554005084659907283
@test_throws Base.Math.AmosException airy(200im)
@test_throws Base.Math.AmosException airybi(200)
@test_throws ArgumentError airy(5,one(Complex128))
z = 1.8 + 1.0im
@test_approx_eq airyx(z) airyx(0,z)
@test_approx_eq airyx(0, z) airy(0, z) * exp(2/3 * z * sqrt(z))
@test_approx_eq airyx(1, z) airy(1, z) * exp(2/3 * z * sqrt(z))
@test_approx_eq airyx(2, z) airy(2, z) * exp(-abs(real(2/3 * z * sqrt(z))))
Expand Down Expand Up @@ -218,6 +232,7 @@ y33 = bessely(3,3.)
@test_throws Base.Math.AmosException bessely(200.5,0.1)
@test_throws DomainError bessely(0.4,-1.0)
@test_throws DomainError bessely(0.4,Float32(-1.0))
@test_throws DomainError bessely(1,Float32(-1.0))

# issue #6653
for f in (besselj,bessely,besseli,besselk,hankelh1,hankelh2)
Expand Down Expand Up @@ -473,6 +488,14 @@ with_bigfloat_precision(10_000) do
@test log(2,big(2)^400) == 400
end

for T in (Float32,Float64)
@test log(zero(T)) == -Inf
@test isnan(log(NaN))
@test_throws DomainError log(-one(T))
@test log1p(-one(T)) == -Inf
@test isnan(log1p(NaN))
@test_throws DomainError log1p(-2*one(T))
end
# test vectorization of 2-arg vectorized functions
binary_math_functions = [
copysign, flipsign, log, atan2, hypot, max, min,
Expand Down