Skip to content

Commit 8b41976

Browse files
lkdvosclaude
andcommitted
fix: correct the out-of-domain DomainError messages
The `lazy"..."` literal introduced for the negative-eigenvalue error does not process the `\` line continuation, so the rendered message contained a literal backslash, newline and indentation. Use ordinary string concatenation, as elsewhere in this file, and report the tolerance under its user-facing name `domain_atol` rather than `atol`. Restore the distinction between the real ("a negative real eigenvalue") and complex ("an eigenvalue on the negative real axis") variants by passing the description to the shared helper, and give `_check_nonzero_eigenvalues` the same `@noinline` throw helper so both domain checks stay free of error-path code on the GPU. Thread `domain_atol` back into the diagonal kernels invoked by the `squareroot` eigendecomposition paths, and let the `DiagonalAlgorithm` kernel reuse `_clamp_domain_eigenvalues!` again so it matches `power!`. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
1 parent f009a95 commit 8b41976

3 files changed

Lines changed: 30 additions & 38 deletions

File tree

src/implementations/matrixfunctions.jl

Lines changed: 24 additions & 26 deletions
Original file line numberDiff line numberDiff line change
@@ -1,20 +1,33 @@
11
# Shared helpers for matrix functions with a restricted domain
22
# -------------------------------------------------------------
33

4+
# The throwing branches live in `@noinline` helpers so that the reductions and broadcasts
5+
# below stay free of error-path code, which keeps them GPU friendly.
6+
@noinline function throw_negative_eigenvalue(λmin, atol, what)
7+
return throw(
8+
DomainError(
9+
λmin,
10+
"The matrix has $what beyond `domain_atol = $atol` and the result of this matrix function is complex. " *
11+
"Pass a complex matrix to obtain the principal value, or increase `domain_atol` if the eigenvalue is a rounding artifact."
12+
)
13+
)
14+
end
15+
16+
@noinline function throw_zero_eigenvalue(amin, atol)
17+
return throw(
18+
DomainError(
19+
amin,
20+
"The matrix has a (numerically) zero eigenvalue within `domain_atol = $atol`, for which this matrix function is not defined."
21+
)
22+
)
23+
end
24+
425
# Clamp real eigenvalues that are negative within `atol` (rounding artifacts) to zero,
526
# and throw a `DomainError` for eigenvalues that are genuinely negative, since then the
627
# result cannot be expressed with the same (real) scalar type.
728
function _clamp_domain_eigenvalues!::AbstractVector{<:Real}, atol::Real)
829
λmin = minimum(λ; init = zero(eltype(λ)))
9-
if λmin < -atol
10-
throw(
11-
DomainError(
12-
λmin,
13-
"The matrix has a negative real eigenvalue beyond `domain_atol = $atol` and the result of this matrix function is complex. " *
14-
"Pass a complex matrix to obtain the principal value, or increase `domain_atol` if the eigenvalue is a rounding artifact."
15-
)
16-
)
17-
end
30+
λmin < -atol && throw_negative_eigenvalue(λmin, atol, "a negative real eigenvalue")
1831
λ .= max.(λ, zero(eltype(λ)))
1932
return λ
2033
end
@@ -32,15 +45,7 @@ end
3245
function _clamp_domain_eigenvalues!::AbstractVector{<:Complex}, atol::Real)
3346
onaxis = x -> abs(imag(x)) <= atol && real(x) < 0
3447
λmin = mapreduce(x -> onaxis(x) ? real(x) : zero(real(x)), min, λ; init = zero(real(eltype(λ))))
35-
if λmin < -atol
36-
throw(
37-
DomainError(
38-
λmin,
39-
"The matrix has an eigenvalue on the negative real axis beyond `domain_atol = $atol` and the result of this matrix function is complex. " *
40-
"Pass a complex matrix to obtain the principal value, or increase `domain_atol` if the eigenvalue is a rounding artifact."
41-
)
42-
)
43-
end
48+
λmin < -atol && throw_negative_eigenvalue(λmin, atol, "an eigenvalue on the negative real axis")
4449
λ .= ifelse.(onaxis.(λ), zero(eltype(λ)), λ)
4550
return λ
4651
end
@@ -49,14 +54,7 @@ end
4954
# e.g. `logarithm` and `power` with a negative fractional power.
5055
function _check_nonzero_eigenvalues(λ, atol::Real)
5156
amin = minimum(abs, λ; init = typemax(real(eltype(λ))))
52-
if amin <= atol
53-
throw(
54-
DomainError(
55-
amin,
56-
"The matrix has a (numerically) zero eigenvalue within `domain_atol = $atol`, for which this matrix function is not defined."
57-
)
58-
)
59-
end
57+
amin <= atol && throw_zero_eigenvalue(amin, atol)
6058
return λ
6159
end
6260

src/implementations/squareroot.jl

Lines changed: 2 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -38,9 +38,8 @@ function squareroot!(A::AbstractMatrix, sqrtA, alg::MatrixFunctionViaLA)
3838
# `LinearAlgebra.sqrt` of a real matrix is real whenever the principal square root is,
3939
# so a complex result with a real output signals a genuine domain violation
4040
sqrtAc = LinearAlgebra.sqrt(A)
41-
if eltype(sqrtAc) <: Complex && !(eltype(sqrtA) <: Complex)
41+
eltype(sqrtAc) <: Complex && !(eltype(sqrtA) <: Complex) &&
4242
throw(_realness_domainerror(squareroot!))
43-
end
4443
copy!(sqrtA, sqrtAc)
4544
return sqrtA
4645
end
@@ -50,8 +49,7 @@ function squareroot!(A::AbstractMatrix, sqrtA, alg::MatrixFunctionViaEigh)
5049
D, V = eigh_full!(A, alg.eigh_alg)
5150
diag_alg = DiagonalAlgorithm(; domain_atol = alg.domain_atol)
5251
# `sqrt(A) = (V * D^(1/4)) * (V * D^(1/4))'` is hermitian by construction
53-
sqrtD = squareroot!(D, D, diag_alg)
54-
Vs = rmul!(V, squareroot!(sqrtD, sqrtD, diag_alg))
52+
Vs = rmul!(V, power!(D, 1 // 4, D, diag_alg))
5553
return _mul_herm!(sqrtA, Vs)
5654
end
5755

src/interface/logarithm.jl

Lines changed: 4 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -11,18 +11,14 @@ Compute the principal logarithm `logA` of the square matrix `A`, i.e. the logari
1111
whose eigenvalues have imaginary part in `(-π, π]`.
1212
1313
The scalar type of the output matches that of the input.
14-
As a consequence, a real matrix with eigenvalues on the negative real axis, for which
15-
the principal logarithm is complex, leads to a `DomainError`; pass a complex matrix
16-
to obtain the principal value.
17-
A matrix with (numerically) zero eigenvalues has no logarithm and also leads to a
18-
`DomainError`.
14+
As a consequence, a real matrix with eigenvalues on the negative real axis, for which the principal logarithm is complex, leads to a `DomainError`; pass a complex matrix to obtain the principal value.
15+
A matrix with (numerically) zero eigenvalues has no logarithm and also leads to a `DomainError`.
1916
Both checks use a tolerance `domain_atol`, which can be specified for the algorithms
2017
that support it and defaults to [`default_domain_atol`](@ref).
2118
2219
!!! note
23-
The bang method `logarithm!` optionally accepts the output structure and
24-
possibly destroys the input matrix `A`. Always use the return value of the function
25-
as it may not always be possible to use the provided `logA` as output.
20+
The bang method `logarithm!` optionally accepts the output structure and possibly destroys the input matrix `A`.
21+
Always use the return value of the function as it may not always be possible to use the provided `logA` as output.
2622
"""
2723
@functiondef logarithm
2824

0 commit comments

Comments
 (0)