Skip to content

Commit

Permalink
Real Quantities
Browse files Browse the repository at this point in the history
  • Loading branch information
aplavin committed Nov 11, 2023
1 parent 026f3d9 commit 8a45d7b
Show file tree
Hide file tree
Showing 8 changed files with 569 additions and 549 deletions.
8 changes: 0 additions & 8 deletions src/complex.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,14 +5,6 @@
#
# It is currently incomplete.

complex(z::Quantity{T,D,U}) where {T<:Complex,D,U} = z
function complex(x::Quantity{T}, y = zero(x)) where {T<:Real}
r, i = promote(x, y)
return Quantity(complex(ustrip(r), ustrip(i)), unit(r))
end
complex(::Type{Quantity{T,D,U}}) where {T,D,U} =
Quantity{complex(T),D,U}

# implement Base.widen for real and complex quantities because Unitful
# does not have an implementation for widen yet
Base.widen(::Type{Quantity{T,D,U}}) where {T,D,U} =
Expand Down
17 changes: 14 additions & 3 deletions src/conversion.jl
Original file line number Diff line number Diff line change
Expand Up @@ -73,7 +73,9 @@ function uconvert(a::Units, x::Quantity{T,D,U}) where {T,D,U}
end
end

function uconvert(a::Units, x::Number)
uconvert(a::Units, x::Complex) = complex(uconvert(a, real(x)), uconvert(a, imag(x)))

function uconvert(a::Units, x::Real)
if dimension(a) == NoDims
Quantity(x * convfact(a, NoUnits), a)
else
Expand All @@ -99,7 +101,15 @@ uconvert(a::Units, x::Missing) = missing
end
end

function convert(::Type{Quantity{T,D,U}}, x::Number) where {T,D,U}
function convert(::Type{Quantity{T,D,U}}, x::Real) where {T,D,U}
if dimension(x) == D
Quantity(T(uconvert(U(),x).val), U())
else
throw(DimensionError(U(),x))
end
end
# disambiguation:
function convert(::Type{Quantity{T,D,U}}, x::Quantity) where {T,D,U}
if dimension(x) == D
Quantity(T(uconvert(U(),x).val), U())
else
Expand Down Expand Up @@ -136,7 +146,8 @@ function convert(::Type{DimensionlessQuantity{T,U}}, x::Quantity) where {T,U}
end

convert(::Type{Number}, y::Quantity) = y
convert(::Type{Real}, y::Quantity) = y
convert(::Type{T}, y::Quantity) where {T <: Real} =
T(uconvert(NoUnits, y))
convert(::Type{T}, y::Quantity) where {T <: Complex} =
T(uconvert(NoUnits, y))
T(convert(real(T), y), zero(real(T)))
7 changes: 5 additions & 2 deletions src/promotion.jl
Original file line number Diff line number Diff line change
Expand Up @@ -94,7 +94,7 @@ function Base.promote_rule(::Type{Quantity{S1,D,U1}},
end

# number, quantity
function Base.promote_rule(::Type{Quantity{S,D,U}}, ::Type{T}) where {S,T <: Number,D,U}
function Base.promote_rule(::Type{Quantity{S,D,U}}, ::Type{T}) where {S,T <: Real,D,U}
if D == NoDims
promote_type(S,T,typeof(convfact(NoUnits,U())))
else
Expand All @@ -105,7 +105,7 @@ end
Base.promote_rule(::Type{S}, ::Type{T}) where {S<:AbstractIrrational,S2,T<:Quantity{S2}} =
promote_type(promote_type(S, real(S2)), T)

Base.promote_rule(::Type{Quantity{S}}, ::Type{T}) where {S,T <: Number} =
Base.promote_rule(::Type{Quantity{S}}, ::Type{T}) where {S,T <: Real} =
Quantity{promote_type(S,T)}

# With only one of these, you can get a segmentation fault because you # fall back to the
Expand All @@ -115,3 +115,6 @@ Base.promote_rule(::Type{Quantity{T}}, ::Type{Quantity{S,D,U}}) where {T,S,D,U}

Base.promote_rule(::Type{Quantity{S,D,U}}, ::Type{Quantity{T}}) where {T,S,D,U} =
Quantity{promote_type(T,S)}

Base.promote_rule(::Type{<:AbstractFloat}, ::Type{<:AbstractQuantity}) = Union{}
Base.promote_rule(::Type{BigFloat}, ::Type{<:AbstractQuantity}) = Union{}
97 changes: 50 additions & 47 deletions src/quantities.jl
Original file line number Diff line number Diff line change
Expand Up @@ -23,18 +23,19 @@ julia> Quantity(5, u"m")
Quantity(x::Number, y::Units) = _Quantity(x, y)
Quantity(x::Number, y::Units{()}) = x

*(x::Number, y::Units, z::Units...) = Quantity(x,*(y,z...))
*(x::Number, y::Units, z::Units...) = Quantity(ustrip(x),*(unit(x),y,z...))
*(x::Number, y::typeof(NoUnits)) = x
*(x::Units, y::Number) = *(y,x)

*(x::AbstractQuantity, y::Units, z::Units...) = Quantity(x.val, *(unit(x),y,z...))
*(x::AbstractQuantity, y::AbstractQuantity) = Quantity(x.val*y.val, unit(x)*unit(y))

function *(x::Number, y::AbstractQuantity)
function *(x::Real, y::AbstractQuantity)
y isa AffineQuantity &&
throw(AffineError("an invalid operation was attempted with affine quantities: $x*$y"))
return Quantity(x*y.val, unit(y))
end
function *(x::AbstractQuantity, y::Number)
function *(x::AbstractQuantity, y::Real)
x isa AffineQuantity &&
throw(AffineError("an invalid operation was attempted with affine quantities: $x*$y"))
return Quantity(x.val*y, unit(x))
Expand All @@ -48,7 +49,7 @@ end
# Division (units)
/(x::AbstractQuantity, y::Units) = Quantity(x.val, unit(x) / y)
/(x::Units, y::AbstractQuantity) = Quantity(1/y.val, x / unit(y))
/(x::Number, y::Units) = Quantity(x,inv(y))
/(x::Number, y::Units) = x * inv(y)
/(x::Units, y::Number) = (1/y) * x

//(x::AbstractQuantity, y::Units) = Quantity(x.val, unit(x) / y)
Expand All @@ -57,38 +58,39 @@ end
//(x::Units, y::Number) = (1//y) * x

/(x::AbstractQuantity, y::AbstractQuantity) = Quantity(/(x.val, y.val), unit(x) / unit(y))
/(x::AbstractQuantity, y::Number) = Quantity(/(x.val, y), unit(x) / unit(y))
/(x::Number, y::AbstractQuantity) = Quantity(/(x, y.val), unit(x) / unit(y))
/(x::AbstractQuantity, y::Real) = Quantity(/(x.val, y), unit(x) / unit(y))
/(x::Real, y::AbstractQuantity) = Quantity(/(x, y.val), unit(x) / unit(y))
//(x::AbstractQuantity, y::AbstractQuantity) = Quantity(//(x.val, y.val), unit(x) / unit(y))
//(x::AbstractQuantity, y::Number) = Quantity(//(x.val, y), unit(x) // unit(y))
//(x::Number, y::AbstractQuantity) = Quantity(//(x, y.val), unit(x) / unit(y))

# ambiguity resolution
//(x::AbstractQuantity, y::Complex) = Quantity(//(x.val, y), unit(x))

for f in (:fld, :cld)
@eval begin
function ($f)(x::AbstractQuantity, y::AbstractQuantity)
z = uconvert(unit(y), x) # TODO: use promote?
($f)(z.val,y.val)
end

($f)(x::Number, y::AbstractQuantity) = Quantity(($f)(x, ustrip(y)), unit(x) / unit(y))
($f)(x::Real, y::AbstractQuantity) = Quantity(($f)(x, ustrip(y)), unit(x) / unit(y))

($f)(x::AbstractQuantity, y::Number) = Quantity(($f)(ustrip(x), y), unit(x))
($f)(x::AbstractQuantity, y::Real) = Quantity(($f)(ustrip(x), y), unit(x))
end
end

function div(x::AbstractQuantity, y::AbstractQuantity, r...)
function div(x::AbstractQuantity, y::AbstractQuantity)
z = uconvert(unit(y), x) # TODO: use promote?
div(z.val,y.val)
end
function div(x::AbstractQuantity, y::AbstractQuantity, r::RoundingMode)
z = uconvert(unit(y), x) # TODO: use promote?
div(z.val,y.val, r...)
div(z.val,y.val, r)
end

function div(x::Number, y::AbstractQuantity, r...)
function div(x::Real, y::AbstractQuantity, r...)
Quantity(div(x, ustrip(y), r...), unit(x) / unit(y))
end

function div(x::AbstractQuantity, y::Number, r...)
function div(x::AbstractQuantity, y::Real, r...)
Quantity(div(ustrip(x), y, r...), unit(x))
end

Expand Down Expand Up @@ -202,8 +204,8 @@ for (_x,_y) in [(:fma, :_fma), (:muladd, :_muladd)]
end
end

sqrt(x::AbstractQuantity) = Quantity(sqrt(x.val), sqrt(unit(x)))
cbrt(x::AbstractQuantity) = Quantity(cbrt(x.val), cbrt(unit(x)))
sqrt(x::AbstractQuantityOrComplex) = Quantity(sqrt(ustrip(x)), sqrt(unit(x)))
cbrt(x::AbstractQuantityOrComplex) = Quantity(cbrt(ustrip(x)), cbrt(unit(x)))

for _y in (:sin, :cos, :tan, :asin, :acos, :atan, :sinh, :cosh, :tanh, :asinh, :acosh, :atanh,
:sinpi, :cospi, :tanpi, :sinc, :cosc, :cis, :cispi, :sincospi)
Expand All @@ -219,20 +221,15 @@ atan(y::AbstractQuantity, x::AbstractQuantity) = throw(DimensionError(x,y))

abs(x::AbstractQuantity) = Quantity(abs(x.val), unit(x))
abs2(x::AbstractQuantity) = Quantity(abs2(x.val), unit(x)*unit(x))
angle(x::AbstractQuantity{<:Complex}) = angle(x.val)

copysign(x::AbstractQuantity, y::Number) = Quantity(copysign(x.val,y/unit(y)), unit(x))
copysign(x::Number, y::AbstractQuantity) = copysign(x,y/unit(y))
copysign(x::AbstractQuantity, y::AbstractQuantity) = Quantity(copysign(x.val,y/unit(y)), unit(x))

flipsign(x::AbstractQuantity, y::Number) = Quantity(flipsign(x.val,y/unit(y)), unit(x))
flipsign(x::Number, y::AbstractQuantity) = flipsign(x,y/unit(y))
flipsign(x::AbstractQuantity, y::AbstractQuantity) = Quantity(flipsign(x.val,y/unit(y)), unit(x))
for T in (:Real, :Signed, :Float32, :Float64) # for disambiguation
@eval flipsign(x::$T, y::AbstractQuantity) = flipsign(x,y/unit(y))
end

for (i,j) in zip((:<, :<=, :isless), (:_lt, :_le, :_isless))
@eval ($i)(x::AbstractQuantity, y::AbstractQuantity) = ($j)(x,y)
@eval ($i)(x::AbstractQuantity, y::Number) = ($i)(promote(x,y)...)
@eval ($i)(x::Number, y::AbstractQuantity) = ($i)(promote(x,y)...)
@eval ($i)(x::AbstractQuantity, y::Real) = ($i)(promote(x,y)...)
@eval ($i)(x::Real, y::AbstractQuantity) = ($i)(promote(x,y)...)

# promotion might not yield Quantity types
@eval @inline ($j)(x::AbstractQuantity{T1}, y::AbstractQuantity{T2}) where {T1,T2} = ($i)(promote(x,y)...)
Expand All @@ -259,12 +256,12 @@ function isapprox(x::AbstractQuantity, y::AbstractQuantity; kwargs...)
return isapprox(promote(x,y)...; kwargs...)
end

isapprox(x::AbstractQuantity, y::Number; kwargs...) = isapprox(promote(x,y)...; kwargs...)
isapprox(x::Number, y::AbstractQuantity; kwargs...) = isapprox(y, x; kwargs...)
isapprox(x::AbstractQuantity, y::Real; kwargs...) = isapprox(promote(x,y)...; kwargs...)
isapprox(x::Real, y::AbstractQuantity; kwargs...) = isapprox(y, x; kwargs...)

function isapprox(
x::AbstractArray{<:AbstractQuantity{T1,D,U1}},
y::AbstractArray{<:AbstractQuantity{T2,D,U2}};
x::AbstractArray{<:AbstractQuantityOrComplex{T1,D,U1}},
y::AbstractArray{<:AbstractQuantityOrComplex{T2,D,U2}};
rtol::Real=Base.rtoldefault(T1,T2,0),
atol=zero(Quantity{real(T1),D,U1}),
norm::Function=norm,
Expand All @@ -280,10 +277,10 @@ function isapprox(
end

isapprox(x::AbstractArray{S}, y::AbstractArray{T};
kwargs...) where {S <: AbstractQuantity,T <: AbstractQuantity} = false
kwargs...) where {S <: AbstractQuantityOrComplex,T <: AbstractQuantityOrComplex} = false

function isapprox(x::AbstractArray{S}, y::AbstractArray{N};
kwargs...) where {S <: AbstractQuantity,N <: Number}
kwargs...) where {S <: AbstractQuantityOrComplex,N <: Number}
if dimension(N) == dimension(S)
isapprox(map(x->uconvert(NoUnits,x),x),y; kwargs...)
else
Expand All @@ -292,7 +289,7 @@ function isapprox(x::AbstractArray{S}, y::AbstractArray{N};
end

isapprox(y::AbstractArray{N}, x::AbstractArray{S};
kwargs...) where {S <: AbstractQuantity,N <: Number} = isapprox(x,y; kwargs...)
kwargs...) where {S <: AbstractQuantityOrComplex,N <: Number} = isapprox(x,y; kwargs...)

for cmp in [:(==), :isequal]
@eval $cmp(x::AbstractQuantity{S,D,U}, y::AbstractQuantity{T,D,U}) where {S,T,D,U} = $cmp(x.val, y.val)
Expand All @@ -301,10 +298,10 @@ for cmp in [:(==), :isequal]
$cmp(promote(x,y)...)
end

@eval function $cmp(x::AbstractQuantity, y::Number)
@eval function $cmp(x::AbstractQuantity, y::Real)
$cmp(promote(x,y)...)
end
@eval $cmp(x::Number, y::AbstractQuantity) = $cmp(y,x)
@eval $cmp(x::Real, y::AbstractQuantity) = $cmp(y,x)
end

_dimerr(f) = error("$f can only be well-defined for dimensionless ",
Expand Down Expand Up @@ -358,12 +355,18 @@ for (f,r) = ((:trunc, :RoundToZero), (:floor, :RoundDown), (:ceil, :RoundUp))
@eval $f(u::Units, x::AbstractQuantity; kwargs...) = round(u, x, $r; kwargs...)
end

zero(x::AbstractQuantity) = Quantity(zero(x.val), unit(x))
# same as Base methods for Any arguments
# unlike Real methods, do not perform promotion
Base.min(x::AbstractQuantity, y::AbstractQuantity) = ifelse(isless(y, x), y, x)
Base.max(x::AbstractQuantity, y::AbstractQuantity) = ifelse(isless(y, x), x, y)
Base.minmax(x::AbstractQuantity, y::AbstractQuantity) = ifelse(isless(y, x), (y, x), (x, y))

zero(x::AbstractQuantityOrComplex) = Quantity(zero(ustrip(x)), unit(x))
zero(x::AffineQuantity) = Quantity(zero(x.val), absoluteunit(x))
zero(x::Type{<:AbstractQuantity{T}}) where {T} = throw(ArgumentError("zero($x) not defined."))
zero(x::Type{<:AbstractQuantity{T,D}}) where {T,D} = zero(T) * upreferred(D)
zero(x::Type{<:AbstractQuantity{T,D,U}}) where {T,D,U<:ScalarUnits} = zero(T)*U()
zero(x::Type{<:AbstractQuantity{T,D,U}}) where {T,D,U<:AffineUnits} = zero(T)*absoluteunit(U())
zero(x::Type{<:AbstractQuantityOrComplex{T}}) where {T} = throw(ArgumentError("zero($x) not defined."))
zero(x::Type{<:AbstractQuantityOrComplex{T,D}}) where {T,D} = zero(T) * upreferred(D)
zero(x::Type{<:AbstractQuantityOrComplex{T,D,U}}) where {T,D,U<:ScalarUnits} = zero(T)*U()
zero(x::Type{<:AbstractQuantityOrComplex{T,D,U}}) where {T,D,U<:AffineUnits} = zero(T)*absoluteunit(U())

function zero(x::AbstractArray{T}) where T<:AbstractQuantity
if isconcretetype(T)
Expand Down Expand Up @@ -479,18 +482,18 @@ typemin(x::AbstractQuantity{T}) where {T} = typemin(T)*unit(x)
typemax(::Type{<:AbstractQuantity{T,D,U}}) where {T,D,U} = typemax(T)*U()
typemax(x::AbstractQuantity{T}) where {T} = typemax(T)*unit(x)

Base.literal_pow(::typeof(^), x::AbstractQuantity, ::Val{v}) where {v} =
Quantity(Base.literal_pow(^, x.val, Val(v)),
Base.literal_pow(::typeof(^), x::AbstractQuantityOrComplex, ::Val{v}) where {v} =
Quantity(Base.literal_pow(^, ustrip(x), Val(v)),
Base.literal_pow(^, unit(x), Val(v)))

# All of these are needed for ambiguity resolution
^(x::AbstractQuantity, y::Integer) = Quantity((x.val)^y, unit(x)^y)
^(x::AbstractQuantity, y::Integer) = Quantity(ustrip(x)^y, unit(x)^y)
@static if VERSION v"1.8.0-DEV.501"
Base.@constprop(:aggressive, ^(x::AbstractQuantity, y::Rational) = Quantity((x.val)^y, unit(x)^y))
Base.@constprop(:aggressive, ^(x::AbstractQuantityOrComplex, y::Rational) = Quantity(ustrip(x)^y, unit(x)^y))
else
^(x::AbstractQuantity, y::Rational) = Quantity((x.val)^y, unit(x)^y)
^(x::AbstractQuantityOrComplex, y::Rational) = Quantity(ustrip(x)^y, unit(x)^y)
end
^(x::AbstractQuantity, y::Real) = Quantity((x.val)^y, unit(x)^y)
^(x::AbstractQuantityOrComplex, y::Real) = Quantity(ustrip(x)^y, unit(x)^y)

Base.rand(r::Random.AbstractRNG, ::Random.SamplerType{<:AbstractQuantity{T,D,U}}) where {T,D,U} =
rand(r, T) * U()
Expand Down
11 changes: 7 additions & 4 deletions src/range.jl
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ Base._range(a::T, step::T, ::Nothing, len::Integer) where {T<:Quantity} =
else
Base._rangestyle(OrderStyle(a), ArithmeticStyle(a), a, step, len)
end
Base._range(a::Quantity{<:Real}, step::Quantity{<:AbstractFloat}, ::Nothing, len::Integer) =
Base._range(a::Quantity, step::Quantity{<:AbstractFloat}, ::Nothing, len::Integer) =
_unitful_start_step_length(float(a), step, len)
Base._range(a::Quantity{<:AbstractFloat}, step::Quantity{<:Real}, ::Nothing, len::Integer) =
_unitful_start_step_length(a, float(step), len)
Expand Down Expand Up @@ -86,7 +86,10 @@ end
colon(start::A, step, stop::C) where {A<:Real,C<:Quantity} = colonstartstop(start,step,stop)
colon(start::A, step, stop::C) where {A<:Quantity,C<:Real} = colonstartstop(start,step,stop)
colon(a::T, b::Quantity, c::T) where {T<:Real} = colon(promote(a,b,c)...)
colon(start::Quantity{<:Real}, step, stop::Quantity{<:Real}) =
colon(a::T, b::Quantity, c::T) where {T<:AbstractFloat} = colon(promote(a,b,c)...) # disambiguation
colon(start::Quantity, step, stop::Quantity) =
colon(promote(start, step, stop)...)
colon(start::Quantity, step::AbstractFloat, stop::Quantity) = # disambiguation
colon(promote(start, step, stop)...)

# promotes start and stop
Expand Down Expand Up @@ -123,8 +126,8 @@ colon(start::T, step::T, stop::T) where {T<:Quantity{<:Base.IEEEFloat}} =
colon(ustrip(start), ustrip(step), ustrip(stop)) * unit(T) # This will always return a StepRangeLen

# two-argument colon
colon(start, stop::Quantity) = _unitful_start_stop(start, stop)
colon(start::Quantity, stop) = _unitful_start_stop(start, stop)
colon(start::Real, stop::Quantity) = _unitful_start_stop(start, stop)
colon(start::Quantity, stop::Real) = _unitful_start_stop(start, stop)
colon(start::Quantity, stop::Quantity) = _unitful_start_stop(start, stop)
function _unitful_start_stop(start, stop)
dimension(start) != dimension(stop) && throw(DimensionError(start, stop))
Expand Down
10 changes: 7 additions & 3 deletions src/types.jl
Original file line number Diff line number Diff line change
Expand Up @@ -148,7 +148,9 @@ The type parameter `T` represents the numeric backing type. The type parameters
Of course, the dimensions follow from the units, but the type parameters are
kept separate to permit convenient dispatch on dimensions.
"""
abstract type AbstractQuantity{T,D,U} <: Number end
abstract type AbstractQuantity{T,D,U} <: Real end

const AbstractQuantityOrComplex{T,D,U} = Union{AbstractQuantity{T,D,U},Complex{<:AbstractQuantity{T,D,U}}}

"""
struct Quantity{T,D,U} <: AbstractQuantity{T,D,U}
Expand All @@ -159,10 +161,12 @@ The type parameter `T` represents the numeric backing type. The type parameters
"""
struct Quantity{T,D,U} <: AbstractQuantity{T,D,U}
val::T
Quantity{T,D,U}(v::Number) where {T,D,U} = new{T,D,U}(v)
Quantity{T,D,U}(v::Real) where {T,D,U} = new{T,D,U}(v)
Quantity{T,D,U}(v::Quantity) where {T,D,U} = convert(Quantity{T,D,U}, v)
end

Quantity{T,D,U}(v::Complex) where {T,D,U} = complex(Quantity{real(T),D,U}(real(v)), Quantity{real(T),D,U}(imag(v)))

# Field-only constructor
Quantity{<:Any,D,U}(val::Number) where {D,U} = Quantity{typeof(val),D,U}(val)

Expand Down Expand Up @@ -233,7 +237,7 @@ struct LogInfo{N,B,P} end
Abstract supertype of [`Unitful.Level`](@ref) and [`Unitful.Gain`](@ref). It is only
used in promotion to put levels and gains onto a common log scale.
"""
abstract type LogScaled{L<:LogInfo} <: Number end
abstract type LogScaled{L<:LogInfo} <: Real end

const RealOrRealQuantity = Union{Real, AbstractQuantity{<:Real}}

Expand Down
4 changes: 4 additions & 0 deletions src/utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -55,6 +55,8 @@ true
```
"""
@inline ustrip(x::Number) = x / unit(x)
@inline ustrip(x::Real) = x
@inline ustrip(x::Complex) = complex(ustrip(real(x)), ustrip(imag(x)))
@inline ustrip(x::Quantity) = ustrip(x.val)
@inline ustrip(x::Missing) = missing

Expand Down Expand Up @@ -116,7 +118,9 @@ true
```
"""
@inline unit(x::AbstractQuantity{T,D,U}) where {T,D,U} = U()
@inline unit(x::Complex{T}) where {T} = unit(T)
@inline unit(::Type{<:AbstractQuantity{T,D,U}}) where {T,D,U} = U()
@inline unit(::Type{<:Complex{T}}) where {T} = unit(T)


"""
Expand Down
Loading

0 comments on commit 8a45d7b

Please sign in to comment.