Skip to content

Commit 4ca74fe

Browse files
ssikdar1KristofferC
authored andcommitted
Have matrix square root account for Rii=Rjj=0 to prevent NaN values (#35758)
(cherry picked from commit 880c731)
1 parent d6bdf50 commit 4ca74fe

File tree

2 files changed

+8
-1
lines changed

2 files changed

+8
-1
lines changed

stdlib/LinearAlgebra/src/triangular.jl

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2607,7 +2607,9 @@ function sqrt(A::UpperTriangular{T},::Val{realmatrix}) where {T,realmatrix}
26072607
@simd for k = i+1:j-1
26082608
r -= R[i,k]*R[k,j]
26092609
end
2610-
iszero(r) || (R[i,j] = sylvester(R[i,i],R[j,j],-r))
2610+
if !(iszero(r) || (iszero(R[i,i]) && iszero(R[j,j])))
2611+
R[i,j] = sylvester(R[i,i],R[j,j],-r)
2612+
end
26112613
end
26122614
end
26132615
return UpperTriangular(R)

stdlib/LinearAlgebra/test/triangular.jl

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -651,4 +651,9 @@ end
651651
@test_throws ArgumentError LinearAlgebra.powm(A, 2.2)
652652
end
653653

654+
# Issue 35058
655+
let A = [0.9999999999999998 4.649058915617843e-16 -1.3149405273715513e-16 9.9959579317056e-17; -8.326672684688674e-16 1.0000000000000004 2.9280733590254494e-16 -2.9993900031619594e-16; 9.43689570931383e-16 -1.339206523454095e-15 1.0000000000000007 -8.550505126287743e-16; -6.245004513516506e-16 -2.0122792321330962e-16 1.183061278035052e-16 1.0000000000000002],
656+
B = [0.09648289218436859 0.023497875751503007 0.0 0.0; 0.023497875751503007 0.045787575150300804 0.0 0.0; 0.0 0.0 0.0 0.0; 0.0 0.0 0.0 0.0]
657+
@test sqrt(A*B*A')^2 A*B*A'
658+
end
654659
end # module TestTriangular

0 commit comments

Comments
 (0)