# Defining Sparse Matrix from Scratch using vectors

Hi

I want to create a custom tridiagonal matrix type whose inputs will be 3 different vectors (`c`, `d` and `e`). These 3 vectors will make the matrix in a format as follows: I tried writing a struct

``````struct tridiag <: AbstractMatrix{<:AbstractFloat}
c::Array{<:AbstractFloat}
d::Array{<:AbstractFloat}
e::Array{<:AbstractFloat}
end
``````

but I’m getting an error

``````invalid subtyping in definition of tridiag

Stacktrace:
 top-level scope at In:1
``````

I’m also not sure how I can define a custom type for tridiagonal (sprase) matrix effectively like this. Any other ideas are welcome and I want to writing this type from scratch and not use something like `Sparce matrix`.

Will `Tridiagonal` do what you want? You may have to do `using LinearAlgebra` to get it, but it’s built in.

``````julia> a=[1,1,1]; b=[2,2,2,2]; c=[3,3,3];

julia> Tridiagonal(a,b,c)
4×4 Tridiagonal{Int64,Array{Int64,1}}:
2  3  ⋅  ⋅
1  2  3  ⋅
⋅  1  2  3
⋅  ⋅  1  2
``````
2 Likes

I’m just starting out with Julia and learning by writing basic numerical methods from scratch without any libraries. I want to do something like this but want to apply from scratch. I think I can create something similar using struct. Maybe I can read the source code for this linear algebra implementation of tridiagonal. Can you suggest me something here?

Here are the docs Linear Algebra · The Julia Language and there’s a link to the source code in the lower right corner (you might need to select or hover for it to show up).

1 Like

You want the syntax:

``````struct Tridiag{T<:Number} <: AbstractMatrix{T}
c::Vector{T}
d::Vector{T}
e::Vector{T}
end
``````

Note that this defines a whole set of types `Tridiag{Float64}`, `Tridiag{Int}`, and so on, where for each of which the fields are concrete types. (Note also that `Array{T}` is abstract because you didn’t specify the dimensionality — you want `Array{T,1}` or equivalently `Vector{T}`.)

5 Likes

You can also see the Alan Edelman’s video on Matrix Structures for inspiration.

I have written (to learn and for fun) the structure as suggested by Steven and some two other functions, and already matrix multiplications work. This is very nice ``````julia> A = Tridiag([1,1],[2,2,2],[3,3])
3×3 Tridiag{Int64}:
2  3  0
1  2  3
0  1  2

julia> A*A
3×3 Array{Int64,2}:
7  12   9
4  10  12
1   4   7

``````

If you want to see what I did so far, here is the code:

Code
``````julia> struct Tridiag{T<:Number} <: AbstractMatrix{T}
a :: Vector{T}
b :: Vector{T}
c :: Vector{T}
end

julia> function Base.getindex(M::Tridiag{T}, i::Int, j::Int) where T
if i == j
M.b[i]
elseif i == j + 1
M.a[j]
elseif i == j - 1
M.c[i]
else
zero(T)
end
end

julia> Base.size(M::Tridiag) = (length(M.b),length(M.b))

julia> A = Tridiag([1,1],[2,2,2],[3,3])
3×3 Tridiag{Int64}:
2  3  0
1  2  3
0  1  2

julia> A*A
3×3 Array{Int64,2}:
7  12   9
4  10  12
1   4   7
``````
4 Likes

@lmiq, thanks for sharing that nice video link and your code too. If one has missed the previous lessons in the series, one might not know that some of the examples shown require using several packages (ex: LinearAlgebra, SparseArrays, Statistics).

1 Like