Skip to content

Commit

Permalink
Merge pull request #12104 from JuliaLang/anj/rectfull
Browse files Browse the repository at this point in the history 
Move rectuangular full packed matrix support from base to LinearAlgebra.jl
  • Loading branch information
@andreasnoack
andreasnoack committed Jul 13, 2015
2 parents d4aabee + 1af3420 commit ab4f684
Show file tree
Hide file tree
Showing 3 changed files with 0 additions and 239 deletions.
1 change: 0 additions & 1 deletion base/linalg.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -224,7 +224,6 @@ include("linalg/symmetric.jl")
include("linalg/diagonal.jl")
include("linalg/bidiag.jl")
include("linalg/uniformscaling.jl")
include("linalg/rectfullpacked.jl")
include("linalg/givens.jl")
include("linalg/special.jl")
include("linalg/bitarray.jl")
Expand Down
198 changes: 0 additions & 198 deletions base/linalg/lapack.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4085,204 +4085,6 @@ for (trsen, tgsen, elty) in
end
end

### Rectangular full packed format

# Symmetric rank-k operation for matrix in RFP format.
for (fn, elty, relty) in ((:dsfrk_, :Float64, :Float64),
(:ssfrk_, :Float32, :Float32),
(:zhfrk_, :Complex128, :Float64),
(:chfrk_, :Complex64, :Float32))
@eval begin
function sfrk!(transr::Char, uplo::Char, trans::Char, alpha::Real, A::StridedMatrix{$elty}, beta::Real, C::StridedVector{$elty})
chkuplo(uplo)
chkstride1(A)
if trans == 'N' || trans == 'n'
n, k = size(A)
elseif trans == 'T' || trans == 't'
k, n = size(A)
end
lda = max(1, stride(A, 2))
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt},
Ptr{BlasInt}, Ptr{$relty}, Ptr{$elty}, Ptr{BlasInt},
Ptr{$relty}, Ptr{$elty}),
&transr, &uplo, &trans, &n,
&k, &alpha, A, &lda,
&beta, C)
C
end
end
end

# Cholesky factorization of a real symmetric positive definite matrix A
for (fn, elty) in ((:dpftrf_, :Float64),
(:spftrf_, :Float32),
(:zpftrf_, :Complex128),
(:cpftrf_, :Complex64))
@eval begin
function pftrf!(transr::Char, uplo::Char, A::StridedVector{$elty})
chkuplo(uplo)
n = round(Int,div(sqrt(8length(A)), 2))
info = Array(BlasInt, 1)
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt}, Ptr{$elty},
Ptr{BlasInt}),
&transr, &uplo, &n, A,
info)
@assertargsok
@assertnonsingular
A
end
end
end

# Computes the inverse of a (real) symmetric positive definite matrix A using the Cholesky factorization
for (fn, elty) in ((:dpftri_, :Float64),
(:spftri_, :Float32),
(:zpftri_, :Complex128),
(:cpftri_, :Complex64))
@eval begin
function pftri!(transr::Char, uplo::Char, A::StridedVector{$elty})
chkuplo(uplo)
n = round(Int,div(sqrt(8length(A)), 2))
info = Array(BlasInt, 1)
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt}, Ptr{$elty},
Ptr{BlasInt}),
&transr, &uplo, &n, A,
info)
@assertargsok
@assertnonsingular
A
end
end
end

# DPFTRS solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization
for (fn, elty) in ((:dpftrs_, :Float64),
(:spftrs_, :Float32),
(:zpftrs_, :Complex128),
(:cpftrs_, :Complex64))
@eval begin
function pftrs!(transr::Char, uplo::Char, A::StridedVector{$elty}, B::StridedVecOrMat{$elty})
chkuplo(uplo)
chkstride1(B)
n = round(Int,div(sqrt(8length(A)), 2))
if n != size(B, 1)
throw(DimensionMismatch("B has first dimension $(size(B,1)) but needs $n"))
end
nhrs = size(B, 2)
ldb = max(1, stride(B, 2))
info = Array(BlasInt, 1)
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt}, Ptr{BlasInt},
Ptr{$elty}, Ptr{$elty}, Ptr{BlasInt}, Ptr{BlasInt}),
&transr, &uplo, &n, &nhrs,
A, B, &ldb, info)
@assertargsok
@assertposdef
B
end
end
end

# Solves a matrix equation (one operand is a triangular matrix in RFP format)
for (fn, elty) in ((:dtfsm_, :Float64),
(:stfsm_, :Float32),
(:ztfsm_, :Complex128),
(:ctfsm_, :Complex64))
@eval begin
function pftrs!(transr::Char, side::Char, uplo::Char, trans::Char, diag::Char, alpha::Real, A::StridedVector{$elty}, B::StridedMatrix{$elty})
chkuplo(uplo)
chkside(side)
chkdiag(diag)
chkstride1(B)
m, n = size(B)
if round(Int, div(sqrt(8length(A)), 2)) != m
throw(DimensionMismatch("First dimension of B must equal $(round(Int, div(sqrt(8length(A)), 2))), got $m"))
end
ldb = max(1, stride(B, 2))
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{UInt8}, Ptr{UInt8},
Ptr{UInt8}, Ptr{BlasInt}, Ptr{BlasInt}, Ptr{$elty},
Ptr{$elty}, Ptr{$elty}, Ptr{BlasInt}),
&transr, &side, &uplo, &trans,
&diag, &m, &n, &alpha,
A, B, &ldb)
B
end
end
end

# Computes the inverse of a triangular matrix A stored in RFP format.
for (fn, elty) in ((:dtftri_, :Float64),
(:stftri_, :Float32),
(:ztftri_, :Complex128),
(:ctftri_, :Complex64))
@eval begin
function tftri!(transr::Char, uplo::Char, diag::Char, A::StridedVector{$elty})
chkuplo(uplo)
chkdiag(diag)
n = round(Int,div(sqrt(8length(A)), 2))
info = Array(BlasInt, 1)
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt},
Ptr{$elty}, Ptr{BlasInt}),
&transr, &uplo, &diag, &n,
A, info)
@assertargsok
@assertnonsingular
A
end
end
end

# Copies a triangular matrix from the rectangular full packed format (TF) to the standard full format (TR)
for (fn, elty) in ((:dtfttr_, :Float64),
(:stfttr_, :Float32),
(:ztfttr_, :Complex128),
(:ctfttr_, :Complex64))
@eval begin
function tfttr!(transr::Char, uplo::Char, Arf::StridedVector{$elty})
chkuplo(uplo)
n = round(Int,div(sqrt(8length(Arf)), 2))
info = Array(BlasInt, 1)
A = similar(Arf, $elty, n, n)
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt}, Ptr{$elty},
Ptr{$elty}, Ptr{BlasInt}, Ptr{BlasInt}),
&transr, &uplo, &n, Arf,
A, &n, info)
@assertargsok
A
end
end
end

# Copies a triangular matrix from the standard full format (TR) to the rectangular full packed format (TF).
for (fn, elty) in ((:dtrttf_, :Float64),
(:strttf_, :Float32),
(:ztrttf_, :Complex128),
(:ctrttf_, :Complex64))
@eval begin
function trttf!(transr::Char, uplo::Char, A::StridedMatrix{$elty})
chkuplo(uplo)
chkstride1(A)
n = size(A, 1)
lda = max(1, stride(A, 2))
info = Array(BlasInt, 1)
Arf = similar(A, $elty, div(n*(n+1), 2))
ccall(($(blasfunc(fn)), liblapack), Void,
(Ptr{UInt8}, Ptr{UInt8}, Ptr{BlasInt}, Ptr{$elty},
Ptr{BlasInt}, Ptr{$elty}, Ptr{BlasInt}),
&transr, &uplo, &n, A,
&lda, Arf, info)
@assertargsok
Arf
end
end
end

# Solves the real Sylvester matrix equation: op(A)*X +- X*op(B) = scale*C and A and B are both upper quasi triangular.
for (fn, elty, relty) in ((:dtrsyl_, :Float64, :Float64),
(:strsyl_, :Float32, :Float32),
Expand Down
40 changes: 0 additions & 40 deletions base/linalg/rectfullpacked.jl
View file

This file was deleted.

0 comments on commit ab4f684

Please sign in to comment.