I am a bit curious about the question in the title. Since it is unanswered in the docstring, I list its definition as a reference.
# Symmetric and Hermitian matrices
struct Symmetric{T,S<:AbstractMatrix{<:T}} <: AbstractMatrix{T}
data::S
uplo::Char
function Symmetric{T,S}(data, uplo::Char) where {T,S<:AbstractMatrix{<:T}}
require_one_based_indexing(data)
(uplo != 'U' && uplo != 'L') && throw_uplo()
new{T,S}(data, uplo)
end
end
Since Symmetric is merely a “shell”, I wonder what is the purpose of T (the design that singles it out). Consider the following example
julia> using LinearAlgebra
julia> function m() return rand(-9:.1:9, 2, 2) end
m (generic function with 1 method)
julia> A = m()
2Ă—2 Matrix{Float64}:
1.9 3.7
8.7 -1.8
julia> LinearAlgebra.symmetric_type(typeof(A[1, 1]))
Float64
julia> LinearAlgebra.symmetric_type(typeof(A))
Symmetric{Float64, Matrix{Float64}}
julia> B = [m() for i in 1:2, j in 1:2]
2Ă—2 Matrix{Matrix{Float64}}:
[8.6 2.1; 0.8 8.9] [0.7 -8.4; -2.3 -2.8]
[6.4 4.1; -5.8 1.7] [4.6 -3.0; 4.5 5.6]
julia> LinearAlgebra.symmetric_type(typeof(B[1, 1]))
Symmetric{Float64, Matrix{Float64}}
julia> LinearAlgebra.symmetric_type(typeof(B))
Symmetric{AbstractMatrix, Matrix{Matrix{Float64}}}
See the last line, AbstractMatrix is hardly instructive.
I mean, what about
The related code is here
(and line 60 therein).
T here is the type of the elements inside of the Symmetric, indicated by <: AbstractMatrix{T} (where that same T designates the element type).
Are you perhaps looking for eltype?
1 Like
Yes. But then, isn’t this T redundant? making the type expression looks cumbersome. (I’ve revised the title of the topic). Is the T indispensable? We could have infer/dispatch conveniently, e.g. with the aid of eltype you suggested.
You need the T so you can write <: AbstractMatrix{T} in the definition. You can’t just do struct Symmetric{S} or struct Symmetric{S<:AbstractMatrix{T}} (since this last one needs T to be defined first).
3 Likes
Why not do this
julia> struct MySymmetric{S <: AbstractMatrix{T} where T <: Union{Number, AbstractMatrix}}
data::S
end
julia> using LinearAlgebra
julia> A = rand(2, 2); B = [rand(2, 2) for i in 1:2, j in 1:2];
julia> MySymmetric(A)
MySymmetric{Matrix{Float64}}([0.26880897363440703 0.6476423313174478; 0.6105024571886544 0.1442311885018166])
julia> LinearAlgebra.Symmetric(A)
2Ă—2 Symmetric{Float64, Matrix{Float64}}:
0.268809 0.647642
0.647642 0.144231
julia> MySymmetric(B)
MySymmetric{Matrix{Matrix{Float64}}}([[0.8054973493291133 0.7597395285524011; 0.14616471160829714 0.42034533671996954] [0.737336878028752 0.8408610455649348; 0.48620788006404236 0.15674775620016002]; [0.8856275147903988 0.7373043349918851; 0.3587458066422289 0.32907006075331047] [0.8815327773307925 0.8533663448400511; 0.15092667149027095 0.18201636320297676]])
julia> LinearAlgebra.Symmetric(B)
2Ă—2 Symmetric{AbstractMatrix, Matrix{Matrix{Float64}}}:
[0.805497 0.75974; 0.75974 0.420345] [0.737337 0.840861; 0.486208 0.156748]
[0.737337 0.486208; 0.840861 0.156748] [0.881533 0.853366; 0.853366 0.182016]
julia> typeof(MySymmetric(A)) <: AbstractMatrix{Float64}
false
julia> typeof(Symmetric(A)) <: AbstractMatrix{Float64}
true
So you want to declare
struct MySymmetric{S <: AbstractMatrix{T}} <: AbstractMatrix{T} end
but that doesn’t work because T is not defined. You can fix it by adding the extra type parameter
struct MySymmetric{T, S <: AbstractMatrix{T}} <: AbstractMatrix{T} end
4 Likes
You can define that with a local T
struct MySymmetric{S<:AbstractMatrix} end
# or equivalently
struct MySymmetric{S<:AbstractMatrix{T where T}} end
but the main problem is then still that MySymmetric{Matrix{Float64} is never a subtype of AbstractMatrix{Float64} as one can see in the first typeof in the post above me. And that is something one wants for nice default dispatches for all eltypes.
3 Likes
Yes, your last line is good. I copy-paste it here
struct MySymmetric{T, S <: AbstractMatrix{T}} <: AbstractMatrix{T} end
But the LinearAlgebra.jl define it as
struct Symmetric{T, S <: AbstractMatrix{<:T}} <: AbstractMatrix{T} end
We see that the latter uses <: T, what is the notivation? so that
Looks a bit inconsistent.
I now understand this point, thank you.
I’m unsure if this is intentional, probably just an oversight that has not been corrected since.
There is a discussion in the relevant PR but it did not clarify the issue for me, I don’t understand the possibility they mention.
1 Like
My understanding is that it’s because of matrices of matrices (transpose of a number has the same type, but transpose of a matrix is of type Transpose).
Example behavior:
julia> A = [1 2;
3 4]
2Ă—2 Matrix{Int64}:
1 2
3 4
julia> B = [[A] [A];
[A] [A]]
2Ă—2 Matrix{Matrix{Int64}}:
[1 2; 3 4] [1 2; 3 4]
[1 2; 3 4] [1 2; 3 4]
julia> C = Symmetric(B)
2Ă—2 Symmetric{AbstractMatrix, Matrix{Matrix{Int64}}}:
[1 2; 2 4] [1 2; 3 4]
[1 3; 2 4] [1 2; 2 4]
julia> C[1, 1]
2Ă—2 Symmetric{Int64, Matrix{Int64}}:
1 2
2 4
julia> C[1, 2]
2Ă—2 Matrix{Int64}:
1 2
3 4
julia> C[2, 1]
2Ă—2 transpose(::Matrix{Int64}) with eltype Int64:
1 3
2 4
(three different eltypes, all of which are <: AbstractMatrix)
2 Likes
Do make sense. But in this case, T is not some “eltype” but a “least common eltype”. And in your example, T is AbstractMatrix, indicating that the dispatch advantage is lost. (But I guess this usage is not practical)
Wait a minute. There is still a question.
The LinearAlgebra.jl defines
struct Symmetric{T, S <: AbstractMatrix{<:T}} <: AbstractMatrix{T} end
And you call C = Symmetric(B), where B isa Matrix{Matrix{Int64}}.
At the time of construction (which is “earlier”), T can surely be Matrix{Int64}.
It’s after some manipulations (at a “later” time) that julia find AbstractMatrix can serve as a “least common supertype” for all the eltypes.
The weird thing is: Did I tell julia that what role I want T to play? Did I stipulate that T should be the “least common supertype” at that “later” time? Seems I hadn’t.
Then why can’t T also be Matrix{Int64}?
