Skip to content

Commit

Permalink
WIP: add oneunit(x) for dimensionful version of one(x) (#20268)
Browse files Browse the repository at this point in the history
add oneunit(x) for dimensionful version of one(x), and change one -> oneunit where appropriate.
  • Loading branch information
stevengj authored Feb 6, 2017
1 parent e3bb870 commit 7fb758a
Show file tree
Hide file tree
Showing 33 changed files with 191 additions and 128 deletions.
3 changes: 3 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -180,6 +180,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
Expand Down
2 changes: 1 addition & 1 deletion base/Enums.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
13 changes: 8 additions & 5 deletions base/array.jl
Original file line number Diff line number Diff line change
Expand Up @@ -237,10 +237,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
Expand All @@ -251,7 +251,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)
Expand Down Expand Up @@ -281,14 +281,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)
Expand Down
2 changes: 1 addition & 1 deletion base/bool.jl
Original file line number Diff line number Diff line change
Expand Up @@ -93,7 +93,7 @@ iszero(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

Expand Down
4 changes: 2 additions & 2 deletions base/deprecated.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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))
Expand All @@ -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))
Expand Down
8 changes: 4 additions & 4 deletions base/dft.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down Expand Up @@ -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...) =
Expand Down
1 change: 1 addition & 0 deletions base/exports.jl
Original file line number Diff line number Diff line change
Expand Up @@ -382,6 +382,7 @@ export
numerator,
num2hex,
one,
oneunit,
powermod,
prevfloat,
prevpow,
Expand Down
2 changes: 1 addition & 1 deletion base/fft/FFTW.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
4 changes: 2 additions & 2 deletions base/intfuncs.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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)))

"""
Expand Down
2 changes: 1 addition & 1 deletion base/iterators.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down
2 changes: 1 addition & 1 deletion base/linalg/dense.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
2 changes: 1 addition & 1 deletion base/linalg/factorization.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)
Expand Down
7 changes: 4 additions & 3 deletions base/linalg/generic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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)

Expand Down Expand Up @@ -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

"""
Expand Down
47 changes: 26 additions & 21 deletions base/linalg/givens.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand All @@ -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
Expand Down Expand Up @@ -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
Expand All @@ -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.

Expand Down Expand Up @@ -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
Expand Down
2 changes: 1 addition & 1 deletion base/linalg/hessenberg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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

Expand Down
2 changes: 1 addition & 1 deletion base/linalg/linalg.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down
2 changes: 1 addition & 1 deletion base/linalg/qr.jl
Original file line number Diff line number Diff line change
Expand Up @@ -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

Expand Down
Loading

0 comments on commit 7fb758a

Please sign in to comment.