Add qr(::AbstractVector), qr!(::AbstractVector)

Compute polar decomposition of vector

`qr[!]` is equivalent to `v->(normalize[!](v), norm(v))` but is
convenient for reusing the norm calculation, e.g. for Krylov subspace
algorithms.

Also refactors the in place normalization code to a separate routine,
`__normalize!()`, so that it can be reused by `qr!`

base/linalg/generic.jl

@@ -547,19 +547,22 @@ Normalize the vector `v` in-place with respect to the `p`-norm.
 
 # See also
 
-`normalize`
+`normalize`, `qr`
 
 """
 function normalize!(v::AbstractVector, p::Real=2)
     nrm = norm(v, p)
+    __normalize!(v, nrm)
+end
 
+@inline function __normalize!{T<:AbstractFloat}(v::AbstractVector{T}, nrm::T)
     #The largest positive floating point number whose inverse is less than
     #infinity
-    const δ = inv(prevfloat(typemax(float(nrm))))
+    const δ = inv(prevfloat(typemax(T)))
 
     if nrm ≥ δ #Safe to multiply with inverse
         invnrm = inv(nrm)
-        scale!(v, invrm)
+        scale!(v, invnrm)
 
     else #Divide by norm; slower but more correct
         #Note 2015-10-19: As of Julia 0.4, @simd does not vectorize floating
@@ -590,7 +593,7 @@ Normalize the vector `v` with respect to the `p`-norm.
 
 # See also
 
-`normalize!`
+`normalize!`, `qr`
 """
 normalize(v::AbstractVector, p::Real=2) = v/norm(v, p)
 

base/linalg/qr.jl

@@ -108,6 +108,49 @@ function _qr(A::Union{Number, AbstractMatrix}, ::Type{Val{true}}; thin::Bool=tru
     full(getq(F), thin=thin), F[:R]::Matrix{eltype(F)}, F[:p]::Vector{BlasInt}
 end
 
+"""
+    qr(v::AbstractVector)
+
+Computes the polar decomposition of a vector.
+
+# Input
+- `v::AbstractVector` - vector to normalize
+
+# Outputs
+- `w` - A unit vector in the direction of `v`
+- `r` - The norm of `v`
+
+# See also
+
+`normalize`, `normalize!`, `qr!`
+"""
+function qr(v::AbstractVector)
+    nrm = norm(v)
+    v/nrm, nrm
+end
+
+"""
+    qr!(v::AbstractVector)
+
+Computes the polar decomposition of a vector.
+
+# Input
+- `v::AbstractVector` - vector to normalize
+
+# Outputs
+- `w` - A unit vector in the direction of `v`
+- `r` - The norm of `v`
+
+# See also
+
+`normalize`, `normalize!`, `qr`
+"""
+function qr!(v::AbstractVector)
+    nrm = norm(v)
+    __normalize!(v, nrm), nrm
+end
+
+
 convert{T}(::Type{QR{T}},A::QR) = QR(convert(AbstractMatrix{T}, A.factors), convert(Vector{T}, A.τ))
 convert{T}(::Type{Factorization{T}}, A::QR) = convert(QR{T}, A)
 convert{T}(::Type{QRCompactWY{T}},A::QRCompactWY) = QRCompactWY(convert(AbstractMatrix{T}, A.factors), convert(AbstractMatrix{T}, A.T))

test/linalg/generic.jl

@@ -157,17 +157,17 @@ let
     w = normalize(v)
     @test w == [0.6, 0.8]
     @test norm(w) === 1.0
-    @test normalize!(v) == w
+    @test norm(normalize!(v) - w, Inf) < eps()
 end
 
 #Test potential overflow in normalize!
 let
-    δ = inv(prevfloat(typemax(float(nrm))))
+    δ = inv(prevfloat(typemax(Float64)))
     v = [δ, -δ]
 
     @test norm(v) === 7.866824069956793e-309
     w = normalize(v)
-    @test norm(w) === 1.0
     @test w ≈ [1/√2, -1/√2]
-    @test normalize!(v) == w
+    @test norm(w) === 1.0
+    @test norm(normalize!(v) - w, Inf) < eps()
 end

test/linalg/qr.jl

@@ -162,3 +162,11 @@ let
     Q=full(qrfact(A)[:Q])
     @test vecnorm(A-Q) < eps()
 end
+
+let
+    debug && println("qr on AbstractVector")
+
+    v = [3.0, 4.0]
+    @test qr(v) == ([0.6, 0.8], 5.0)
+end
+