Merge pull request #21 from SymbolicML/fast-fractions

Create `FixedRational` type for faster dimensions

Author: MilesCranmer

src/DynamicQuantities.jl

@@ -3,6 +3,7 @@ module DynamicQuantities
 export Quantity, Dimensions, DimensionError, ustrip, dimension, valid
 export ulength, umass, utime, ucurrent, utemperature, uluminosity, uamount
 
+include("fixed_rational.jl")
 include("types.jl")
 include("utils.jl")
 include("math.jl")

src/fixed_rational.jl

@@ -0,0 +1,51 @@
+"""
+    FixedRational{T,den}
+
+A rational number with a fixed denominator. Significantly
+faster than `Rational{T}`, as it never needs to compute
+the `gcd` apart from when printing.
+"""
+struct FixedRational{T<:Integer,den} <: Real
+    num::T
+    global unsafe_fixed_rational(num::Integer, ::Type{T}, ::Val{den}) where {T,den} = new{T,den}(num)
+end
+
+"""
+    denom(F::FixedRational)
+
+Since `den` can be a different type than `T`, this function
+is used to get the denominator as a `T`.
+"""
+denom(::Type{F}) where {T,den,F<:FixedRational{T,den}} = convert(T, den)
+
+# But, for Val(den), we need to use the same type as at init.
+# Otherwise, we would have type instability.
+val_denom(::Type{F}) where {T,den,F<:FixedRational{T,den}} = Val(den)
+
+Base.eltype(::Type{F}) where {T,den,F<:FixedRational{T,den}} = T
+
+(::Type{F})(x::Integer) where {F<:FixedRational} = unsafe_fixed_rational(x * denom(F), eltype(F), val_denom(F))
+(::Type{F})(x::Rational) where {F<:FixedRational} = unsafe_fixed_rational(widemul(x.num, denom(F)) ÷ x.den, eltype(F), val_denom(F))
+
+Base.:*(l::F, r::F) where {F<:FixedRational} = unsafe_fixed_rational(widemul(l.num, r.num) ÷ denom(F), eltype(F), val_denom(F))
+Base.:+(l::F, r::F) where {F<:FixedRational} = unsafe_fixed_rational(l.num + r.num, eltype(F), val_denom(F))
+Base.:-(l::F, r::F) where {F<:FixedRational} = unsafe_fixed_rational(l.num - r.num, eltype(F), val_denom(F))
+Base.:-(x::F) where {F<:FixedRational} = unsafe_fixed_rational(-x.num, eltype(F), val_denom(F))
+Base.inv(x::F) where {F<:FixedRational} = unsafe_fixed_rational(widemul(denom(F), denom(F)) ÷ x.num, eltype(F), val_denom(F))
+
+Base.:(==)(x::F, y::F) where {F<:FixedRational} = x.num == y.num
+Base.iszero(x::FixedRational) = iszero(x.num)
+Base.isinteger(x::F) where {F<:FixedRational} = iszero(x.num % denom(F))
+Base.convert(::Type{F}, x::Integer) where {F<:FixedRational} = unsafe_fixed_rational(x * denom(F), eltype(F), val_denom(F))
+Base.convert(::Type{F}, x::Rational) where {F<:FixedRational} = F(x)
+Base.convert(::Type{Rational}, x::F) where {F<:FixedRational} = Rational{eltype(F)}(x.num, denom(F))
+Base.convert(::Type{AF}, x::F) where {AF<:AbstractFloat,F<:FixedRational} = convert(AF, x.num) / convert(AF, denom(F))
+Base.round(::Type{T}, x::F) where {T,F<:FixedRational} = div(convert(T, x.num), convert(T, denom(F)), RoundNearest)
+Base.promote(x, y::F) where {F<:FixedRational} = promote(x, convert(Rational, y))
+Base.promote(x::F, y) where {F<:FixedRational} = promote(convert(Rational, x), y)
+Base.show(io::IO, x::F) where {F<:FixedRational} = show(io, convert(Rational, x))
+Base.zero(::Type{F}) where {F<:FixedRational} = unsafe_fixed_rational(0, eltype(F), val_denom(F))
+
+tryrationalize(::Type{F}, x::F) where {F<:FixedRational} = x
+tryrationalize(::Type{F}, x::Union{Rational,Integer}) where {F<:FixedRational} = convert(F, x)
+tryrationalize(::Type{F}, x) where {F<:FixedRational} = unsafe_fixed_rational(round(eltype(F), x * denom(F)), eltype(F), val_denom(F))

src/types.jl

@@ -1,4 +1,4 @@
-const DEFAULT_DIM_TYPE = Rational{Int16}
+const DEFAULT_DIM_TYPE = FixedRational{Int32, 2^4 * 3^2 * 5^2 * 7}
 const DEFAULT_VALUE_TYPE = Float64
 
 """

test/test_unitful.jl

@@ -1,4 +1,5 @@
 import DynamicQuantities
+using DynamicQuantities: DEFAULT_DIM_TYPE, DEFAULT_VALUE_TYPE
 import Unitful
 import Unitful: @u_str
 import Ratios: SimpleRatio
@@ -15,7 +16,7 @@ risapprox(x::Unitful.Quantity, y::Unitful.Quantity; kws...) =
 
 factor_for_preferred_units = 1e-3
 
-for T in [Float16, Float32, Float64], R in [Rational{Int16}, Rational{Int32}, SimpleRatio{Int}, SimpleRatio{SafeInt16}]
+for T in [DEFAULT_VALUE_TYPE, Float16, Float32, Float64], R in [DEFAULT_DIM_TYPE, Rational{Int16}, Rational{Int32}, SimpleRatio{Int}, SimpleRatio{SafeInt16}]
     x = DynamicQuantities.Quantity(T(0.2*factor_for_preferred_units), R, length=1, amount=2, current=-1 // 2, luminosity=2 // 5)
     x_unitful = T(0.2)u"m*mol^2*A^(-1//2)*cd^(2//5)"
 

test/unittests.jl

@@ -6,7 +6,7 @@ using Test
 
 @testset "Basic utilities" begin
 
-    for T in [Float16, Float32, Float64], R in [Rational{Int16}, Rational{Int32}, SimpleRatio{Int}, SimpleRatio{SafeInt16}]
+    for T in [DEFAULT_VALUE_TYPE, Float16, Float32, Float64], R in [DEFAULT_DIM_TYPE, Rational{Int16}, Rational{Int32}, SimpleRatio{Int}, SimpleRatio{SafeInt16}]
         x = Quantity(T(0.2), R, length=1, mass=2.5)
 
         @test typeof(x).parameters[1] == T