|
| 1 | +""" |
| 2 | +Wilkinson Matrix |
| 3 | +================ |
| 4 | +The Wilkinson matrix is a symmetric tridiagonal matrix with pairs |
| 5 | +of nearly equal eigenvalues. The most frequently used case |
| 6 | +is `matrixdepot("wilkinson", 21)`. The result is of type `Tridiagonal`. |
| 7 | +
|
| 8 | +*Input options:* |
| 9 | +
|
| 10 | ++ dim: the dimension of the matrix. |
| 11 | +
|
| 12 | +*References:* |
| 13 | +
|
| 14 | +**J. H. Wilkinson**, Error analysis of direct methods |
| 15 | +of matrix inversion, J. Assoc. Comput. Mach., 8 (1961), pp. 281-330. |
| 16 | +""" |
| 17 | +struct Wilkinson{T<:Number} <: AbstractMatrix{T} |
| 18 | + n::Integer |
| 19 | + |
| 20 | + function Wilkinson{T}(n::Integer) where {T<:Number} |
| 21 | + n >= 0 || throw(ArgumentError("$n < 0")) |
| 22 | + return new{T}(n) |
| 23 | + end |
| 24 | +end |
| 25 | + |
| 26 | +# constructors |
| 27 | +Wilkinson(n::Integer) = Wilkinson{Float64}(n) |
| 28 | + |
| 29 | +# metadata |
| 30 | +@properties Wilkinson [:symmetric, :eigen] |
| 31 | + |
| 32 | +# properties |
| 33 | +size(A::Wilkinson) = (A.n, A.n) |
| 34 | +LinearAlgebra.isdiag(A::Wilkinson) = A.n <= 1 ? true : false |
| 35 | +LinearAlgebra.issymmetric(::Wilkinson) = true |
| 36 | + |
| 37 | +# functions |
| 38 | +@inline Base.@propagate_inbounds function getindex(A::Wilkinson{T}, i::Integer, j::Integer) where {T} |
| 39 | + @boundscheck checkbounds(A, i, j) |
| 40 | + n = A.n |
| 41 | + if n == 1 |
| 42 | + return one(T) |
| 43 | + elseif i == j |
| 44 | + m = (n - 1) / 2 |
| 45 | + return T(abs(m - (i - 1))) |
| 46 | + elseif i == j + 1 || i == j - 1 |
| 47 | + return one(T) |
| 48 | + else |
| 49 | + return zero(T) |
| 50 | + end |
| 51 | +end |
| 52 | + |
| 53 | +function Base.replace_in_print_matrix(A::Wilkinson, i::Integer, j::Integer, s::AbstractString) |
| 54 | + i == j - 1 || i == j || i == j + 1 ? s : Base.replace_with_centered_mark(s) |
| 55 | +end |
0 commit comments