Inverse of a view of a symmetric matrix

Specific Domains Numerics
1

Hi,

I am trying to compute the inverse of a view of a symmetric matrix via:

using LinearAlgebra
V=Symmetric(rand(10,10)+50*I);
V_view = @view V[1:5,1:5];
V_view_inv = inv(V_view);

However, it does not work and I get the following error message:

ERROR: MethodError: no method matching factorize(::SubArray{Float64,2,Symmetric{Float64,Array{Float64,2}},Tuple{UnitRange{Int64},UnitRange{Int64}},false})

I wonder whether I should use a different (specialised) inverse function to perform this operation. What would you suggest to do?

2

You could use V_view = Symmetric(@view parent(V)[1:5, 1:5]), so that the Symmetric wrapper is on the outside rather than on the inside. There should probably be a fallback that allows you to do what you tried, but I think I think this is better anyway, since what you’re doing causes there to be a loss of compile time information regarding the structure of the matrix: a view of a symmetric matrix is not necessarily symmetric (for arbitrary UnitRanges), so Julia can’t dispatch to a specialized method for computing the inverse of a symmetric matrix.

2 Likes
3

Thank you.

I am not familiar with the parent(.) syntax. Would you mind giving me more details on that? Is there any difference between V.data and parent(V)?

4

It’s the same:

https://github.com/JuliaLang/julia/blob/8d4f6d24c0fa87509fcf40788ca2c83cdfecd9af/stdlib/LinearAlgebra/src/symmetric.jl#L270

But directly accessing a field of the Symmetric type could be seen as relying on an implementation detail. It’s not likely that the data field will be renamed in the future, but if there’s a documented function you can use to achieve the same, it’s generally preferable to use it.

3 Likes