Cauchy{T,X,Y} <: AbstractMatrix{T}
Cauchy(x [,y])

Construct lazy Cauchy matrix where

• Cauchy(x,y)[i,j] = 1 / (x[i] + y[j])
• Cauchy(x) = Cauchy(x,x)
• Cauchy(k::Int) = Cauchy(1:k)

Both x and y can be any AbstractArray or NTuple of Numbers, but all elements of x must have the same type; likewise for y. The elements of x and of y should be distinct.

The element type T corresponds to the reciprocal of the sum of pairs of elements of x and y.

julia> Cauchy([2.0 1], (0, 1, 2))
2×3 Cauchy{Float64, Matrix{Float64}, Tuple{Int64, Int64, Int64}}:
 0.5  0.333333  0.25
 1.0  0.5       0.333333

julia> Cauchy(1:3)
3×3 Cauchy{Rational{Int64}, UnitRange{Int64}, UnitRange{Int64}}:
 1//2  1//3  1//4
 1//3  1//4  1//5
 1//4  1//5  1//6

julia> Cauchy(3)
3×3 Cauchy{Rational{Int64}, UnitRange{Int64}, UnitRange{Int64}}:
 1//2  1//3  1//4
 1//3  1//4  1//5
 1//4  1//5  1//6

Companion(v::Union{AbstractVector,Polynomial})::AbstractMatrix

Construct a lazy companion matrix from the coefficients of its characteristic (monic) polynomial.

The matrix is n × n for a vector input of length n or an input Polynomial of degree n, but the storage here is only O(n). This version puts the coefficients along the last column of the matrix. Some texts put the coefficients along the first row, the transpose of the convention used here.

This type has efficient methods for mul! and inv.

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

Also, directly from a polynomial:

julia> using Polynomials

julia> P = Polynomial(2:5)
Polynomials.Polynomial(2 + 3*x + 4*x^2 + 5*x^3)

julia> C = Companion(P)
3×3 Companion{Float64}:
 0.0  0.0  -0.4
 1.0  0.0  -0.6
 0.0  1.0  -0.8

[`Frobenius` matrix](http://en.wikipedia.org/wiki/Frobenius_matrix)

Frobenius matrices or Gaussian elimination matrices of the form

[ 1 0 ...     0 ]
[ 0 1 ...     0 ]
[ .........     ]
[ ... 1 ...     ]
[ ... c1 1 ...  ]
[ ... c2 0 1 ...]
[ ............. ]
[ ... ck ...   1]

i.e., an identity matrix with possibly nonzero subdiagonal elements along a single column.

In this implementation, the subdiagonal of the nonzero column is stored as a dense vector, so that the size can be inferred automatically as j+k, where j is the column index and k is the number of subdiagonal elements.

julia> using SpecialMatrices

julia> F = Frobenius(3, 4:6) # Specify subdiagonals of column 3
6×6 Frobenius{Int64}:
 1  0  0  0  0  0
 0  1  0  0  0  0
 0  0  1  0  0  0
 0  0  4  1  0  0
 0  0  5  0  1  0
 0  0  6  0  0  1

julia> inv(F) # Special form of inverse
6×6 Frobenius{Int64}:
 1  0   0  0  0  0
 0  1   0  0  0  0
 0  0   1  0  0  0
 0  0  -4  1  0  0
 0  0  -5  0  1  0
 0  0  -6  0  0  1

julia> F*F # Special form preserved if the same column has the subdiagonals
6×6 Frobenius{Int64}:
 1  0   0  0  0  0
 0  1   0  0  0  0
 0  0   1  0  0  0
 0  0   8  1  0  0
 0  0  10  0  1  0
 0  0  12  0  0  1

julia> F*Frobenius(2, 2:5) # Promotes to Matrix
6×6 Matrix{Int64}:
 1   0  0  0  0  0
 0   1  0  0  0  0
 0   2  1  0  0  0
 0  11  4  1  0  0
 0  14  5  0  1  0
 0  17  6  0  0  1

julia> F * [10.0,20,30,40,50,60.0]
6-element Vector{Float64}:
  10.0
  20.0
  30.0
 160.0
 200.0
 240.0

Hilbert(m [,n])

Construct m × m or m × n Hilbert matrix from its specified dimensions, where element i,j equal to 1 / (i+j-1).

julia> A = Hilbert(5)
5×5 Hilbert{Rational{Int64}}:
 1//1  1//2  1//3  1//4  1//5
 1//2  1//3  1//4  1//5  1//6
 1//3  1//4  1//5  1//6  1//7
 1//4  1//5  1//6  1//7  1//8
 1//5  1//6  1//7  1//8  1//9

julia> Matrix(A)
5×5 Matrix{Rational{Int64}}:
 1//1  1//2  1//3  1//4  1//5
 1//2  1//3  1//4  1//5  1//6
 1//3  1//4  1//5  1//6  1//7
 1//4  1//5  1//6  1//7  1//8
 1//5  1//6  1//7  1//8  1//9

Inverses are also integer matrices:

julia> inv(A)
5×5 InverseHilbert{Rational{Int64}}:
    25//1    -300//1     1050//1    -1400//1     630//1
  -300//1    4800//1   -18900//1    26880//1  -12600//1
  1050//1  -18900//1    79380//1  -117600//1   56700//1
 -1400//1   26880//1  -117600//1   179200//1  -88200//1
   630//1  -12600//1    56700//1   -88200//1   44100//1

Kahan matrix

julia> A = Kahan(5,5,1,35)
5×5 Kahan{Int64, Int64}:
 1.0  -0.540302  -0.540302  -0.540302  -0.540302
 0.0   0.841471  -0.454649  -0.454649  -0.454649
 0.0   0.0        0.708073  -0.382574  -0.382574
 0.0   0.0        0.0        0.595823  -0.321925
 0.0   0.0        0.0        0.0        0.501368

julia> A = Kahan(5,3,0.5,0)
5×3 Kahan{Float64, Int64}:
 1.0  -0.877583  -0.877583
 0.0   0.479426  -0.420735
 0.0   0.0        0.229849
 0.0   0.0        0.0
 0.0   0.0        0.0

julia> A = Kahan(3,5,0.5,1e-3)
3×5 Kahan{Float64, Float64}:
 1.0  -0.877583  -0.877583  -0.877583  -0.877583
 0.0   0.479426  -0.420735  -0.420735  -0.420735
 0.0   0.0        0.229849  -0.201711  -0.201711

For more details see: N. Higham, "A Survey of Condition Number Estimation for Triangular Matrices," SIMAX, Vol. 29, No. 4, pp. 588, 1987, https://doi.org/10.1137/1029112

Riemann(N::Int)

Construct N × N Riemann matrix, defined as A = B[2:N+1, 2:N+1], where B[i,j] = i-1 if i divides j, and -1 otherwise. The Riemann hypothesis holds if and only if det(A) = O(N! N^(-1/2+ϵ)) for every ϵ > 0.

julia> Riemann(7)
7×7 Riemann{Int64}:
  1  -1   1  -1   1  -1   1
 -1   2  -1  -1   2  -1  -1
 -1  -1   3  -1  -1  -1   3
 -1  -1  -1   4  -1  -1  -1
 -1  -1  -1  -1   5  -1  -1
 -1  -1  -1  -1  -1   6  -1
 -1  -1  -1  -1  -1  -1   7

For more details see Friedrich Roesler, "Riemann's hypothesis as an eigenvalue problem," Linear Algebra and its Applications, Vol. 81, p.153-198, Sep. 1986. https://doi.org/10.1016/0024-3795(86)90255-7

Strang([T=Int,] n::Int)

Construct Strang matrix with elements of type T. This special matrix named after Gilbert Strang is symmetric, tridiagonal, and Toeplitz. See Fig. 1 of Julia paper.

julia> Strang(4)
4×4 Strang{Int64}:
  2  -1   0   0
 -1   2  -1   0
  0  -1   2  -1
  0   0  -1   2

julia> Strang(Int16, 3)
3×3 Strang{Int16}:
  2  -1   0
 -1   2  -1
  0  -1   2

Vandermonde(c::AbstractVector)

Create a "lazy" n × n Vandermonde matrix where $A_{ij} = c_i^{j-1}$, requiring only O(n) storage for the vector c.

The transpose and adjoint operations are also lazy.

julia> a = 1:5; A = Vandermonde(a)
5×5 Vandermonde{Int64}:
 1  1   1    1    1
 1  2   4    8   16
 1  3   9   27   81
 1  4  16   64  256
 1  5  25  125  625

julia> A'
5×5 adjoint(::Vandermonde{Int64}) with eltype Int64:
 1   1   1    1    1
 1   2   3    4    5
 1   4   9   16   25
 1   8  27   64  125
 1  16  81  256  625

The backslash operator \ is overloaded to solve Vandermonde and adjoint Vandermonde systems in O(n^2) time using the algorithm of Björck & Pereyra (1970), https://doi.org/10.2307/2004623.

julia> A \ a
5-element Vector{Float64}:
 0.0
 1.0
 0.0
 0.0
 0.0

julia> A' \ A[2,:]
5-element Vector{Float64}:
 0.0
 1.0
 0.0
 0.0
 0.0

det(A::InverseHilbert)

Explicit formula for the determinant. Caution: this function overflows easily.

dvand!(a, b) -> b

Solve system $A*x = b$ in-place.

A is Vandermonde matrix $A_{ij} = a_i^{j-1}$.

Algorithm by Bjorck & Pereyra, Mathematics of Computation, Vol. 24, No. 112 (1970), pp. 893-903, https://doi.org/10.2307/2004623

pvand!(a, b) -> b

Solve system $A^T*x = b$ in-place.

$A^T$ is transpose of Vandermonde matrix $A_{ij} = a_i^{j-1}$.

Algorithm by Bjorck & Pereyra, Mathematics of Computation, Vol. 24, No. 112 (1970), pp. 893-903, https://doi.org/10.2307/2004623

vandtype(T1::Type, T2::Type)

Determine the return type of Vandermonde{T1} \ Vector{T2}.