diff --git a/NEWS.md b/NEWS.md index 19dd81a117fed..52a4503ba44d6 100644 --- a/NEWS.md +++ b/NEWS.md @@ -355,6 +355,8 @@ Deprecated or removed full path if you need access to executables or libraries in the `JULIA_HOME` directory, e.g. `joinpath(JULIA_HOME, "7z.exe")` for `7z.exe` ([#21540]). + * `sqrtm` has been deprecated in favor of `sqrt` ([#23504]). + * `expm` has been deprecated in favor of `exp` ([#23233]). * `logm` has been deprecated in favor of `log` ([#23505]). diff --git a/base/deprecated.jl b/base/deprecated.jl index d20e243f921ba..2dc1cec006b6f 100644 --- a/base/deprecated.jl +++ b/base/deprecated.jl @@ -253,7 +253,7 @@ for f in ( # base/math.jl :cbrt, :sinh, :cosh, :tanh, :atan, :asinh, :exp2, :expm1, :exp10, :sin, :cos, :tan, :asin, :acos, :acosh, :atanh, - :log2, :log10, :lgamma, #=:log1p,=# :sqrt, + :log2, :log10, :lgamma, #=:log1p,=# # base/floatfuncs.jl :abs, :abs2, :angle, :isnan, :isinf, :isfinite, # base/complex.jl @@ -1620,6 +1620,9 @@ function Tridiagonal(dl::AbstractVector{Tl}, d::AbstractVector{Td}, du::Abstract Tridiagonal(map(v->convert(Vector{promote_type(Tl,Td,Tu)}, v), (dl, d, du))...) end +# deprecate sqrtm in favor of sqrt +@deprecate sqrtm sqrt + # deprecate expm in favor of exp @deprecate expm! exp! @deprecate expm exp diff --git a/base/exports.jl b/base/exports.jl index 75607641a48f8..addd53d810f7f 100644 --- a/base/exports.jl +++ b/base/exports.jl @@ -609,7 +609,6 @@ export schur, schurfact!, schurfact, - sqrtm, svd, svdfact!, svdfact, diff --git a/base/linalg/dense.jl b/base/linalg/dense.jl index 3298f0c42c737..dcdd94b60ffd6 100644 --- a/base/linalg/dense.jl +++ b/base/linalg/dense.jl @@ -374,8 +374,8 @@ function schurpow(A::AbstractMatrix, p) retmat = A ^ floor(p) # Real part if p - floor(p) == 0.5 - # special case: A^0.5 === sqrtm(A) - retmat = retmat * sqrtm(A) + # special case: A^0.5 === sqrt(A) + retmat = retmat * sqrt(A) else retmat = retmat * powm!(UpperTriangular(float.(A)), real(p - floor(p))) end @@ -385,8 +385,8 @@ function schurpow(A::AbstractMatrix, p) R = S ^ floor(p) # Real part if p - floor(p) == 0.5 - # special case: A^0.5 === sqrtm(A) - R = R * sqrtm(S) + # special case: A^0.5 === sqrt(A) + R = R * sqrt(S) else R = R * powm!(UpperTriangular(float.(S)), real(p - floor(p))) end @@ -618,7 +618,7 @@ function log(A::StridedMatrix{T}) where T end """ - sqrtm(A) + sqrt(A::AbstractMatrix) If `A` has no negative real eigenvalues, compute the principal matrix square root of `A`, that is the unique matrix ``X`` with eigenvalues having positive real part such that @@ -642,40 +642,38 @@ julia> A = [4 0; 0 4] 4 0 0 4 -julia> sqrtm(A) +julia> sqrt(A) 2×2 Array{Float64,2}: 2.0 0.0 0.0 2.0 ``` """ -function sqrtm(A::StridedMatrix{<:Real}) +function sqrt(A::StridedMatrix{<:Real}) if issymmetric(A) - return full(sqrtm(Symmetric(A))) + return full(sqrt(Symmetric(A))) end n = checksquare(A) if istriu(A) - return full(sqrtm(UpperTriangular(A))) + return full(sqrt(UpperTriangular(A))) else SchurF = schurfact(complex(A)) - R = full(sqrtm(UpperTriangular(SchurF[:T]))) + R = full(sqrt(UpperTriangular(SchurF[:T]))) return SchurF[:vectors] * R * SchurF[:vectors]' end end -function sqrtm(A::StridedMatrix{<:Complex}) +function sqrt(A::StridedMatrix{<:Complex}) if ishermitian(A) - return full(sqrtm(Hermitian(A))) + return full(sqrt(Hermitian(A))) end n = checksquare(A) if istriu(A) - return full(sqrtm(UpperTriangular(A))) + return full(sqrt(UpperTriangular(A))) else SchurF = schurfact(A) - R = full(sqrtm(UpperTriangular(SchurF[:T]))) + R = full(sqrt(UpperTriangular(SchurF[:T]))) return SchurF[:vectors] * R * SchurF[:vectors]' end end -sqrtm(a::Number) = (b = sqrt(complex(a)); imag(b) == 0 ? real(b) : b) -sqrtm(a::Complex) = sqrt(a) function inv(A::StridedMatrix{T}) where T checksquare(A) diff --git a/base/linalg/diagonal.jl b/base/linalg/diagonal.jl index beac1cf844cbd..dfa8408594780 100644 --- a/base/linalg/diagonal.jl +++ b/base/linalg/diagonal.jl @@ -328,8 +328,7 @@ eye(::Type{Diagonal{T}}, n::Int) where {T} = Diagonal(ones(T,n)) # Matrix functions exp(D::Diagonal) = Diagonal(exp.(D.diag)) log(D::Diagonal) = Diagonal(log.(D.diag)) -sqrtm(D::Diagonal) = Diagonal(sqrt.(D.diag)) -sqrtm(D::Diagonal{<:AbstractMatrix}) = Diagonal(sqrtm.(D.diag)) +sqrt(D::Diagonal) = Diagonal(sqrt.(D.diag)) #Linear solver function A_ldiv_B!(D::Diagonal, B::StridedVecOrMat) diff --git a/base/linalg/linalg.jl b/base/linalg/linalg.jl index 46025dd8d2999..09816b9e0a73c 100644 --- a/base/linalg/linalg.jl +++ b/base/linalg/linalg.jl @@ -9,7 +9,7 @@ import Base: USE_BLAS64, abs, big, broadcast, ceil, conj, convert, copy, copy!, adjoint, eltype, exp, eye, findmax, findmin, fill!, floor, full, getindex, hcat, imag, indices, inv, isapprox, isone, IndexStyle, kron, length, log, map, ndims, oneunit, parent, power_by_squaring, print_matrix, promote_rule, real, round, - setindex!, show, similar, size, transpose, trunc, typed_hcat + setindex!, show, similar, size, sqrt, transpose, trunc, typed_hcat using Base: hvcat_fill, iszero, IndexLinear, _length, promote_op, promote_typeof, @propagate_inbounds, @pure, reduce, typed_vcat # We use `_length` because of non-1 indices; releases after julia 0.5 @@ -124,7 +124,6 @@ export schur, schurfact!, schurfact, - sqrtm, svd, svdfact!, svdfact, diff --git a/base/linalg/symmetric.jl b/base/linalg/symmetric.jl index 99d95ba112ab6..adf1c0762f71d 100644 --- a/base/linalg/symmetric.jl +++ b/base/linalg/symmetric.jl @@ -607,40 +607,7 @@ function exp(A::Hermitian{T}) where T end end -for (funm, func) in ([:sqrtm,:sqrt],) - @eval begin - function ($funm)(A::Symmetric{T}) where T<:Real - F = eigfact(A) - if all(λ -> λ ≥ 0, F.values) - retmat = (F.vectors * Diagonal(($func).(F.values))) * F.vectors' - else - retmat = (F.vectors * Diagonal(($func).(complex.(F.values)))) * F.vectors' - end - return Symmetric(retmat) - end - - function ($funm)(A::Hermitian{T}) where T - n = checksquare(A) - F = eigfact(A) - if all(λ -> λ ≥ 0, F.values) - retmat = (F.vectors * Diagonal(($func).(F.values))) * F.vectors' - if T <: Real - return Hermitian(retmat) - else - for i = 1:n - retmat[i,i] = real(retmat[i,i]) - end - return Hermitian(retmat) - end - else - retmat = (F.vectors * Diagonal(($func).(complex(F.values)))) * F.vectors' - return retmat - end - end - end -end - -for func in (:log, #=:sqrtm=#) +for func in (:log, :sqrt) @eval begin function ($func)(A::Symmetric{T}) where T<:Real F = eigfact(A) diff --git a/base/linalg/triangular.jl b/base/linalg/triangular.jl index eb2af9ab869d5..df46c9239931d 100644 --- a/base/linalg/triangular.jl +++ b/base/linalg/triangular.jl @@ -1815,7 +1815,7 @@ function log(A0::UpperTriangular{T}) where T<:Union{Float64,Complex{Float64}} end s0 = s for k = 1:min(s, maxsqrt) - A = sqrtm(A) + A = sqrt(A) end AmI = A - I @@ -1871,7 +1871,7 @@ function log(A0::UpperTriangular{T}) where T<:Union{Float64,Complex{Float64}} m = tmax break end - A = sqrtm(A) + A = sqrt(A) AmI = A - I s = s + 1 end @@ -2015,7 +2015,7 @@ function invsquaring(A0::UpperTriangular, theta) end s0 = s for k = 1:min(s, maxsqrt) - A = sqrtm(A) + A = sqrt(A) end AmI = A - I @@ -2073,7 +2073,7 @@ function invsquaring(A0::UpperTriangular, theta) m = tmax break end - A = sqrtm(A) + A = sqrt(A) AmI = A - I s = s + 1 end @@ -2119,7 +2119,7 @@ unw(x::Number) = ceil((imag(x) - pi) / (2 * pi)) # End of auxiliary functions for matrix logarithm and matrix power -function sqrtm(A::UpperTriangular) +function sqrt(A::UpperTriangular) realmatrix = false if isreal(A) realmatrix = true @@ -2131,9 +2131,9 @@ function sqrtm(A::UpperTriangular) end end end - sqrtm(A,Val(realmatrix)) + sqrt(A,Val(realmatrix)) end -function sqrtm(A::UpperTriangular{T},::Val{realmatrix}) where {T,realmatrix} +function sqrt(A::UpperTriangular{T},::Val{realmatrix}) where {T,realmatrix} B = A.data n = checksquare(B) t = realmatrix ? typeof(sqrt(zero(T))) : typeof(sqrt(complex(zero(T)))) @@ -2151,7 +2151,7 @@ function sqrtm(A::UpperTriangular{T},::Val{realmatrix}) where {T,realmatrix} end return UpperTriangular(R) end -function sqrtm(A::UnitUpperTriangular{T}) where T +function sqrt(A::UnitUpperTriangular{T}) where T B = A.data n = checksquare(B) t = typeof(sqrt(zero(T))) @@ -2169,8 +2169,8 @@ function sqrtm(A::UnitUpperTriangular{T}) where T end return UnitUpperTriangular(R) end -sqrtm(A::LowerTriangular) = sqrtm(A.').' -sqrtm(A::UnitLowerTriangular) = sqrtm(A.').' +sqrt(A::LowerTriangular) = sqrt(A.').' +sqrt(A::UnitLowerTriangular) = sqrt(A.').' # Generic eigensystems eigvals(A::AbstractTriangular) = diag(A) diff --git a/doc/src/manual/linear-algebra.md b/doc/src/manual/linear-algebra.md index 3a0fb0fb4e6cc..7c86c389c0cb4 100644 --- a/doc/src/manual/linear-algebra.md +++ b/doc/src/manual/linear-algebra.md @@ -177,8 +177,8 @@ as well as whether hooks to various optimized methods for them in LAPACK are ava | Matrix type | `+` | `-` | `*` | `\` | Other functions with optimized methods | |:------------------------- |:--- |:--- |:--- |:--- |:------------------------------------------------------------------- | -| [`Symmetric`](@ref) | | | | MV | [`inv()`](@ref), [`sqrtm()`](@ref), [`exp()`](@ref) | -| [`Hermitian`](@ref) | | | | MV | [`inv()`](@ref), [`sqrtm()`](@ref), [`exp()`](@ref) | +| [`Symmetric`](@ref) | | | | MV | [`inv()`](@ref), [`sqrt()`](@ref), [`exp()`](@ref) | +| [`Hermitian`](@ref) | | | | MV | [`inv()`](@ref), [`sqrt()`](@ref), [`exp()`](@ref) | | [`UpperTriangular`](@ref) | | | MV | MV | [`inv()`](@ref), [`det()`](@ref) | | [`LowerTriangular`](@ref) | | | MV | MV | [`inv()`](@ref), [`det()`](@ref) | | [`SymTridiagonal`](@ref) | M | M | MS | MV | [`eigmax()`](@ref), [`eigmin()`](@ref) | diff --git a/doc/src/stdlib/linalg.md b/doc/src/stdlib/linalg.md index 99c8c2f3d3412..126d74954d4d2 100644 --- a/doc/src/stdlib/linalg.md +++ b/doc/src/stdlib/linalg.md @@ -92,7 +92,6 @@ Base.repmat Base.kron Base.SparseArrays.blkdiag Base.LinAlg.linreg -Base.LinAlg.sqrtm Base.LinAlg.lyap Base.LinAlg.sylvester Base.LinAlg.issuccess diff --git a/test/linalg/dense.jl b/test/linalg/dense.jl index ba43ee24345b5..ea313bc8f7579 100644 --- a/test/linalg/dense.jl +++ b/test/linalg/dense.jl @@ -115,10 +115,10 @@ bimg = randn(n,2)/2 end @testset "Matrix square root" begin - asq = sqrtm(a) + asq = sqrt(a) @test asq*asq ≈ a asym = a'+a # symmetric indefinite - asymsq = sqrtm(asym) + asymsq = sqrt(asym) @test asymsq*asymsq ≈ asym end @@ -370,10 +370,10 @@ end @testset "issue #2246" begin A = [1 2 0 0; 0 1 0 0; 0 0 0 0; 0 0 0 0] - Asq = sqrtm(A) + Asq = sqrt(A) @test Asq*Asq ≈ A A2 = view(A, 1:2, 1:2) - A2sq = sqrtm(A2) + A2sq = sqrt(A2) @test A2sq*A2sq ≈ A2 N = 3 @@ -569,12 +569,12 @@ end end for A in (Aa, Ab, Ac, Ad, Ah, ADi) - @test A^(1/2) ≈ sqrtm(A) - @test A^(-1/2) ≈ inv(sqrtm(A)) - @test A^(3/4) ≈ sqrtm(A) * sqrtm(sqrtm(A)) - @test A^(-3/4) ≈ inv(A) * sqrtm(sqrtm(A)) - @test A^(17/8) ≈ A^2 * sqrtm(sqrtm(sqrtm(A))) - @test A^(-17/8) ≈ inv(A^2 * sqrtm(sqrtm(sqrtm(A)))) + @test A^(1/2) ≈ sqrt(A) + @test A^(-1/2) ≈ inv(sqrt(A)) + @test A^(3/4) ≈ sqrt(A) * sqrt(sqrt(A)) + @test A^(-3/4) ≈ inv(A) * sqrt(sqrt(A)) + @test A^(17/8) ≈ A^2 * sqrt(sqrt(sqrt(A))) + @test A^(-17/8) ≈ inv(A^2 * sqrt(sqrt(sqrt(A)))) @test (A^0.2)^5 ≈ A @test (A^(2/3))*(A^(1/3)) ≈ A @test (A^im)^(-im) ≈ A @@ -680,7 +680,6 @@ end @testset "test ops on Numbers for $elty" for elty in [Float32,Float64,Complex64,Complex128] a = rand(elty) @test isposdef(one(elty)) - @test sqrtm(a) == sqrt(a) @test lyap(one(elty),a) == -a/2 end diff --git a/test/linalg/diagonal.jl b/test/linalg/diagonal.jl index 47b7eaf509194..4d490c2d6e281 100644 --- a/test/linalg/diagonal.jl +++ b/test/linalg/diagonal.jl @@ -66,7 +66,7 @@ srand(1) @test log(Diagonal(abs.(D.diag))) ≈ log(abs.(DM)) atol=n^3*eps(relty) end if elty <: BlasComplex - for func in (logdet, sqrtm) + for func in (logdet, sqrt) @test func(D) ≈ func(DM) atol=n^2*eps(relty)*2 end end @@ -383,7 +383,7 @@ end @test exp(D) == Diagonal([exp([1 2; 3 4]), exp([1 2; 3 4])]) @test log(D) == Diagonal([log([1 2; 3 4]), log([1 2; 3 4])]) - @test sqrtm(D) == Diagonal([sqrtm([1 2; 3 4]), sqrtm([1 2; 3 4])]) + @test sqrt(D) == Diagonal([sqrt([1 2; 3 4]), sqrt([1 2; 3 4])]) end @testset "multiplication with Symmetric/Hermitian" begin diff --git a/test/linalg/triangular.jl b/test/linalg/triangular.jl index 8589b530160ff..ff03cbcd2f26c 100644 --- a/test/linalg/triangular.jl +++ b/test/linalg/triangular.jl @@ -237,7 +237,7 @@ for elty1 in (Float32, Float64, BigFloat, Complex64, Complex128, Complex{BigFloa @test ladb ≈ fladb atol=sqrt(eps(real(float(one(elty1)))))*n*n # Matrix square root - @test sqrtm(A1) |> t -> t*t ≈ A1 + @test sqrt(A1) |> t -> t*t ≈ A1 # naivesub errors @test_throws DimensionMismatch naivesub!(A1,ones(elty1,n+1)) @@ -391,10 +391,10 @@ end # Matrix square root Atn = UpperTriangular([-1 1 2; 0 -2 2; 0 0 -3]) Atp = UpperTriangular([1 1 2; 0 2 2; 0 0 3]) -@test sqrtm(Atn) |> t->t*t ≈ Atn -@test typeof(sqrtm(Atn)[1,1]) <: Complex -@test sqrtm(Atp) |> t->t*t ≈ Atp -@test typeof(sqrtm(Atp)[1,1]) <: Real +@test sqrt(Atn) |> t->t*t ≈ Atn +@test typeof(sqrt(Atn)[1,1]) <: Complex +@test sqrt(Atp) |> t->t*t ≈ Atp +@test typeof(sqrt(Atp)[1,1]) <: Real Areal = randn(n, n)/2 Aimg = randn(n, n)/2 @@ -511,5 +511,5 @@ end isdefined(Main, :TestHelpers) || @eval Main include("../TestHelpers.jl") using TestHelpers.Furlong let A = UpperTriangular([Furlong(1) Furlong(4); Furlong(0) Furlong(1)]) - @test sqrtm(A) == Furlong{1//2}.(UpperTriangular([1 2; 0 1])) + @test sqrt(A) == Furlong{1//2}.(UpperTriangular([1 2; 0 1])) end