Skip to content
7 changes: 4 additions & 3 deletions src/ToeplitzMatrices.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,9 @@ module ToeplitzMatrices
# import StatsBase: levinson!, levinson
import DSP: conv

import Base: adjoint, convert, transpose, size, getindex, similar, copy, getproperty, inv, sqrt, copyto!, reverse, conj, zero, fill!, checkbounds, real, imag, isfinite, DimsInteger, iszero
import Base: adjoint, convert, transpose, size, getindex, similar, copy, getproperty, inv, sqrt, copyto!, reverse, conj, zero, fill!, checkbounds, real, isfinite, DimsInteger, _all
import Base: ==, +, -, *, \
import Base.Broadcast: broadcasted, DefaultMatrixStyle
import LinearAlgebra: Cholesky, Factorization
import LinearAlgebra: ldiv!, factorize, lmul!, pinv, eigvals, cholesky!, cholesky, tril!, triu!, checksquare, rmul!, dot, mul!, tril, triu
import LinearAlgebra: UpperTriangular, LowerTriangular, Symmetric, Adjoint
Expand All @@ -16,7 +17,7 @@ export durbin, trench, levinson
include("iterativeLinearSolvers.jl")

# Abstract
abstract type AbstractToeplitz{T<:Number} <: AbstractMatrix{T} end
abstract type AbstractToeplitz{T} <: AbstractMatrix{T} end

size(A::AbstractToeplitz) = (length(A.vc),length(A.vr))
@inline _vr(A::AbstractToeplitz) = A.vr
Expand All @@ -28,7 +29,7 @@ AbstractArray{T}(A::AbstractToeplitz) where T = AbstractToeplitz{T}(A)
convert(::Type{AbstractToeplitz{T}}, A::AbstractToeplitz) where T = AbstractToeplitz{T}(A)

isconcrete(A::AbstractToeplitz) = isconcretetype(typeof(A.vc)) && isconcretetype(typeof(A.vr))
iszero(A::AbstractToeplitz) = iszero(A.vc) && iszero(A.vr)
_all(f, A::AbstractToeplitz, ::Colon) = _all(f, A.vc, :) && _all(f, A.vr, :)

"""
ToeplitzFactorization
Expand Down
10 changes: 4 additions & 6 deletions src/special.jl
Original file line number Diff line number Diff line change
Expand Up @@ -14,9 +14,10 @@ for TYPE in (:SymmetricToeplitz, :Circulant, :LowerTriangularToeplitz, :UpperTri

size(A::$TYPE) = (length(A.v),length(A.v))

broadcasted(::DefaultMatrixStyle, f, A::$TYPE) = $TYPE(f.(A.v))
broadcasted(::DefaultMatrixStyle, f, x::Number, A::$TYPE) = $TYPE(f.(x, A.v))
broadcasted(::DefaultMatrixStyle, f, A::$TYPE, x::Number) = $TYPE(f.(A.v, x))
adjoint(A::$TYPE) = transpose(conj(A))
(*)(scalar::Number, C::$TYPE) = $TYPE(scalar * C.v)
(*)(C::$TYPE, scalar::Number) = $TYPE(C.v * scalar)
(==)(A::$TYPE,B::$TYPE) = A.v==B.v
function zero!(A::$TYPE)
if isconcrete(A)
Expand All @@ -40,12 +41,9 @@ for TYPE in (:SymmetricToeplitz, :Circulant, :LowerTriangularToeplitz, :UpperTri
A
end
end
for fun in (:zero, :conj, :copy, :-, :real, :imag)
for fun in (:zero, :copy)
@eval $fun(A::$TYPE) = $TYPE($fun(A.v))
end
for fun in (:iszero,)
@eval $fun(A::$TYPE) = $fun(A.v)
end
for op in (:+, :-)
@eval $op(A::$TYPE,B::$TYPE) = $TYPE($op(A.v,B.v))
end
Expand Down
31 changes: 21 additions & 10 deletions src/toeplitz.jl
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ struct Toeplitz{T, VC<:AbstractVector{T}, VR<:AbstractVector{T}} <: AbstractToep
vr::VR

function Toeplitz{T, VC, VR}(vc::VC, vr::VR) where {T, VC<:AbstractVector{T}, VR<:AbstractVector{T}}
if first(vc) != first(vr)
if !isequal(first(vc), first(vr))
error("First element of the vectors must be the same")
end
return new{T,VC,VR}(vc, vr)
Expand Down Expand Up @@ -42,6 +42,22 @@ AbstractToeplitz{T}(A::Toeplitz) where T = Toeplitz{T}(A)
convert(::Type{Toeplitz{T}}, A::AbstractToeplitz) where {T} = Toeplitz{T}(A)
convert(::Type{Toeplitz}, A::AbstractToeplitz) = Toeplitz(A)

# Dims for disambiguity
for DIM in (:DimsInteger, :Dims)
@eval function similar(::Type{Toeplitz{T, VC, VR}}, dims::$DIM{2}) where {T, VC, VR}
vc = similar(VC, dims[1])
vr = similar(VR, dims[2])
vr[1] = vc[1]
Toeplitz(vc, vr)
end
@eval function similar(A::AbstractToeplitz, T::Type = eltype(A), dims::$DIM{2} = size(A))
vc = similar(A.vc, T, dims[1])
vr = similar(A.vr, T, dims[2])
vr[1] = vc[1]
Toeplitz{T}(vc, vr)
end
end

# Retrieve an entry
Base.@propagate_inbounds function getindex(A::AbstractToeplitz, i::Integer, j::Integer)
@boundscheck checkbounds(A,i,j)
Expand All @@ -53,6 +69,9 @@ Base.@propagate_inbounds function getindex(A::AbstractToeplitz, i::Integer, j::I
end
end

broadcasted(::DefaultMatrixStyle, f, A::AbstractToeplitz) = Toeplitz(f.(A.vc), f.(A.vr))
broadcasted(::DefaultMatrixStyle, f, x::Number, A::AbstractToeplitz) = Toeplitz(f.(x, A.vc), f.(x, A.vr))

checknonaliased(A::Toeplitz) = Base.mightalias(A.vc, A.vr) && throw(ArgumentError("Cannot modify Toeplitz matrices in place with aliased data"))

function tril!(A::Toeplitz, k::Integer=0)
Expand Down Expand Up @@ -104,13 +123,7 @@ end

adjoint(A::AbstractToeplitz) = transpose(conj(A))
transpose(A::AbstractToeplitz) = Toeplitz(A.vr, A.vc)
function similar(A::AbstractToeplitz, T::Type, dims::Dims{2})
vc=similar(A.vc, T, dims[1])
vr=similar(A.vr, T, dims[2])
vr[1]=vc[1]
Toeplitz{T}(vc,vr)
end
for fun in (:zero, :conj, :copy, :-, :real, :imag)
for fun in (:zero, :copy)
@eval $fun(A::AbstractToeplitz)=Toeplitz($fun(A.vc),$fun(A.vr))
end
for op in (:+, :-)
Expand All @@ -129,8 +142,6 @@ function fill!(A::Toeplitz, x::Number)
fill!(A.vr,x)
A
end
(*)(scalar::Number, C::AbstractToeplitz) = Toeplitz(scalar * C.vc, scalar * C.vr)
(*)(C::AbstractToeplitz,scalar::Number) = Toeplitz(C.vc * scalar, C.vr * scalar)

function lmul!(x::Number, A::Toeplitz)
checknonaliased(A)
Expand Down