Skip to content

Commit

Permalink
Rename QRPDense to QRDensePivoted. Use qr and qrpivot instead of qrd …
Browse files Browse the repository at this point in the history
…and qrpd.

Update docs for qr.
  • Loading branch information
ViralBShah committed Feb 7, 2013
1 parent c5b0fe0 commit 2fe3c04
Show file tree
Hide file tree
Showing 4 changed files with 40 additions and 38 deletions.
10 changes: 4 additions & 6 deletions base/exports.jl
Original file line number Diff line number Diff line change
Expand Up @@ -109,7 +109,7 @@ export
LUTridiagonal,
LDLTTridiagonal,
QRDense,
QRPDense,
QRDensePivoted,
InsertionSort,
QuickSort,
MergeSort,
Expand Down Expand Up @@ -559,11 +559,9 @@ export
null,
pinv,
qr,
qrd!,
qrd,
qrp,
qrpd!,
qrpd,
qr!,
qrpivot,
qrpivot!,
randsym,
rank,
rref,
Expand Down
42 changes: 21 additions & 21 deletions base/linalg_dense.jl
Original file line number Diff line number Diff line change
Expand Up @@ -573,17 +573,16 @@ end
size(A::QRDense) = size(A.hh)
size(A::QRDense,n) = size(A.hh,n)

qrd!{T<:BlasFloat}(A::StridedMatrix{T}) = QRDense{T}(LAPACK.geqrf!(A)...)
qrd{T<:BlasFloat}(A::StridedMatrix{T}) = qrd!(copy(A))
qrd{T<:Real}(A::StridedMatrix{T}) = qrd(float64(A))
qr!{T<:BlasFloat}(A::StridedMatrix{T}) = QRDense{T}(LAPACK.geqrf!(A)...)
qr{T<:BlasFloat}(A::StridedMatrix{T}) = qr!(copy(A))
qr{T<:Real}(A::StridedMatrix{T}) = qr(float64(A))

function factors{T<:BlasFloat}(qrd::QRDense{T})
aa = copy(qrd.hh)
function factors{T<:BlasFloat}(qr::QRDense{T})
aa = copy(qr.hh)
R = triu(aa[1:min(size(aa)),:]) # must be *before* call to orgqr!
LAPACK.orgqr!(aa, qrd.tau, size(aa,2)), R
LAPACK.orgqr!(aa, qr.tau, size(aa,2)), R
end

qr{T<:Number}(x::StridedMatrix{T}) = factors(qrd(x))
qr(x::Number) = (one(x), x)

## Multiplication by Q from the QR decomposition
Expand All @@ -602,41 +601,42 @@ function (\){T<:BlasFloat}(A::QRDense{T}, B::StridedVecOrMat{T})
return ans
end

type QRPDense{T} <: Factorization{T}
## Pivoted QR decomposition

type QRDensePivoted{T} <: Factorization{T}
hh::Matrix{T}
tau::Vector{T}
jpvt::Vector{BlasInt}
function QRPDense(hh::Matrix{T}, tau::Vector{T}, jpvt::Vector{BlasInt})
function QRDensePivoted(hh::Matrix{T}, tau::Vector{T}, jpvt::Vector{BlasInt})
m, n = size(hh)
if length(tau) != min(m,n) || length(jpvt) != n
error("QRPDense: mismatched dimensions")
error("QRDensePivoted: mismatched dimensions")
end
new(hh,tau,jpvt)
end
end
size(x::QRPDense) = size(x.hh)
size(x::QRPDense,d) = size(x.hh,d)
size(x::QRDensePivoted) = size(x.hh)
size(x::QRDensePivoted,d) = size(x.hh,d)
## Multiplication by Q from the QR decomposition
(*){T<:BlasFloat}(A::QRPDense{T}, B::StridedVecOrMat{T}) =
(*){T<:BlasFloat}(A::QRDensePivoted{T}, B::StridedVecOrMat{T}) =
LAPACK.ormqr!('L', 'N', A.hh, size(A,2), A.tau, copy(B))
## Multiplication by Q' from the QR decomposition
Ac_mul_B{T<:BlasFloat}(A::QRPDense{T}, B::StridedVecOrMat{T}) =
Ac_mul_B{T<:BlasFloat}(A::QRDensePivoted{T}, B::StridedVecOrMat{T}) =
LAPACK.ormqr!('L', iscomplex(A.tau)?'C':'T', A.hh, size(A,2), A.tau, copy(B))

qrpd!{T<:BlasFloat}(A::StridedMatrix{T}) = QRPDense{T}(LAPACK.geqp3!(A)...)
qrpd{T<:BlasFloat}(A::StridedMatrix{T}) = qrpd!(copy(A))
qrpd{T<:Real}(x::StridedMatrix{T}) = qrpd(float64(x))
qrpivot!{T<:BlasFloat}(A::StridedMatrix{T}) = QRDensePivoted{T}(LAPACK.geqp3!(A)...)
qrpivot{T<:BlasFloat}(A::StridedMatrix{T}) = qrpivot!(copy(A))
qrpivot{T<:Real}(x::StridedMatrix{T}) = qrpivot(float64(x))

function factors{T<:BlasFloat}(x::QRPDense{T})
function factors{T<:BlasFloat}(x::QRDensePivoted{T})
aa = copy(x.hh)
R = triu(aa[1:min(size(aa)),:])
LAPACK.orgqr!(aa, x.tau, size(aa,2)), R, x.jpvt
end

qrp{T<:BlasFloat}(x::StridedMatrix{T}) = factors(qrpd(x))
qrp{T<:Real}(x::StridedMatrix{T}) = qrp(float64(x))
qrpivot{T<:Real}(x::StridedMatrix{T}) = qrpivot(float64(x))

function (\){T<:BlasFloat}(A::QRPDense{T}, B::StridedVecOrMat{T})
function (\){T<:BlasFloat}(A::QRDensePivoted{T}, B::StridedVecOrMat{T})
n = length(A.tau)
x, info = LAPACK.trtrs!('U','N','N',A.hh[1:n,:],(A'*B)[1:n,:])
if info > 0; throw(LAPACK.SingularException(info)); end
Expand Down
14 changes: 11 additions & 3 deletions doc/stdlib/base.rst
Original file line number Diff line number Diff line change
Expand Up @@ -1742,11 +1742,19 @@ Linear algebra functions in Julia are largely implemented by calling functions f

.. function:: qr(A)

Compute QR factorization
Compute the QR factorization of ``A`` and return a ``QRDense`` object. ``factors(qr(A))`` returns ``Q`` and ``R``. The following functions are available for ``QRDense`` objects: ``size``, ``factors``, ``*``, ``Ac_mul_B``, ``\``.

.. function:: qrp(A)
.. function:: qr!(A)

Compute QR factorization with pivoting
``qr!`` is the same as ``qr``, but overwrites the input matrix A with the factorization.

.. function:: qrpivot(A)

Compute the QR factorization of ``A`` and return a ``QRDensePivoted`` object. ``factors(qr(A))`` returns ``Q`` and ``R``. The following functions are available for ``QRDensePivoted`` objects: ``size``, ``factors``, ``*``, ``Ac_mul_B``, ``\``.

.. function:: qrpivot!(A)

``qrpivot!`` is the same as ``qrpivot``, but overwrites the input matrix A with the factorization.

.. function:: factors(F)

Expand Down
12 changes: 4 additions & 8 deletions test/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -39,23 +39,19 @@ for elty in (Float32, Float64, Complex64, Complex128)
@assert_approx_eq a * inv(lua) eye(elty, n)
@assert_approx_eq a*(lua\b) b

qra = qrd(a) # QR decomposition
qra = qr(a) # QR decomposition
q,r = factors(qra)
@assert_approx_eq q'*q eye(elty, n)
@assert_approx_eq q*q' eye(elty, n)
Q,R = qr(a)
@test q == Q && r == R
@assert_approx_eq q*r a
@assert_approx_eq qra*b Q*b
@assert_approx_eq qra'*b Q'*b
@assert_approx_eq qra*b q*b
@assert_approx_eq qra'*b q'*b
@assert_approx_eq a*(qra\b) b

qrpa = qrpd(a) # pivoted QR decomposition
qrpa = qrpivot(a) # pivoted QR decomposition
q,r,p = factors(qrpa)
@assert_approx_eq q'*q eye(elty, n)
@assert_approx_eq q*q' eye(elty, n)
Q,R,P = qrp(a)
@test q == Q && r == R && p == P
@assert_approx_eq q*r a[:,p]
@assert_approx_eq q*r[:,invperm(p)] a
@assert_approx_eq a*(qrpa\b) b
Expand Down

0 comments on commit 2fe3c04

Please sign in to comment.