Error in Sqrt for julia > 1.10.2

[CodeLenz]

Reference thread

https://discourse.julialang.org/t/question-regarding-numerical-stability-among-julia-versions/112741

sqrt(A) is failing when A is non-symmetric with repeated eigenvalues.

[CodeLenz]

The culprit seems to be schur(A), computed inside sqrt(A).

[dkarrasch]

@aravindh-krishnamoorthy You may be interested in this. I saw you're working on efficient square roots of matrices.

[aravindh-krishnamoorthy]

Interesting problem! I see from the Discourse thread that @CodeLenz has studied it quite a bit. I'll take a closer look at the issue.

[CodeLenz]

@aravindh-krishnamoorthy

I tried the code in your repository, but the same pattern remains. I managed to reduce the test Matrix to

function Test_sqrt(N=1000)

    # Initialize matrix A
    A = zeros(N,N)

    # A is tridiagonal
    E = 210E9
    p = 7.3E-4
    l = 20.0/N

    # Main diagonal
    cte = E/(p*l^2)
    A[diagind(A,0)].= 2.0*cte

    # Upper and lower first diagonals
    cte = -E/(p*l^2)
    A[diagind(A,1)].= cte
    A[diagind(A,-1)].= cte

    # Sqrt
    sA = sqrt(A)

    # Test sqrt
    return A, norm(sA*sA .- A) 

end

[aravindh-krishnamoorthy]

@aravindh-krishnamoorthy

I tried the code in your repository, but the same pattern remains. I managed to reduce the test Matrix to

[snip|

Hello @CodeLenz. I just started with your function Test_sqrt provided above. But, as seen from below, there's no difference in the error values produced between versions 1.9.4 and 1.11.0-alpha2 (which is > 1.10.2) on my PC. Could you please help me identify the exact problem to be solved? Do you think it's better for me to look at the original problem reported in the Discourse thread with the square root of A = -4*Kb having an issue between these versions?

julia> VERSION
v"1.9.4"
julia> using LinearAlgebra
julia> include("test_sqrt.jl")
Test_sqrt (generic function with 2 methods)
julia> A, e = Test_sqrt()
([1.4383561643835615e18 -7.191780821917807e17 … 0.0 0.0; -7.191780821917807e17 1.4383561643835615e18 … 0.0 0.0; … ; 0.0 0.0 … 1.4383561643835615e18 -7.191780821917807e17; 0.0 0.0 … -7.191780821917807e17 1.4383561643835615e18], 5.647486380061708e6)

julia> VERSION
v"1.11.0-alpha2"
julia> using LinearAlgebra
julia> include("test_sqrt.jl")
Test_sqrt (generic function with 2 methods)
julia> A, e = Test_sqrt()
([1.4383561643835615e18 -7.191780821917807e17 … 0.0 0.0; -7.191780821917807e17 1.4383561643835615e18 … 0.0 0.0; … ; 0.0 0.0 … 1.4383561643835615e18 -7.191780821917807e17; 0.0 0.0 … -7.191780821917807e17 1.4383561643835615e18], 5.647486380061708e6)