WIP: add oneunit(x) for dimensionful version of one(x)

NEWS.md

@@ -177,6 +177,9 @@ Library improvements
   * `max`, `min`, and related functions (`minmax`, `maximum`, `minimum`,
     `extrema`) now return `NaN` for `NaN` arguments ([#12563]).
 
+  * `oneunit(x)` function to return a dimensionful version of `one(x)` (which
+    is clarified to mean a dimensionless quantity if `x` is dimensionful) ([#20268]).
+
   * The `chop` and `chomp` functions now return a `SubString` ([#18339]).
 
   * Numbered stackframes printed in stacktraces can be opened in an editor by

base/Enums.jl

@@ -95,7 +95,7 @@ macro enum(T,syms...)
             lo = min(lo, i)
             hi = max(hi, i)
         end
-        i += one(i)
+        i += oneunit(i)
     end
     values = basetype[i[2] for i in vals]
     if hasexpr && values != unique(values)

base/array.jl

@@ -234,10 +234,10 @@ end
 `m`-by-`n` identity matrix.
 The default element type is `Float64`.
 """
-function eye(T::Type, m::Integer, n::Integer)
+function eye{T}(::Type{T}, m::Integer, n::Integer)
     a = zeros(T,m,n)
     for i = 1:min(m,n)
-        a[i,i] = one(T)
+        a[i,i] = oneunit(T)
     end
     return a
 end
@@ -248,7 +248,7 @@ end
 `m`-by-`n` identity matrix.
 """
 eye(m::Integer, n::Integer) = eye(Float64, m, n)
-eye(T::Type, n::Integer) = eye(T, n, n)
+eye{T}(::Type{T}, n::Integer) = eye(T, n, n)
 """
     eye([T::Type=Float64,] n::Integer)
 
@@ -278,14 +278,17 @@ julia> eye(A)
 
 Note the difference from [`ones`](@ref).
 """
-eye{T}(x::AbstractMatrix{T}) = eye(T, size(x, 1), size(x, 2))
+eye{T}(x::AbstractMatrix{T}) = eye(typeof(one(T)), size(x, 1), size(x, 2))
 
-function one{T}(x::AbstractMatrix{T})
+function _one{T}(unit::T, x::AbstractMatrix)
     m,n = size(x)
     m==n || throw(DimensionMismatch("multiplicative identity defined only for square matrices"))
     eye(T, m)
 end
 
+one{T}(x::AbstractMatrix{T}) = _one(one(T), x)
+oneunit{T}(x::AbstractMatrix{T}) = _one(oneunit(T), x)
+
 ## Conversions ##
 
 convert{T}(::Type{Vector}, x::AbstractVector{T}) = convert(Vector{T}, x)

base/bool.jl

@@ -92,7 +92,7 @@ abs2(x::Bool) = x
 ^(x::Integer, y::Bool) = ifelse(y, x, one(x))
 
 function +{T<:AbstractFloat}(x::Bool, y::T)::promote_type(Bool,T)
-    return ifelse(x, one(y) + y, y)
+    return ifelse(x, oneunit(y) + y, y)
 end
 +(y::AbstractFloat, x::Bool) = x + y
 

base/deprecated.jl

@@ -141,7 +141,7 @@ function convert(::Type{Base.LinAlg.UnitUpperTriangular}, A::Diagonal)
         "that convert `Diagonal`/`Bidiagonal` to `<:AbstractTriangular` are deprecated. ",
         "Consider calling the `UnitUpperTriangular` constructor directly ",
         "(`Base.LinAlg.UnitUpperTriangular(A)`) instead."), :convert)
-    if !all(x -> x == one(x), A.diag)
+    if !all(x -> x == oneunit(x), A.diag)
         throw(ArgumentError("matrix cannot be represented as UnitUpperTriangular"))
     end
     Base.LinAlg.UnitUpperTriangular(Array(A))
@@ -151,7 +151,7 @@ function convert(::Type{Base.LinAlg.UnitLowerTriangular}, A::Diagonal)
         "that convert `Diagonal`/`Bidiagonal` to `<:AbstractTriangular` are deprecated. ",
         "Consider calling the `UnitLowerTriangular` constructor directly ",
         "(`Base.LinAlg.UnitLowerTriangular(A)`) instead."), :convert)
-    if !all(x -> x == one(x), A.diag)
+    if !all(x -> x == oneunit(x), A.diag)
         throw(ArgumentError("matrix cannot be represented as UnitLowerTriangular"))
     end
     Base.LinAlg.UnitLowerTriangular(Array(A))

base/dft.jl

@@ -33,10 +33,10 @@ realfloat{T<:FFTWFloat}(x::StridedArray{T}) = x
 
 # return an Array, rather than similar(x), to avoid an extra copy for FFTW
 # (which only works on StridedArray types).
-complexfloat{T<:Complex}(x::AbstractArray{T}) = copy1(typeof(fftwfloat(one(T))), x)
-complexfloat{T<:Real}(x::AbstractArray{T}) = copy1(typeof(complex(fftwfloat(one(T)))), x)
+complexfloat{T<:Complex}(x::AbstractArray{T}) = copy1(typeof(fftwfloat(zero(T))), x)
+complexfloat{T<:Real}(x::AbstractArray{T}) = copy1(typeof(complex(fftwfloat(zero(T)))), x)
 
-realfloat{T<:Real}(x::AbstractArray{T}) = copy1(typeof(fftwfloat(one(T))), x)
+realfloat{T<:Real}(x::AbstractArray{T}) = copy1(typeof(fftwfloat(zero(T))), x)
 
 # copy to a 1-based array, using circular permutation
 function copy1{T}(::Type{T}, x)
@@ -262,7 +262,7 @@ summary(p::ScaledPlan) = string(p.scale, " * ", summary(p.p))
 *(p::Plan, I::UniformScaling) = ScaledPlan(p, I.λ)
 
 # Normalization for ifft, given unscaled bfft, is 1/prod(dimensions)
-normalization(T, sz, region) = (one(T) / prod([sz...][[region...]]))::T
+normalization(T, sz, region) = one(T) / Int(prod([sz...][[region...]]))
 normalization(X, region) = normalization(real(eltype(X)), size(X), region)
 
 plan_ifft(x::AbstractArray, region; kws...) =

base/exports.jl

@@ -386,6 +386,7 @@ export
     numerator,
     num2hex,
     one,
+    oneunit,
     powermod,
     prevfloat,
     prevpow,

base/fft/FFTW.jl

@@ -568,7 +568,7 @@ end
 # (FIXME: is there a way to use the Julia promotion rules more cleverly here?)
 fftwcomplex{T<:fftwComplex}(X::StridedArray{T}) = X
 fftwcomplex{T<:fftwReal}(X::AbstractArray{T}) =
-    copy!(Array{typeof(complex(one(T)))}(size(X)), X)
+    copy!(Array{typeof(complex(zero(T)))}(size(X)), X)
 fftwcomplex{T<:Real}(X::AbstractArray{T}) = copy!(Array{Complex128}(size(X)),X)
 fftwcomplex{T<:Complex}(X::AbstractArray{T}) =
     copy!(Array{Complex128}(size(X)), X)

base/intfuncs.jl

@@ -109,7 +109,7 @@ julia> gcdx(240, 46)
 """
 function gcdx{T<:Integer}(a::T, b::T)
     # a0, b0 = a, b
-    s0, s1 = one(T), zero(T)
+    s0, s1 = oneunit(T), zero(T)
     t0, t1 = s1, s0
     # The loop invariant is: s0*a0 + t0*b0 == a
     while b != 0
@@ -242,7 +242,7 @@ julia> nextpow2(17)
 32
 ```
 """
-nextpow2(x::Unsigned) = one(x)<<((sizeof(x)<<3)-leading_zeros(x-one(x)))
+nextpow2(x::Unsigned) = oneunit(x)<<((sizeof(x)<<3)-leading_zeros(x-oneunit(x)))
 nextpow2(x::Integer) = reinterpret(typeof(x),x < 0 ? -nextpow2(unsigned(-x)) : nextpow2(unsigned(x)))
 
 """

base/iterators.jl

@@ -248,7 +248,7 @@ end
 An iterator that counts forever, starting at `start` and incrementing by `step`.
 """
 countfrom(start::Number, step::Number) = Count(promote(start, step)...)
-countfrom(start::Number)               = Count(start, one(start))
+countfrom(start::Number)               = Count(start, oneunit(start))
 countfrom()                            = Count(1, 1)
 
 eltype{S}(::Type{Count{S}}) = S

base/linalg/dense.jl

@@ -165,7 +165,7 @@ tril(M::Matrix, k::Integer) = tril!(copy(M), k)
 
 function gradient(F::AbstractVector, h::Vector)
     n = length(F)
-    T = typeof(one(eltype(F))/one(eltype(h)))
+    T = typeof(oneunit(eltype(F))/oneunit(eltype(h)))
     g = similar(F, T)
     if n == 1
         g[1] = zero(T)

base/linalg/factorization.jl

@@ -42,7 +42,7 @@ for (f1, f2) in ((:\, :A_ldiv_B!),
                  (:Ac_ldiv_B, :Ac_ldiv_B!))
     @eval begin
         function $f1(F::Factorization, B::AbstractVecOrMat)
-            TFB = typeof(one(eltype(F)) / one(eltype(B)))
+            TFB = typeof(oneunit(eltype(B)) / oneunit(eltype(F)))
             BB = similar(B, TFB, size(B))
             copy!(BB, B)
             $f2(convert(Factorization{TFB}, F), BB)

base/linalg/generic.jl

@@ -447,7 +447,7 @@ function vecnorm(itr, p::Real=2)
         vecnormp(itr,p)
     end
 end
-@inline vecnorm(x::Number, p::Real=2) = p == 0 ? (x==0 ? zero(abs(x)) : one(abs(x))) : abs(x)
+@inline vecnorm(x::Number, p::Real=2) = p == 0 ? (x==0 ? zero(abs(x)) : oneunit(abs(x))) : abs(x)
 
 norm(x::AbstractVector, p::Real=2) = vecnorm(x, p)
 
@@ -719,8 +719,9 @@ true
 ```
 """
 function inv{T}(A::AbstractMatrix{T})
-    S = typeof(zero(T)/one(T))
-    A_ldiv_B!(factorize(convert(AbstractMatrix{S}, A)), eye(S, checksquare(A)))
+    S = typeof(zero(T)/one(T))      # dimensionful
+    S0 = typeof(zero(T)/oneunit(T)) # dimensionless
+    A_ldiv_B!(factorize(convert(AbstractMatrix{S}, A)), eye(S0, checksquare(A)))
 end
 
 """

base/linalg/givens.jl

@@ -55,14 +55,16 @@ realmin2{T}(::Type{T}) = (twopar = 2one(T); twopar^trunc(Integer,log(realmin(T)/
 # Univ. of Colorado Denver
 # NAG Ltd.
 function givensAlgorithm{T<:AbstractFloat}(f::T, g::T)
-    zeropar = zero(T)
     onepar = one(T)
     twopar = 2one(T)
+    T0 = typeof(onepar) # dimensionless
+    zeropar = T0(zero(T)) # must be dimensionless
 
-    safmin = realmin(T)
-    epspar = eps(T)
-    safmn2 = realmin2(T)
-    safmx2 = onepar/safmn2
+    # need both dimensionful and dimensionless versions of these:
+    safmn2 = realmin2(T0)
+    safmn2u = realmin2(T)
+    safmx2 = one(T)/safmn2
+    safmx2u = oneunit(T)/safmn2
 
     if g == 0
         cs = onepar
@@ -76,29 +78,29 @@ function givensAlgorithm{T<:AbstractFloat}(f::T, g::T)
         f1 = f
         g1 = g
         scalepar = max(abs(f1), abs(g1))
-        if scalepar >= safmx2
+        if scalepar >= safmx2u
             count = 0
             while true
                 count += 1
                 f1 *= safmn2
                 g1 *= safmn2
                 scalepar = max(abs(f1), abs(g1))
-                if scalepar < safmx2 break end
+                if scalepar < safmx2u break end
             end
             r = sqrt(f1*f1 + g1*g1)
             cs = f1/r
             sn = g1/r
             for i = 1:count
                 r *= safmx2
             end
-        elseif scalepar <= safmn2
+        elseif scalepar <= safmn2u
             count = 0
             while true
                 count += 1
                 f1 *= safmx2
                 g1 *= safmx2
                 scalepar = max(abs(f1), abs(g1))
-                if scalepar > safmn2 break end
+                if scalepar > safmn2u break end
             end
             r = sqrt(f1*f1 + g1*g1)
             cs = f1/r
@@ -127,27 +129,30 @@ end
 # Univ. of Colorado Denver
 # NAG Ltd.
 function givensAlgorithm{T<:AbstractFloat}(f::Complex{T}, g::Complex{T})
-    twopar, onepar, zeropar = 2one(T), one(T), zero(T)
-    czero = zero(Complex{T})
+    twopar, onepar = 2one(T), one(T)
+    T0 = typeof(onepar) # dimensionless
+    zeropar = T0(zero(T)) # must be dimensionless
+    czero = complex(zeropar)
 
     abs1(ff) = max(abs(real(ff)), abs(imag(ff)))
-    safmin = realmin(T)
-    epspar = eps(T)
-    safmn2 = realmin2(T)
-    safmx2 = onepar/safmn2
+    safmin = realmin(T0)
+    safmn2 = realmin2(T0)
+    safmn2u = realmin2(T)
+    safmx2 = one(T)/safmn2
+    safmx2u = oneunit(T)/safmn2
     scalepar = max(abs1(f), abs1(g))
     fs = f
     gs = g
     count = 0
-    if scalepar >= safmx2
+    if scalepar >= safmx2u
         while true
             count += 1
             fs *= safmn2
             gs *= safmn2
             scalepar *= safmn2
-            if scalepar < safmx2 break end
+            if scalepar < safmx2u break end
         end
-    elseif scalepar <= safmn2
+    elseif scalepar <= safmn2u
         if g == 0
             cs = onepar
             sn = czero
@@ -159,12 +164,12 @@ function givensAlgorithm{T<:AbstractFloat}(f::Complex{T}, g::Complex{T})
             fs *= safmx2
             gs *= safmx2
             scalepar *= safmx2
-            if scalepar > safmn2 break end
+            if scalepar > safmn2u break end
         end
     end
     f2 = abs2(fs)
     g2 = abs2(gs)
-    if f2 <= max(g2, onepar)*safmin
+    if f2 <= max(g2, oneunit(T))*safmin
 
      # This is a rare case: F is very small.
 
@@ -305,7 +310,7 @@ function getindex(G::Givens, i::Integer, j::Integer)
         if i == G.i1 || i == G.i2
             G.c
         else
-            one(G.c)
+            oneunit(G.c)
         end
     elseif i == G.i1 && j == G.i2
         G.s

base/linalg/hessenberg.jl

@@ -48,7 +48,7 @@ julia> F[:Q] * F[:H] * F[:Q]'
 ```
 """
 function hessfact{T}(A::StridedMatrix{T})
-    S = promote_type(Float32, typeof(one(T)/norm(one(T))))
+    S = promote_type(Float32, typeof(zero(T)/norm(one(T))))
     return hessfact!(copy_oftype(A, S))
 end
 

base/linalg/linalg.jl

@@ -8,7 +8,7 @@ import Base: A_mul_Bt, At_ldiv_Bt, A_rdiv_Bc, At_ldiv_B, Ac_mul_Bc, A_mul_Bc, Ac
 import Base: USE_BLAS64, abs, big, broadcast, ceil, conj, convert, copy, copy!,
     ctranspose, eltype, eye, findmax, findmin, fill!, floor, full, getindex,
     hcat, imag, indices, inv, isapprox, kron, length, linearindexing, map,
-    ndims, parent, power_by_squaring, print_matrix, promote_rule, real, round,
+    ndims, oneunit, parent, power_by_squaring, print_matrix, promote_rule, real, round,
     setindex!, show, similar, size, transpose, trunc, typed_hcat
 using Base: promote_op, _length, iszero, @pure, @propagate_inbounds, LinearFast,
     reduce, hvcat_fill, typed_vcat, promote_typeof

base/linalg/qr.jl

@@ -246,7 +246,7 @@ function qr(v::AbstractVector)
         return __normalize!(vv, nrm), nrm
     else
         T = typeof(zero(eltype(v))/nrm)
-        return T[], one(T)
+        return T[], oneunit(T)
     end
 end