diff --git a/Project.toml b/Project.toml index c4f1857..d1b6578 100644 --- a/Project.toml +++ b/Project.toml @@ -3,18 +3,30 @@ uuid = "e1ba4f0e-776d-440f-acd9-e1d2e9742647" authors = ["Sheehan Olver "] version = "0.1.13" +[weakdeps] +ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" +Static = "aedffcd0-7271-4cad-89d0-dc628f76c6d3" + +[extensions] +InfinitiesForwardDiffExt = "ForwardDiff" +InfinitiesStaticExt = "Static" + [compat] Aqua = "0.8" Base64 = "1" +ForwardDiff = "1" JET = "0.9, 0.10, 0.11, 0.12" +Static = "1.4" Test = "1" julia = "1.10" [extras] Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595" Base64 = "2a0f44e3-6c83-55bd-87e4-b1978d98bd5f" +ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b" +Static = "aedffcd0-7271-4cad-89d0-dc628f76c6d3" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" [targets] -test = ["Aqua", "Base64", "JET", "Test"] +test = ["Aqua", "Base64", "ForwardDiff", "JET", "Static", "Test"] diff --git a/README.md b/README.md index 2fca520..10d445e 100644 --- a/README.md +++ b/README.md @@ -16,6 +16,86 @@ This Julia package is used to represent infinities, including: Note that we subtype based on interfaces, rather than strict mathematical definitions. For example, `ℵ₀ isa Integer` as `Integer` is often used to represent the size of a set or vector. Similarly, `∞ isa Real`. +## Extending `RealInfinity` + +To add another representation of positive or negative infinity, subtype `RealInfinity` +and define `Base.signbit`: return `false` for positive infinity and `true` for negative infinity. + +```julia +struct SignedInfinity <: RealInfinity + negative::Bool +end +Base.signbit(inf::SignedInfinity) = inf.negative + +SignedInfinity(true) == -∞ # true +SignedInfinity(true)^2 === +∞ # true +``` + +Read the sign from your representation in `signbit`. Do not define it as `x < 0`, +because the inherited comparisons call `signbit` and would cause infinite recursion. +You do not need to implement arithmetic, comparisons, floating-point conversion, or +hashing: your subtype inherits these operations from `RealInfinity`. Any additional +fields you store are ignored when comparing or hashing values. Calling `zero` on your +type or an instance returns `0.0`. Calling `one` or `oneunit` returns `1.0`. + +A subtype may represent only one sign. Results need not retain its concrete type or +metadata: negation and powers may return `PositiveInfinity` or `NegativeInfinity`. +Constructors and representation-preserving conversions are the subtype's responsibility. +Specialize standard Base operations when representation preservation is needed. + +## Static.jl integration + +Loading [Static.jl](https://github.com/SciML/Static.jl) enables an optional package extension. +The values `∞`, `+∞`, `-∞`, `InfiniteCardinal{k}()`, and `NotANumber()` are already +fully determined by their types, so they need no separate static representation: + +```julia +using Infinities, Static + +static(∞) === ∞ # true +is_static(typeof(-∞)) === True() # true +known(typeof(ℵ₀)) === ℵ₀ # true +Static.lt(static(2), ∞) === True() # true +static(2) + ℵ₀ === ℵ₀ # true +``` + +Mixed arithmetic and comparisons use the same rules as the corresponding ordinary +numbers, including undefined results and NaN propagation. Base operations retain their +ordinary return types: `one(∞)` is `1` and `isinf(∞)` is `true`. Use `static` on a result +or Static.jl's comparison functions when a static result type is needed. + +`ComplexInfinity` stores its direction as a value, not in its type. Its arithmetic +accepts static operands, but `is_static(ComplexInfinity)` is `False()` and +`static(im*∞)` is unsupported. + +## ForwardDiff.jl integration + +Loading [ForwardDiff.jl](https://github.com/JuliaDiff/ForwardDiff.jl) enables an optional +extension for mixed `+`, `-`, `*`, `/`, `mod`, `rem`, `min`, `max`, and comparisons +of dual numbers with `∞`, `+∞`, `-∞`, and `NotANumber()`. Primal values retain +Infinities' scalar results and exceptions, without converting infinities to floats. +Derivative rules use the same scalar arithmetic, so a zero tangent multiplied by an +infinity becomes `NotANumber()`. Bounded `mod` and `rem` preserve the dividend's +tangents, while an undefined remainder has undefined tangents. + +```julia +using Infinities, ForwardDiff + +ForwardDiff.derivative(input -> input * ∞, 2.0) # +∞ +ForwardDiff.derivative(input -> input + ∞, 2.0) # 1.0 +ForwardDiff.derivative(input -> input / ∞, 2.0) # 0.0 +``` + +`Dual(∞)` preserves the infinity. Explicitly requesting a floating scalar type, such +as `Dual{Nothing, Float64}(∞)`, converts it. Mixed symbolic results can use `Real` +as the dual's scalar parameter, with a corresponding loss of type specialization. + +This is not general support for arbitrary symbolic differentiation: ForwardDiff +operations such as unary negation can reject heterogeneous symbolic partials, and +powers remain subject to the existing scalar and ForwardDiff limitations. There is +no floating-point fallback. The extension does not define mixed dual-number operations +for `InfiniteCardinal` or `ComplexInfinity`. + ## Similar packages This package is meant to eventually replace [Infinity.jl](https://github.com/cjdoris/Infinity.jl) and the definitions of `∞` in [InfiniteArrays.jl](https://github.com/JuliaArrays/InfiniteArrays.jl). We do not yet support Infinity.jl's notions of `InfExtendedReal` but we hope to add this soon. diff --git a/ext/InfinitiesForwardDiffExt.jl b/ext/InfinitiesForwardDiffExt.jl new file mode 100644 index 0000000..971ec2e --- /dev/null +++ b/ext/InfinitiesForwardDiffExt.jl @@ -0,0 +1,51 @@ +module InfinitiesForwardDiffExt + +using Infinities: Infinity, RealInfinity, NotANumber +import ForwardDiff: Dual, Partials, value, partials + +function symbolic_dual(::Dual{Tag}, primal, tangents::NTuple{N, Any}) where {Tag, N} + Dual{Tag, Real, N}(primal, Partials{N, Real}(tangents)) +end + +for Typ in (Infinity, RealInfinity, NotANumber) + for op in (:^, :min, :max, :(==), :isequal, :<, :<=, :isless) + @eval Base.$op(dual::Dual, inf::$Typ) = invoke($op, Tuple{Dual, Real}, dual, inf) + @eval Base.$op(inf::$Typ, dual::Dual) = invoke($op, Tuple{Real, Dual}, inf, dual) + end + + for (op, forward, reverse) in ((:+, :tangent, :tangent), + (:-, :tangent, :(-tangent)), (:*, :(tangent * inf), :(tangent * inf)), + (:/, :(tangent / inf), :(-(primal / value(dual)) * tangent))) + @eval function Base.$op(dual::Dual, inf::$Typ) + primal = $op(value(dual), inf) + symbolic_dual(dual, primal, map(tangent -> $forward, Tuple(partials(dual)))) + end + @eval function Base.$op(inf::$Typ, dual::Dual) + primal = $op(inf, value(dual)) + symbolic_dual(dual, primal, map(tangent -> $reverse, Tuple(partials(dual)))) + end + end + + for op in (:mod, :rem) + @eval function Base.$op(dual::Dual, inf::$Typ) + primal = $op(value(dual), inf) + symbolic_dual(dual, primal, map(tangent -> isnan(primal) ? NotANumber() : tangent, Tuple(partials(dual)))) + end + @eval Base.$op(inf::$Typ, dual::Dual) = + symbolic_dual(dual, $op(inf, value(dual)), map(_ -> NotANumber(), Tuple(partials(dual)))) + end + + @eval begin + Dual(inf::$Typ) = Dual{Nothing}(inf, ()) + Dual{Tag}(inf::$Typ) where {Tag} = Dual{Tag}(inf, ()) + Dual{Tag, Value}(inf::$Typ) where {Tag, Value} = Dual{Tag, Value, 0}(inf) + Dual{Tag, Value, N}(inf::$Typ) where {Tag, Value, N} = + Dual{Tag, Value, N}(convert(Value, inf), zero(Partials{N, Value})) + Dual{Tag}(inf::$Typ, tangents::Partials{N, Value}) where {Tag, N, Value} = + Dual{Tag, Real, N}(inf, Partials{N, Real}(Tuple(tangents))) + Dual{Tag}(inf::Value, tangents::Partials{N, Value}) where {Tag, N, Value<:$Typ} = + Dual{Tag, Real, N}(inf, Partials{N, Real}(Tuple(tangents))) + end +end + +end \ No newline at end of file diff --git a/ext/InfinitiesStaticExt.jl b/ext/InfinitiesStaticExt.jl new file mode 100644 index 0000000..32acfb0 --- /dev/null +++ b/ext/InfinitiesStaticExt.jl @@ -0,0 +1,37 @@ +module InfinitiesStaticExt + +using Infinities: Infinity, PositiveInfinity, NegativeInfinity, InfiniteCardinal, NotANumber +using Infinities: AllInfinities, AllRealInfinities, IntegerInfinities, OrderedInfinities +using Infinities: ComplexInfinity, RealInfinity +using Static: Static, dynamic, StaticNumber, StaticInteger, StaticFloat64, StaticInt, True + +for Typ in (Infinity, PositiveInfinity, NegativeInfinity, NotANumber) + @eval begin + Static.static(x::$Typ) = x + Static.is_static(::Type{$Typ}) = True() + Static.known(::Type{$Typ}) = $Typ() + end +end + +Static.static(x::InfiniteCardinal) = x +Static.is_static(::Type{InfiniteCardinal{N}}) where {N} = True() +Static.known(::Type{InfiniteCardinal{N}}) where {N} = InfiniteCardinal{N}() + +for (ops, Types, StaticTypes) in ( + ((:+, :-, :*, :/, :(==)), (AllInfinities,), (StaticNumber,)), + ((:isequal,), (NotANumber,), (StaticNumber,)), + ((:div, :fld, :cld, :divrem), (IntegerInfinities,), (StaticNumber,)), + ((:<, :<=, :>, :>=), (OrderedInfinities,), (StaticInteger, StaticFloat64)), + ((:isless,), (AllRealInfinities, InfiniteCardinal, NotANumber), (StaticInteger, StaticFloat64)), + ((:min, :max), (OrderedInfinities, NotANumber), (StaticInteger, StaticFloat64)), + ((:mod,), (IntegerInfinities, NotANumber), (StaticNumber,)), + ((:rem,), (IntegerInfinities, NotANumber), (StaticInteger, StaticFloat64)), + ((:*,), (InfiniteCardinal,), (StaticInt{0},)), +), op in ops, Typ in Types, StaticTyp in StaticTypes + @eval Base.$op(x::$Typ, y::$StaticTyp) = $op(x, dynamic(y)) + @eval Base.$op(x::$StaticTyp, y::$Typ) = $op(dynamic(x), y) +end + +Base.:^(x::RealInfinity, y::StaticNumber) = x^dynamic(y) + +end \ No newline at end of file diff --git a/src/Infinities.jl b/src/Infinities.jl index 27bc2d8..6179c26 100644 --- a/src/Infinities.jl +++ b/src/Infinities.jl @@ -3,7 +3,7 @@ module Infinities import Base: angle, isone, iszero, isinf, isfinite, isnan, isreal, abs, one, oneunit, zero, isless, isequal, inv, +, -, *, /, ^, ==, <, ≤, >, ≥, fld, cld, div, mod, rem, divrem, min, max, sign, signbit, isapprox, - string, show, promote_rule, convert, getindex, tryparse, conj, + string, show, promote_rule, convert, getindex, tryparse, conj, complex, isinteger, round, floor, ceil, trunc, float, Bool, Integer @@ -44,7 +44,7 @@ _convert(::Type{T}, ::Infinity) where {T<:Real} = convert(T, Inf)::T (::Type{T})(x::Infinity) where {T<:Real} = _convert(T, x) sign(y::Infinity) = 1 -angle(x::Infinity) = 0 +angle(x::Infinity) = 0.0 signbit(::Infinity) = false one(::Type{Infinity}) = 1 @@ -53,10 +53,29 @@ oneunit(::Infinity) = 1 zero(::Infinity) = 0 zero(::Type{Infinity}) = 0 +""" + RealInfinity <: Real + +Represent a signed real infinity by subtyping `RealInfinity` and implementing `Base.signbit`. + +Every instance must represent exactly positive or negative infinity. Define +`Base.signbit(x::YourInfinity)::Bool` directly, without relying on comparisons that +use `signbit`. Additional fields do not affect numeric equality or hashing. + +Inherited operations follow the built-in signed infinities' value semantics, but +need not preserve the concrete type or its metadata. A subtype need not represent +both signs. Finite identities are `0.0` and `1.0`, including `zero`, `one`, and +`oneunit` called on the type. Define constructors and specialize Base operations +separately when representation preservation is needed. + +Use `RealInfinity(negative::Bool)` to construct from a sign bit. This is not a numeric +conversion: `convert(RealInfinity, negative)` throws `InexactError`. +""" abstract type RealInfinity <: Real end struct PositiveInfinity <: RealInfinity end struct NegativeInfinity <: RealInfinity end +signbit(x::RealInfinity) = throw(ArgumentError("$(typeof(x)) must implement Base.signbit")) signbit(::PositiveInfinity) = false signbit(::NegativeInfinity) = true one(::RealInfinity) = 1.0 @@ -65,6 +84,7 @@ RealInfinity() = PositiveInfinity() RealInfinity(::Infinity) = PositiveInfinity() RealInfinity(x::RealInfinity) = x RealInfinity(x::Bool) = ifelse(x, NegativeInfinity(), PositiveInfinity()) +convert(::Type{RealInfinity}, x::Bool) = throw(InexactError(:convert, RealInfinity, x)) PositiveInfinity(::Infinity) = PositiveInfinity() # otherwise the generic `(::Type{T})(::Infinity) where T<:Real` would route through `Inf` _convert(::Type{Float16}, x::RealInfinity) = sign(x)*Inf16 @@ -85,83 +105,131 @@ show(io::IO, y::RealInfinity) = print(io, string(y)) Base.to_index(i::RealInfinity) = convert(Integer, i) -one(::Type{RealInfinity}) = 1.0 -oneunit(::Type{RealInfinity}) = 1.0 +one(::Type{<:RealInfinity}) = 1.0 +oneunit(::Type{<:RealInfinity}) = 1.0 oneunit(::RealInfinity) = 1.0 zero(::RealInfinity) = 0.0 -zero(::Type{RealInfinity}) = 0.0 +zero(::Type{<:RealInfinity}) = 0.0 ####### # ComplexInfinity ####### -# angle is π*a where a is (false==0) and (true==1) - """ -ComplexInfinity(signbit) + ComplexInfinity(turns::UInt64) + ComplexInfinity(; halfturns::Real = 0) + +Construct an infinity in the complex plane, pointing in a direction held as a count of +`2^-64` turns. -represents an infinity in the complex plane with the angle -specified by `π * signbit`. The use of the name `signbit` is -for consistency with `RealInfinity`. +The count wraps at a full turn, so the stored `UInt64` and the directions are bijective. +`0x0` points along the positive real axis. Values increase counterclockwise. +`0x8000000000000000` points along the negative real axis. Use `reinterpret(UInt64, x)` to +read out the exact value. + +Pass a `UInt64` directly to construct from a direction count. `convert(ComplexInfinity, count)` +throws `InexactError`, because the finite numeric value of the count is not an infinity. + +The complex plane carries no order, so `isless`, `<`, `≤`, `min` and `max` have no method +here, just as they have none for `Complex`. A direction along the real axis is no exception. +Convert it with `RealInfinity` to compare it. + +Multiplying by `∞` takes the direction from the other operand, which usually reads better +than naming an angle: + + im*∞ # 0 + ∞*im + (1+im)*∞ # ∞ + ∞*im + exp(im*π/4)*∞ # the same direction again + +Those forms and the `halfturns` keyword go through `angle`, so they round. It is exact on +the axes and at a quarter turn, but `exp(im*π/8)*∞` lands 256 counts past an eighth turn. +Provide the `UInt64` when you have an off-axis value where accuracy matters. """ -struct ComplexInfinity{T<:Real} <: Number - signbit::T +struct ComplexInfinity <: Number + turns::UInt64 + ComplexInfinity(turns::UInt64) = new(turns) end -ComplexInfinity{T}() where T = ComplexInfinity(zero(T)) -ComplexInfinity() = ComplexInfinity{Bool}() -ComplexInfinity{T}(::Infinity) where T<:Real = ComplexInfinity{T}() +# Half of the `typemax(UInt64) + 1` counts of a full turn, a number that would overflow. +const _HALFTURN = typemax(UInt64) ÷ 2 + 1 # the negative real axis +# The signs of the parts along each of the eight rays, the axes and diagonals. +const _EIGHTH = _HALFTURN ÷ 4 +const _RAYPARTS = ((1, 0), (1, 1), (0, 1), (-1, 1), (-1, 0), (-1, -1), (0, -1), (1, -1)) +# `_turns` and `_halfturns` are inverse: half turns in, count out, and back again. +# `mod` returns 2 itself for a tiny negative angle, since 2 + x rounds back to 2. A full +# turn is the direction zero. Testing `== 2` rather than `< 2` still lets `NaN` throw. +@inline _turns(halfturns::Real) = round(UInt64, (h = mod(halfturns, 2); h == 2 ? zero(h) : h) * 0x1p63) +# Scaling by 2^63 needs 63 bits beyond the numerator, so `Int128` has room for any `Int64`. +_turns(halfturns::Rational) = round(BigInt, mod(halfturns, 2) * big(2)^63) % UInt64 +_turns(halfturns::Rational{<:Base.BitInteger64}) = + round(Int128, mod(halfturns, 2) * Int128(2)^63) % UInt64 +# `Base` puts an angle in `(-π, π]`, so past the half turn the count reads as negative. +@inline _halfturns(x::ComplexInfinity) = + x.turns == _HALFTURN ? 1.0 : reinterpret(Int64, x.turns) / 0x1p63 + +ComplexInfinity(; halfturns::Real = 0) = ComplexInfinity(_turns(halfturns)) ComplexInfinity(::Infinity) = ComplexInfinity() -ComplexInfinity{T}(x::RealInfinity) where T<:Real = ComplexInfinity{T}(signbit(x)) -ComplexInfinity(x::RealInfinity) = ComplexInfinity(signbit(x)) -ComplexInfinity{T}(x::ComplexInfinity) where T<:Real = ComplexInfinity(T(x.signbit)) # ambiguity fix +ComplexInfinity(x::RealInfinity) = ComplexInfinity(_directionof(x)) +ComplexInfinity(x::ComplexInfinity) = x + +signbit(y::ComplexInfinity) = y.turns == _HALFTURN +isreal(y::ComplexInfinity) = iszero(y.turns) || signbit(y) -signbit(y::ComplexInfinity) = mod(y.signbit, 2) == 1 +# `Base` converts a `Complex` to a `Real` the same way, and throws the same error off the axis. +RealInfinity(x::ComplexInfinity) = isreal(x) ? RealInfinity(signbit(x)) : + throw(InexactError(:RealInfinity, RealInfinity, x)) -convert(::Type{ComplexInfinity{T}}, ::Infinity) where T = ComplexInfinity{T}() convert(::Type{ComplexInfinity}, ::Infinity) = ComplexInfinity() -convert(::Type{ComplexInfinity{T}}, x::RealInfinity) where T = ComplexInfinity{T}(x) convert(::Type{ComplexInfinity}, x::RealInfinity) = ComplexInfinity(x) +convert(::Type{ComplexInfinity}, x::UInt64) = throw(InexactError(:convert, ComplexInfinity, x)) -sign(y::ComplexInfinity{<:Integer}) = mod(y.signbit, 2) == 0 ? 1 : -1 -sign(y::ComplexInfinity) = cispi(y.signbit) -angle(x::ComplexInfinity) = π*x.signbit +sign(y::ComplexInfinity) = cispi(_halfturns(y)) +angle(x::ComplexInfinity) = _halfturns(x) * π abs(::ComplexInfinity) = ∞ -conj(y::ComplexInfinity{<:Integer}) = y # an integer factor points along the real axis -conj(y::ComplexInfinity) = ComplexInfinity(mod(-y.signbit, 2)) +conj(y::ComplexInfinity) = ComplexInfinity(-y.turns) # An exact zero has to stay finite, `Inf * 0` being a `NaN`. @inline _ray(c) = iszero(c) ? c : copysign(Inf, c) # `Complex` reaches only the eight rays of its two saturating parts, so the direction lands on the nearest of them. function float(x::ComplexInfinity) - s, c = sincospi(x.signbit) + s, c = sincospi(_halfturns(x)) complex(_ray(c), _ray(s)) end -show(io::IO, x::ComplexInfinity) = print(io, "exp($(x.signbit)*im*π)∞") +# The rays print as `Base` prints an infinite `Complex`. Off the rays the readable form +# names an angle, which recovers most counts but not all, so it is used only where reading +# it back gives the same direction. +function show(io::IO, x::ComplexInfinity) + k, offset = divrem(x.turns, _EIGHTH) + if iszero(offset) + r, i = _RAYPARTS[k + 1] + return print(io, ("-∞", "0", "∞")[r + 2], (" - ∞*im", " + 0im", " + ∞*im")[i + 2]) + end + h = _halfturns(x) + _directionof(cispi(h)) == x.turns ? print(io, "cispi($h)∞") : + print(io, "ComplexInfinity(", repr(x.turns), ")") +end -one(::Type{<:ComplexInfinity}) = one(ComplexF64) -oneunit(::Type{<:ComplexInfinity}) = oneunit(ComplexF64) +one(::Type{ComplexInfinity}) = one(ComplexF64) +oneunit(::Type{ComplexInfinity}) = oneunit(ComplexF64) oneunit(::ComplexInfinity) = oneunit(ComplexF64) zero(::ComplexInfinity) = zero(ComplexF64) -zero(::Type{<:ComplexInfinity}) = zero(ComplexF64) +zero(::Type{ComplexInfinity}) = zero(ComplexF64) # `isequal` implies equal hashes, so the infinities have to hash like the float # infinities they compare equal to. The interface requires implementing `hash(x, h::UInt)`. Base.hash(::Infinity, h::UInt)::UInt = hash(Inf, h) -Base.hash(::PositiveInfinity, h::UInt)::UInt = hash(Inf, h) -Base.hash(::NegativeInfinity, h::UInt)::UInt = hash(-Inf, h) +Base.hash(x::RealInfinity, h::UInt)::UInt = hash(signbit(x) ? -Inf : Inf, h) -# Equality of ComplexInfinity is equality of the angle, hence so is the hash. +# The two real directions have to hash like the real infinities they compare equal to. function Base.hash(x::ComplexInfinity, h::UInt)::UInt - θ = angle(x) - θ == angle(PositiveInfinity()) && return hash(Inf, h) - θ == angle(NegativeInfinity()) && return hash(-Inf, h) - hash(ComplexInfinity, hash(θ, h)) + iszero(x.turns) && return hash(Inf, h) + x.turns == _HALFTURN && return hash(-Inf, h) + hash(ComplexInfinity, hash(x.turns, h)) end diff --git a/src/algebra.jl b/src/algebra.jl index ba07fe1..59f37a0 100644 --- a/src/algebra.jl +++ b/src/algebra.jl @@ -1,4 +1,5 @@ @inline infpromote(x, y) = Base._promote(x, y) +@inline infpromote(x::Bool, y::Union{Infinity, ComplexInfinity}) = (x, y) @inline infpromote(x::ExtendedComplex, y::AllInfinities) = (x, ComplexInfinity(y)) @inline infpromote(x::ExtendedComplex, y::ComplexInfinity) = Base._promote(x, y) @inline infpromote(x::Real, ::InfiniteCardinal) = (x, ∞) @@ -15,29 +16,36 @@ +(::Infinity) = RealInfinity() -(::Infinity) = RealInfinity(true) -(y::RealInfinity) = RealInfinity(!signbit(y)) --(y::ComplexInfinity{B}) where B<:Integer = sign(y) == 1 ? ComplexInfinity(one(B)) : ComplexInfinity(zero(B)) --(y::ComplexInfinity) = ComplexInfinity(mod(y.signbit + 1, 2)) +-(y::ComplexInfinity) = ComplexInfinity(y.turns ⊻ _HALFTURN) +(x::InfiniteCardinal) = x -(::InfiniteCardinal) = -∞ # addition -@inline _sb(x) = signbit(x) -@inline _sb(x::Complex) = angle(x)/π # overloading `signbit` causes type piracy -@inline _sb(x::ComplexInfinity) = x.signbit # the whole angle, not just its sign - @inline toinf(x) = RealInfinity(signbit(x)) -# The field counts half turns, so the radians of `angle` have to be scaled. -@inline toinf(x::Complex) = ComplexInfinity(_sb(x)) +@inline toinf(x::Complex) = ComplexInfinity(_directionof(x)) @inline toinf(x::ComplexInfinity) = x -@inline _infadd(x, y) = angle(x) == angle(y) ? y : NotANumber() +# The undefined value that matches the operands, as `Base` returns `NaN` or `NaN + NaN*im`. +@inline _undefined(x, y) = + x isa ExtendedComplex || y isa ExtendedComplex ? ComplexNotANumber : NotANumber() + +@inline _infadd(x, y) = _directionof(x) == _directionof(y) ? y : + x isa ComplexInfinity || y isa ComplexInfinity ? _rayadd(x, y) : _undefined(x, y) + +# On the eight rays each part is infinite or exactly zero, so `Base` adds part by part. +@inline function _rayadd(x, y) + (kx, ox), (ky, oy) = divrem(_directionof(x), _EIGHTH), divrem(_directionof(y), _EIGHTH) + iszero(ox | oy) || return _undefined(x, y) + (rx, ix), (ry, iy) = _RAYPARTS[kx + 1], _RAYPARTS[ky + 1] + rx * ry < 0 || ix * iy < 0 ? _undefined(x, y) : toinf(complex(sign(rx + ry), sign(ix + iy))) +end @inline __add(x, y::AllInfinities) = isinf(x) ? _infadd(toinf(x), y) : y @inline __add(x::Integer, y::InfiniteCardinal) = max(x, y) # A `NaN` argument makes the result undefined. Types with no `NaN` fold the test away. -@inline _add(x, y) = isnan(x) ? NotANumber() : __add(infpromote(x, y)...) +@inline _add(x, y) = isnan(x) ? _undefined(x, y) : __add(infpromote(x, y)...) +(x::Number, y::AllInfinities) = _add(x, y) +(x::AllInfinities, y::Number) = _add(y, x) @@ -54,15 +62,19 @@ # multiplication -@inline __mul(x, y::AllInfinities) = RealInfinity(_sb(x) ⊻ _sb(y)) -@inline __mul(x, y::ComplexInfinity) = ComplexInfinity(_sb(x) + _sb(y)) -@inline __mul(x, y::ComplexInfinity{Bool}) = ComplexInfinity(_sb(x) ⊻ _sb(y)) -@inline __mul(x::Complex, y::ComplexInfinity{Bool}) = ComplexInfinity(_sb(x) + _sb(y)) +# The count of the direction a value points in. `_turns` instead reads its argument as a +# number of half turns. +@inline _directionof(x::Real) = signbit(x) ? _HALFTURN : zero(UInt64) +@inline _directionof(x::Complex) = _turns(angle(x) / π) # overloading `signbit` causes type piracy +@inline _directionof(x::ComplexInfinity) = x.turns + +@inline __mul(x, y::AllInfinities) = RealInfinity(signbit(x) ⊻ signbit(y)) +@inline __mul(x, y::ComplexInfinity) = ComplexInfinity(_directionof(x) + _directionof(y)) @inline __mul(x::Integer, y::InfiniteCardinal) = x > 0 ? y : throw(ArgumentError("Cannot multiply $x * $y")) @inline function _mul(x, y) - isnan(x) && return NotANumber() - iszero(x) && return NotANumber() + isnan(x) && return _undefined(x, y) + iszero(x) && return _undefined(x, y) __mul(infpromote(x, y)...) end @@ -80,10 +92,12 @@ end # division # `\` needs nothing of its own, `Base` defining it as `y / x`. @inline _div(x, y) = x * inv(y) +# `Base` divides two `Integer`s in floating point. +@inline _div(x::Integer, y::InfiniteCardinal) = float(x) * inv(y) /(x::AllInfinities, y::Number) = _div(x, y) /(x::Number, y::AllInfinities) = _div(x, y) -/(x::AllInfinities, y::AllInfinities) = NotANumber() +/(x::AllInfinities, y::AllInfinities) = _undefined(x, y) # mod @inline function _mod(x::Real, y::IntegerInfinities) @@ -105,15 +119,14 @@ divrem(x::Real, y::IntegerInfinities) = (div(x, y), rem(x, y)) divrem(x::IntegerInfinities, y::Real) = (div(x, y), rem(x, y)) divrem(x::IntegerInfinities, y::IntegerInfinities) = (div(x, y), rem(x, y)) -# fld, cld, div -_divinf(x) = isnan(x) ? NotANumber() : zero(x) -_fldinf(x) = isnan(x) ? NotANumber() : signbit(x) ? -one(x) : zero(x) -_cldinf(x) = isnan(x) ? NotANumber() : signbit(x) ? zero(x) : one(x) -div(x::Real, ::IntegerInfinities) = _divinf(x) -fld(x::Real, ::IntegerInfinities) = _fldinf(x) -cld(x::Real, ::IntegerInfinities) = _cldinf(x) +# A finite numerator over an infinite divisor rounds to zero in every mode. As in `Base`, that +# zero takes its sign from the divisor and its type from the numerator. +_fcdinf(x, y) = isnan(x) || isinf(x) ? NotANumber() : signbit(y) ? -zero(x) : zero(x) +div(x::Real, y::IntegerInfinities) = _fcdinf(x, y) +fld(x::Real, y::IntegerInfinities) = _fcdinf(x, y) +cld(x::Real, y::IntegerInfinities) = _fcdinf(x, y) -_inffcd(x, y) = isnan(y) ? NotANumber() : signbit(y) ? -x : x +_inffcd(x, y) = isnan(y) || isinf(y) ? NotANumber() : signbit(y) ? -x : x for OP in (:fld,:cld,:div) @eval begin $OP(x::IntegerInfinities, y::Real) = _inffcd(x, y) @@ -124,6 +137,7 @@ end # power # Although the base implementation can cover these cases, it can change overtime and yield inconsistent results. # ref: https://github.com/JuliaMath/Infinities.jl/actions/runs/19993302836/ +_infpow(x::RealInfinity, p) = _infpow(RealInfinity(signbit(x)), p) _infpow(::PositiveInfinity, p) = isnan(p) ? NotANumber() : ifelse(iszero(p), one(p), ifelse(p > 0, +∞, +zero(p))) function _infpow(x::NegativeInfinity, p) isnan(p) && return NotANumber() @@ -149,8 +163,8 @@ for op in (:+, :-, :*, :/, :^, :div, :fld, :cld, :mod, :rem, :min, :max) @eval $op(::$Typ, y::NotANumber) = y end for Typ in NotANumberComplexRivals - @eval $op(::NotANumber, ::$Typ) = complex(NotANumber(), NotANumber()) - @eval $op(::$Typ, ::NotANumber) = complex(NotANumber(), NotANumber()) + @eval $op(::NotANumber, ::$Typ) = ComplexNotANumber + @eval $op(::$Typ, ::NotANumber) = ComplexNotANumber end @eval $op(x::NotANumber, ::NotANumber) = x end @@ -158,6 +172,10 @@ for Typ in NotANumberRivals @eval divrem(x::NotANumber, ::$Typ) = (x, x) @eval divrem(::$Typ, y::NotANumber) = (y, y) end +for Typ in NotANumberComplexRivals + @eval divrem(::NotANumber, ::$Typ) = (ComplexNotANumber, ComplexNotANumber) + @eval divrem(::$Typ, ::NotANumber) = (ComplexNotANumber, ComplexNotANumber) +end divrem(x::NotANumber, ::NotANumber) = (x, x) # `Base` has its own `^(::Number, ::Integer)`, which a literal exponent also routes through. ^(x::NotANumber, ::Integer) = x diff --git a/src/ambiguities.jl b/src/ambiguities.jl index 4a80a6b..a81a8f9 100644 --- a/src/ambiguities.jl +++ b/src/ambiguities.jl @@ -1,10 +1,9 @@ for Typ in (Base.TwicePrecision, AbstractChar, Complex) - @eval begin - RealInfinity(x::$Typ) = throw(MethodError(RealInfinity, x)) - ComplexInfinity{T}(x::$Typ) where T<:Real = ComplexInfinity(T(x)) - end + @eval RealInfinity(x::$Typ) = throw(MethodError(RealInfinity, x)) end -ComplexInfinity{T}(x::ComplexInfinity{T}) where T<:Real = x +# `Base` builds a `TwicePrecision` of any type by adding its two halves, which for two opposed +# directions is undefined. +ComplexInfinity(x::Base.TwicePrecision) = ComplexInfinity(halfturns = Float64(x)) for Typ in (Rational, BigInt, BigFloat) for (op, fop) in ((:<, :_lt), (:≤, :_le)) @@ -39,10 +38,8 @@ for Typ in (Rational, ) @eval rem(x::$Typ, ::IntegerInfinities) = x for op in (:fld, :cld, :div) @eval $op(x::InfiniteCardinal, y::$Typ) = _inffcd(x, y) + @eval $op(x::$Typ, y::IntegerInfinities) = _fcdinf(x, y) end - @eval div(x::$Typ, ::IntegerInfinities) = _divinf(x) - @eval fld(x::$Typ, ::IntegerInfinities) = _fldinf(x) - @eval cld(x::$Typ, ::IntegerInfinities) = _cldinf(x) end divrem(x::BigInt, y::IntegerInfinities) = (div(x, y), rem(x, y)) diff --git a/src/cardinality.jl b/src/cardinality.jl index bd71de2..acad061 100644 --- a/src/cardinality.jl +++ b/src/cardinality.jl @@ -23,7 +23,7 @@ isone(::InfiniteCardinal) = false signbit(::InfiniteCardinal) = false sign(::InfiniteCardinal) = 1 -angle(::InfiniteCardinal) = 0 +angle(::InfiniteCardinal) = 0.0 abs(a::InfiniteCardinal) = a zero(::InfiniteCardinal) = 0 zero(::Type{<:InfiniteCardinal}) = 0 @@ -41,6 +41,9 @@ function Integer(x::ComplexInfinity) ℵ₀ end +# Every cardinal is a positive real infinity, so it points along the positive real axis. +ComplexInfinity(::InfiniteCardinal) = ComplexInfinity() + Base.to_index(::Union{Infinity,InfiniteCardinal{0}}) = ℵ₀ Base.to_shape(::Union{Infinity,InfiniteCardinal{0}}) = ℵ₀ diff --git a/src/compare.jl b/src/compare.jl index a3bb812..601ccb5 100644 --- a/src/compare.jl +++ b/src/compare.jl @@ -16,6 +16,10 @@ _isinf(x::Number, y::AllInfinities) = isinf(x) && _angle(x) == angle(y) # On the real line the direction is a comparison against zero. # `signbit(y)` is constant, so the branch folds away and the check becomes a single instruction. _isinf(x::Real, y::AllRealInfinities) = isinf(x) && (signbit(y) ? x < zero(x) : x > zero(x)) +# A direction in the plane is decided by the count, which is exact where an angle in a +# `Float64` is not: the count has 64 bits and the angle has 53. +_isinf(x::Number, y::ComplexInfinity) = isinf(x) && _directionof(x) == y.turns +_isinf(x::ComplexInfinity, y::AllRealInfinities) = x.turns == _directionof(y) # NotANumber # Undefined compares false against everything, itself included, as `NaN` does. @@ -67,14 +71,13 @@ isapprox(::NotANumber, ::NotANumber; kwargs...) = false # `isless` is the sort order. `NaN` sorts after every other value, infinities included. isless(x::AllRealInfinities, y::AllRealInfinities) = signbit(x) && !signbit(y) @generated isless(::InfiniteCardinal{N}, ::InfiniteCardinal{M}) where {N,M} = :($(isless(N, M))) -# The leading `signbit` call discards its result. It is there to reject a non-real `Number`. -for Typ in (Number, Real, AbstractFloat) +for Typ in (Real, AbstractFloat) @eval begin - isless(x::AllRealInfinities, y::$Typ) = (signbit(y); isnan(y) || signbit(x) && y ≠ -∞) - isless(x::$Typ, y::AllRealInfinities) = (signbit(x); !isnan(x) && !signbit(y) && x ≠ ∞) + isless(x::AllRealInfinities, y::$Typ) = isnan(y) || signbit(x) && y ≠ -∞ + isless(x::$Typ, y::AllRealInfinities) = !isnan(x) && !signbit(y) && x ≠ ∞ end end -for Typ in (Number, Real, AbstractFloat, AllRealInfinities) +for Typ in (Real, AbstractFloat, AllRealInfinities) @eval begin isless(::InfiniteCardinal, x::$Typ) = isnan(x) isless(x::$Typ, y::InfiniteCardinal) = isless(x, ∞) || isless(ℵ₀, y) @@ -92,9 +95,9 @@ isless(::InfiniteCardinal{0}, ::InfiniteCardinal{0}) = false for (op, fop) in ((:max, :_max), (:min, :_min), (:<, :_lt), (:≤, :_le)) for Typ in (Real, ) @eval begin - $op(x::AllInfinities, y::$Typ) = $fop(x, y) - $op(x::$Typ, y::AllInfinities) = $fop(x, y) + $op(x::OrderedInfinities, y::$Typ) = $fop(x, y) + $op(x::$Typ, y::OrderedInfinities) = $fop(x, y) end end - @eval $op(x::AllInfinities, y::AllInfinities) = $fop(x, y) + @eval $op(x::OrderedInfinities, y::OrderedInfinities) = $fop(x, y) end diff --git a/src/interface.jl b/src/interface.jl index 7bcabbf..1d670b7 100644 --- a/src/interface.jl +++ b/src/interface.jl @@ -1,7 +1,9 @@ const AllInfinities = Union{Infinity, RealInfinity, ComplexInfinity, InfiniteCardinal} -const AllRealInfinities = Union{Infinity, RealInfinity, ComplexInfinity{<:Integer}} -const IntegerInfinities = Union{Infinity, RealInfinity, ComplexInfinity{<:Integer}, InfiniteCardinal} -const ExtendedComplex{T} = Union{Complex{T}, ComplexInfinity{T}} +const AllRealInfinities = Union{Infinity, RealInfinity} +const IntegerInfinities = Union{Infinity, RealInfinity, InfiniteCardinal} +# The infinities that lie on the real line and so have a place in the numeric ordering. +const OrderedInfinities = Union{Infinity, RealInfinity, InfiniteCardinal} +const ExtendedComplex = Union{Complex, ComplexInfinity} iszero(::AllInfinities) = false isinf(::AllInfinities) = true @@ -25,6 +27,7 @@ const NotANumberRivals = (Number, Real, AbstractFloat, AbstractIrrational, AllIn InfiniteCardinal) # A complex operand makes the undefined result complex, as it does over the floats. const NotANumberComplexRivals = (Complex, Complex{Bool}, ComplexInfinity) +const ComplexNotANumber = complex(NotANumber(), NotANumber()) # `InfiniteCardinal` is absent because `Base` already returns `true` for it through `Integer`. isinteger(::Union{Infinity, RealInfinity, ComplexInfinity}) = false @@ -35,10 +38,17 @@ round(x::Union{AllInfinities, NotANumber}, ::RoundingMode; kwargs...) = x # `Infinity` is positive, so it has no common type with `NegativeInfinity` (as is already the case for `PositiveInfinity`). promote_rule(::Type{Infinity}, ::Type{PositiveInfinity}) = PositiveInfinity -promote_rule(::Type{Infinity}, ::Type{ComplexInfinity{T}}) where T = ComplexInfinity{T} -promote_rule(::Type{<:RealInfinity}, ::Type{ComplexInfinity{T}}) where T = ComplexInfinity{T} -promote_rule(::Type{ComplexInfinity{T}}, ::Type{<:RealInfinity}) where T<:Integer = ComplexInfinity{T} -promote_rule(::Type{ComplexInfinity{T}}, ::Type{ComplexInfinity{S}}) where {T, S} = ComplexInfinity{promote_type(T, S)} +promote_rule(::Type{Infinity}, ::Type{ComplexInfinity}) = ComplexInfinity +promote_rule(::Type{<:RealInfinity}, ::Type{ComplexInfinity}) = ComplexInfinity + +# An infinite part makes the whole number infinite, so `complex` gives the direction it points in. +complex(x::IntegerInfinities) = ComplexInfinity(x) +complex(x::ComplexInfinity) = x +complex(::Type{<:IntegerInfinities}) = ComplexInfinity +complex(::Type{ComplexInfinity}) = ComplexInfinity +complex(x::IntegerInfinities, y::Real) = x + im*y +complex(x::Real, y::IntegerInfinities) = x + im*y +complex(x::IntegerInfinities, y::IntegerInfinities) = x + im*y function tryparse(::Type{NegativeInfinity}, s::AbstractString) i = findfirst(!isspace, s) diff --git a/test/runtests.jl b/test/runtests.jl index d71a2ce..f397823 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -2,6 +2,7 @@ using Infinities, Base64, Test import Infinities: Infinity, AllInfinities, _isinf using Aqua, JET +using Static: Static "An `AbstractString` indexed by character position, so that byte arithmetic on indices is invalid." struct CharString <: AbstractString @@ -15,6 +16,14 @@ Base.isvalid(s::CharString, i::Integer) = 1 ≤ i ≤ ncodeunits(s) Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], i + 1) : nothing @testset "∞" begin + @testset "Boolean arithmetic" begin + for inf in (∞, +∞, -∞, ℵ₀, ℵ₁, NotANumber(), ComplexInfinity(), im*∞), value in (false, true), + (args, expected) in (((inf, value), (inf, Int(value))), ((value, inf), (Int(value), inf))), + op in (+, -, *) + + @test isequal(op(args...), op(expected...)) + end + end @testset "∞" begin @test ∞ ≠ 1 @test 1 ≠ ∞ @@ -54,7 +63,7 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test !signbit(∞) @test sign(∞) ≡ 1 - @test angle(∞) ≡ 0 + @test angle(∞) ≡ 0.0 @test string(∞) == stringmime("text/plain", ∞) == "∞" @@ -79,12 +88,8 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test div(∞, 2) ≡ ∞ @test fld(∞, 2) ≡ ∞ @test cld(∞, 2) ≡ ∞ - @test div(2, ∞) ≡ 0 - @test fld(2, ∞) ≡ 0 - @test cld(2, ∞) ≡ 1 - @test div(-2, ∞) ≡ 0 - @test fld(-2, ∞) ≡ -1 - @test cld(-2, ∞) ≡ 0 + @test div(2, ∞) ≡ fld(2, ∞) ≡ cld(2, ∞) ≡ 0 + @test div(-2, ∞) ≡ fld(-2, ∞) ≡ cld(-2, ∞) ≡ 0 @test mod(2,∞) ≡ 2 @test div(∞,∞) isa NotANumber @test fld(∞,∞) isa NotANumber @@ -92,6 +97,31 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test mod(∞,∞) isa NotANumber @test mod(∞,2) isa NotANumber @test_throws ArgumentError mod(-2,∞) + + for op in (div, fld, cld), x in (0, 2, -2, 0.0, -0.0, 1.5, -1.5, 2//3, -2//3) + @test op(x, ∞) ≡ op(x, +∞) ≡ op(x, ℵ₀) ≡ zero(x) + @test op(x, -∞) ≡ -zero(x) + end + + for op in (div, fld, cld), x in (big(3), big(-3.0)) + @test op(x, ∞) == zero(x) && op(x, -∞) == zero(x) + end + + for op in (div, fld, cld), x in (Inf, -Inf, NaN), inf in (∞, +∞, -∞, ℵ₀) + @test op(x, inf) ≡ NotANumber() + @test op(inf, x) ≡ NotANumber() + end + + # before 1.13, `Base` itself gives e.g. `fld(-1.5, Inf) === NaN` + VERSION ≥ v"1.13" && for op in (div, fld, cld), x in (0.0, -0.0, 1.5, -1.5), + (inf, flt) in ((∞, Inf), (+∞, Inf), (-∞, -Inf), (ℵ₀, Inf)) + + @test op(x, inf) ≡ op(x, flt) + end + + for op in (div, fld, cld), x in (0, 2, -2), inf in (∞, +∞, -∞, ℵ₀) + @test op(x, inf) ≡ 0 + end end @testset "convert" begin @@ -234,23 +264,70 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test convert(Float32, -∞) ≡ Float32(-∞) ≡ -Inf32 @test convert(Float16, -∞) ≡ Float16(-∞) ≡ -Inf16 @test convert(BigFloat, -∞)::BigFloat == BigFloat(-∞)::BigFloat == -BigFloat(Inf) + for negative in (false, true) + inf = RealInfinity(negative) + @test signbit(inf) === negative + @test convert(RealInfinity, inf) === inf + values = RealInfinity[inf] + @test_throws InexactError convert(RealInfinity, negative) + @test_throws InexactError values[1] = negative + @test values[1] === inf + end end @test Base.to_index(RealInfinity()) ≡ ℵ₀ end @testset "ComplexInfinity" begin + # every spelling of the positive real axis is the same value, a cardinal included @test ComplexInfinity(∞) ≡ convert(ComplexInfinity, ∞) ≡ ComplexInfinity() ≡ - ComplexInfinity(false) ≡ ComplexInfinity{Bool}(∞) ≡ ComplexInfinity{Bool}(RealInfinity()) ≡ ComplexInfinity{Bool}(ComplexInfinity()) - - @test convert(ComplexInfinity{Bool}, ∞) ≡ convert(ComplexInfinity, ∞) ≡ ComplexInfinity() - @test convert(ComplexInfinity{Bool}, -∞) ≡ convert(ComplexInfinity, -∞) ≡ -ComplexInfinity() + ComplexInfinity(0x0000000000000000) ≡ ComplexInfinity(RealInfinity()) ≡ + ComplexInfinity(ComplexInfinity()) ≡ ComplexInfinity(ℵ₀) + + @test convert(ComplexInfinity, -∞) ≡ -ComplexInfinity() + for turns in (UInt64(0), UInt64(1), 0x8000000000000000, typemax(UInt64)) + inf = ComplexInfinity(turns) + @test reinterpret(UInt64, inf) === turns + @test convert(ComplexInfinity, inf) === inf + values = ComplexInfinity[inf] + @test_throws InexactError convert(ComplexInfinity, turns) + @test_throws InexactError values[1] = turns + @test values[1] === inf + end + # one direction is one value, however it is spelled + @test ComplexInfinity(halfturns = -0.5) ≡ ComplexInfinity(halfturns = 1.5) ≡ -im*∞ + @test ComplexInfinity(halfturns = 1) ≡ ComplexInfinity(halfturns = 3) ≡ ComplexInfinity(-∞) + @test isreal(ComplexInfinity()) && isreal(-ComplexInfinity()) && !isreal((1+im)*∞) + # a rational direction converts without a float step + @test reinterpret(UInt64, ComplexInfinity(halfturns = 2//3)) ≡ 0x5555555555555555 + # a numerator too wide for `Int128` takes the `BigInt` route to the same count + @test ComplexInfinity(halfturns = big(1)//3) ≡ ComplexInfinity(halfturns = 1//3) + @test ComplexInfinity(0x4000000000000000) ≡ ComplexInfinity(halfturns = 1//2) ≡ im*∞ + # `mod` rounds a hair below the axis up to a full turn, which has no count of its own + @test ComplexInfinity(halfturns = -1e-300) ≡ ComplexInfinity(halfturns = 2.0) ≡ ComplexInfinity() + @test complex(1.0, -1e-17)*∞ ≡ ComplexInfinity() + # no count names a direction that is not one, so the conversion has to refuse + for h in (NaN, Inf, -Inf) + @test_throws InexactError ComplexInfinity(halfturns = h) + end + # the count runs forwards, `angle` reports it on `Base`'s branch of `(-π, π]` + for h in (0.0, 0.25, 0.5, 1.0, -0.25, -0.5, -0.75) + @test angle(ComplexInfinity(halfturns = h)) ≡ h*π + @test complex(cospi(h), sinpi(h))*∞ == ComplexInfinity(halfturns = h) + end + # off the axes the angle has to be rounded to reach a count + @test reinterpret(UInt64, exp(im*π/8)*∞) - reinterpret(UInt64, ComplexInfinity(halfturns = 1//8)) ≡ + 0x0000000000000100 + @test angle(exp(im*0.3)*∞) ≈ angle(∞*exp(im*0.3)) ≈ 0.3 + # the count is finer than an angle in a `Float64`, so equality has to read the count + @test ComplexInfinity(0x7fffffffffffffff) ≠ -ComplexInfinity() + @test im*∞ * ComplexInfinity(0x0000000000000001) ≠ im*∞ @test isinf(ComplexInfinity()) @test !isfinite(ComplexInfinity()) @test promote(∞, RealInfinity(), ComplexInfinity()) ≡ ntuple(_ -> ComplexInfinity(), 3) - @test promote_type(Infinity, ComplexInfinity{Bool}) == promote_type(RealInfinity, ComplexInfinity{Bool}) == ComplexInfinity{Bool} + @test promote_type(Infinity, ComplexInfinity) == promote_type(RealInfinity, ComplexInfinity) == ComplexInfinity @test ComplexInfinity(∞) == ∞ @@ -268,27 +345,54 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test ComplexInfinity() + ∞ ≡ ComplexInfinity() + RealInfinity() ≡ ∞ + ComplexInfinity() ≡ RealInfinity() + ComplexInfinity() ≡ ComplexInfinity() - @test ComplexInfinity(true) + ComplexInfinity(true) == ComplexInfinity(true) - @test ComplexInfinity(false) + ComplexInfinity(false) == ComplexInfinity(false) - @test ComplexInfinity(true)+1 == ComplexInfinity(true) - @test ComplexInfinity(false)+1 == ComplexInfinity(false) + @test ComplexInfinity(-∞) + ComplexInfinity(-∞) == ComplexInfinity(-∞) + @test ComplexInfinity() + ComplexInfinity() == ComplexInfinity() + @test ComplexInfinity(-∞)+1 == ComplexInfinity(-∞) + @test ComplexInfinity()+1 == ComplexInfinity() # An infinite summand reaches `_infadd` through `toinf`, which has to give half turns @test complex(Inf, 0.0) + ∞ ≡ ComplexInfinity() @test complex(-Inf, 0.0) + (-∞) ≡ -ComplexInfinity() @test complex(0.0, Inf) + im*∞ ≡ im*∞ @test complex(0.0, -Inf) + (-im*∞) ≡ -im*∞ - @test complex(0.0, Inf) + ∞ ≡ NotANumber() + @test complex(0.0, Inf) + ∞ ≡ (1+im)*∞ + # on the eight rays the sum goes part by part, as the `Complex` sum does + rays = (∞, (1+im)*∞, im*∞, (-1+im)*∞, -∞, (-1-im)*∞, -im*∞, (1-im)*∞) + for x in rays, y in rays + z = float(x) + float(y) + @test isequal(x + y, !isnan(z) ? z*∞ : + z isa Complex ? complex(NotANumber(), NotANumber()) : NotANumber()) + end + @test ∞ + im*∞ ≡ im*∞ + ∞ ≡ complex(∞, ∞) ≡ (1+im)*∞ + @test exp(0.1im)*∞ + ∞ ≡ ComplexInfinity(0x7fffffffffffffff) + (-ComplexInfinity()) ≡ + complex(NotANumber(), NotANumber()) # two infinite parts are the only way an infinite `Complex` points off the axes for (z, inf) in ((complex(Inf, Inf), (1+im)*∞), (complex(-Inf, Inf), (-1+im)*∞), (complex(-Inf, -Inf), (-1-im)*∞), (complex(Inf, -Inf), (1-im)*∞)) @test z + inf ≡ inf end + # an infinite part makes the whole number infinite, pointing where the `Complex` would + @test complex(∞) ≡ complex(+∞) ≡ complex(ℵ₀) ≡ complex(∞, 0) ≡ complex(+∞, -2.5) ≡ ComplexInfinity() + @test complex(-∞) ≡ complex(-∞, 0) ≡ complex(-∞, -0.0) ≡ -ComplexInfinity() + @test complex(0, ∞) ≡ complex(-3, +∞) ≡ im*∞ + @test complex(0.0, -∞) ≡ -im*∞ + for (x, y, inf) in ((∞, ∞, (1+im)*∞), (0.5∞, ∞, (1+im)*∞), (-∞, +∞, (-1+im)*∞), + (-∞, -∞, (-1-im)*∞), (ℵ₀, -∞, (1-im)*∞), (∞, Inf, (1+im)*∞), + (-Inf, 0.5∞, (-1+im)*∞)) + @test complex(x, y) ≡ inf + end + for nan in (NaN, NotANumber()) + @test complex(∞, nan) ≡ complex(nan, -∞) ≡ complex(NotANumber(), NotANumber()) + end + @test complex(ComplexInfinity()) ≡ ComplexInfinity() && complex(im*∞) ≡ im*∞ + @test complex(Infinity) ≡ complex(RealInfinity) ≡ complex(PositiveInfinity) ≡ + complex(InfiniteCardinal{0}) ≡ complex(ComplexInfinity) ≡ ComplexInfinity + @test ∞ * ComplexInfinity() ≡ RealInfinity() * ComplexInfinity() ≡ ComplexInfinity() * ∞ ≡ ComplexInfinity() * RealInfinity() ≡ ComplexInfinity() - @test 2.0im*∞ ≡ ∞*2.0im ≡ 2.0im * RealInfinity() ≡ RealInfinity() * 2.0im ≡ ComplexInfinity(1/2) + @test 2.0im*∞ ≡ ∞*2.0im ≡ 2.0im * RealInfinity() ≡ RealInfinity() * 2.0im ≡ im*∞ @test 2ComplexInfinity() ≡ ComplexInfinity()*2 ≡ ComplexInfinity() # a factor gives the direction it actually has, so rescaling moves it once it rounds @test 4*(0.3+0.1im)*∞ ≡ (0.3+0.1im)*∞ @@ -301,68 +405,74 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test Inf == ComplexInfinity() @test ComplexInfinity() == Inf - @test isless(-ComplexInfinity(), ComplexInfinity()) - @test isless(5, ComplexInfinity()) - @test !isless(ComplexInfinity(), 5) + # the complex plane carries no order, so these are undefined as they are for `Complex` + for op in (isless, <, ≤, >, ≥, min, max), + y in (5, Inf, -Inf, NaN, ∞, +∞, -∞, ℵ₀, ComplexInfinity(), -ComplexInfinity(), (1+im)*∞) - @test 5 < ComplexInfinity() && 5 ≤ ComplexInfinity() - @test !(ComplexInfinity() < 5) && !(ComplexInfinity() ≤ 5) - @test 5 > -ComplexInfinity() && 5 ≥ -ComplexInfinity() - @test ComplexInfinity() > 5 && ComplexInfinity() ≥ 5 + @test_throws MethodError op(ComplexInfinity(), y) + @test_throws MethodError op(y, ComplexInfinity()) + end + # a direction on the axis converts, as `Real(::Complex)` does + @test RealInfinity(ComplexInfinity()) ≡ +∞ + @test RealInfinity(-ComplexInfinity()) ≡ -∞ + @test_throws InexactError RealInfinity((1+im)*∞) + @test 5 < RealInfinity(ComplexInfinity()) @test 1 + ComplexInfinity() ≡ 1.0 + ComplexInfinity() ≡ ComplexInfinity() + 1 ≡ ComplexInfinity() + 1.0 ≡ ComplexInfinity() @test 5 * ComplexInfinity() ≡ ComplexInfinity() @test (-5) * ComplexInfinity() ≡ -ComplexInfinity() - @test ComplexInfinity(0.25) * ComplexInfinity(0.5) ≡ ComplexInfinity(0.75) - @test ComplexInfinity(0.0) + ComplexInfinity() ≡ ComplexInfinity() + ComplexInfinity(0.0) ≡ ComplexInfinity(0.0) - - @test mod(ComplexInfinity(), 5) ≡ NotANumber() - - @test stringmime("text/plain", ComplexInfinity()) == "exp(false*im*π)∞" - - @testset "min/max" begin - @test min(ComplexInfinity(), -ComplexInfinity()) ≡ -ComplexInfinity() - @test max(ComplexInfinity(), -ComplexInfinity()) ≡ ComplexInfinity() - @test min(ComplexInfinity(), 5) ≡ min(5,ComplexInfinity()) ≡ 5 - @test max(ComplexInfinity(), 5) ≡ max(5,ComplexInfinity()) ≡ ComplexInfinity() - end - - @testset "fld/cld/div" begin - @test div(ComplexInfinity(), 5) ≡ fld(ComplexInfinity(), 5) ≡ ComplexInfinity() - @test div(-ComplexInfinity(),2) ≡ -ComplexInfinity() + @test (1+im)*∞ * (im*∞) ≡ (-1+im)*∞ + @test (2.0+0.0im)*∞ + ComplexInfinity() ≡ ComplexInfinity() + (2.0+0.0im)*∞ ≡ ComplexInfinity() + + @test stringmime("text/plain", ComplexInfinity()) == "∞ + 0im" + @test map(x -> sprint(show, ComplexInfinity(x)), rays) == + ("∞ + 0im", "∞ + ∞*im", "0 + ∞*im", "-∞ + ∞*im", + "-∞ + 0im", "-∞ - ∞*im", "0 - ∞*im", "∞ - ∞*im") + @test sprint(show, ComplexInfinity(halfturns = 0.1)) == "cispi(0.1)∞" + # a count an angle cannot name is shown as itself, so every form reads back + @test sprint(show, ComplexInfinity(0x5555555555555555)) == "ComplexInfinity(0x5555555555555555)" + for x in (ComplexInfinity.(rays)..., ComplexInfinity(0x0ccccccccccccd00), + ComplexInfinity(0x5555555555555555), ComplexInfinity(0xdeadbeefdeadbeef)) + @test Core.eval(@__MODULE__, Meta.parse(sprint(show, x))) ≡ x + end + + @testset "integer operations" begin + # an integer operation needs a real, and `Base` defines none of these for a `Complex` + for op in (div, fld, cld, mod, rem), x in (ComplexInfinity(), (1+im)*∞) + @test_throws MethodError op(x, 5) + @test_throws MethodError op(5, x) + end + @test div(RealInfinity(ComplexInfinity()), 5) ≡ +∞ end - @test signbit(ComplexInfinity(3)) - @test !signbit(ComplexInfinity(100)) + @test signbit(ComplexInfinity(halfturns = 3)) + @test !signbit(ComplexInfinity(halfturns = 100)) # `signbit` returns a `Bool` for every angle, as it does over the reals - @test signbit(ComplexInfinity(1.0)) === signbit(-ComplexInfinity()) === true - @test signbit(ComplexInfinity(0.5)) === signbit(ComplexInfinity()) === false + @test signbit(ComplexInfinity(-∞)) === signbit(-ComplexInfinity()) === true + @test signbit(im*∞) === signbit(ComplexInfinity()) === false @testset "abs/sign/conj/-" begin - @test -ComplexInfinity(0.5) ≡ ComplexInfinity(1.5) - @test -(-ComplexInfinity(0.5)) ≡ ComplexInfinity(0.5) - @test -ComplexInfinity() ≡ ComplexInfinity(true) - @test abs(ComplexInfinity()) ≡ abs(ComplexInfinity(0.5)) ≡ ∞ - @test sign(ComplexInfinity(0.5)) ≡ complex(0.0, 1.0) - @test sign(ComplexInfinity(0.0)) ≡ complex(1.0, 0.0) - @test sign(ComplexInfinity(1.0)) ≡ complex(-1.0, 0.0) - # an integer angle stays on the real line, where the sign is a real ±1 - @test sign(ComplexInfinity(false)) ≡ 1 - @test sign(ComplexInfinity(true)) ≡ -1 - # off the axes conjugation and negation part company: `-ComplexInfinity(0.25)` is `ComplexInfinity(1.25)` - @test conj(ComplexInfinity(0.25)) ≡ ComplexInfinity(1.75) - @test conj(conj(ComplexInfinity(0.25))) ≡ ComplexInfinity(0.25) - @test conj(ComplexInfinity(true)) ≡ ComplexInfinity(true) # the narrow type survives + @test -(im*∞) ≡ -im*∞ + @test -(-(im*∞)) ≡ im*∞ + @test -ComplexInfinity() ≡ ComplexInfinity(-∞) + @test abs(ComplexInfinity()) ≡ abs(im*∞) ≡ ∞ + @test sign(im*∞) ≡ complex(0.0, 1.0) + @test sign(ComplexInfinity()) ≡ complex(1.0, 0.0) + @test sign(ComplexInfinity(-∞)) ≡ complex(-1.0, 0.0) + # conjugation negates the direction, so on the real axis it changes nothing + @test conj((1+im)*∞) ≡ (1-im)*∞ + @test conj(conj((1+im)*∞)) ≡ (1+im)*∞ + @test conj(ComplexInfinity(-∞)) ≡ ComplexInfinity(-∞) end @testset "float" begin - @test float(ComplexInfinity()) ≡ float(ComplexInfinity(0.0)) ≡ complex(Inf, 0.0) - @test float(ComplexInfinity(1/2)) ≡ complex(0.0, Inf) - @test float(ComplexInfinity(1.0)) ≡ float(ComplexInfinity(true)) ≡ complex(-Inf, 0.0) - @test float(ComplexInfinity(-1/2)) ≡ float(ComplexInfinity(3/2)) ≡ complex(0.0, -Inf) + @test float(ComplexInfinity()) ≡ float((2.0+0.0im)*∞) ≡ complex(Inf, 0.0) + @test float(im*∞) ≡ complex(0.0, Inf) + @test float(ComplexInfinity(-∞)) ≡ complex(-Inf, 0.0) + @test float(-im*∞) ≡ float(ComplexInfinity(halfturns = 3/2)) ≡ complex(0.0, -Inf) # `Complex` points along eight rays only, so every other angle collapses onto the nearest - @test float(ComplexInfinity(1/4)) ≡ float(ComplexInfinity(0.3)) ≡ complex(Inf, Inf) + @test float((1+im)*∞) ≡ float(ComplexInfinity(0x1000000000000000)) ≡ complex(Inf, Inf) end end @@ -375,8 +485,11 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @testset "hash" begin infinities = (∞, +∞, -∞, Inf, -Inf, Inf32, -Inf32, Inf16, -Inf16, big(Inf), -big(Inf), - InfiniteCardinal{0}(), ComplexInfinity(false), - ComplexInfinity(true), ComplexInfinity(0.1)) + InfiniteCardinal{0}(), ComplexInfinity(), + ComplexInfinity(-∞), ComplexInfinity(0x1000000000000000), + # counts that an angle in a `Float64` cannot tell apart + ComplexInfinity(0x7fffffffffffffff), + im*∞ * ComplexInfinity(0x0000000000000001)) # isequal must imply equal hashes for a in infinities, b in infinities @@ -417,9 +530,9 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test Base.literal_pow(^, -∞, Val(2)) ≡ (-∞)^2 ≡ +∞ @test Base.literal_pow(^, -∞, Val(-2)) ≡ (-∞)^(-2) ≡ 0.0 - @test Base.literal_pow(^, ComplexInfinity(0.1), Val(0)) ≡ ComplexInfinity(0.1)^0 ≡ 1.0+0.0im - @test Base.literal_pow(^, ComplexInfinity(0.1), Val(1)) ≡ (ComplexInfinity(0.1))^1 ≡ ComplexInfinity(0.1) - @test Base.literal_pow(^, ComplexInfinity(0.1), Val(-1)) ≡ (ComplexInfinity(0.1))^(-1) ≡ 0.0+0.0im + @test Base.literal_pow(^, ComplexInfinity(0x1000000000000000), Val(0)) ≡ ComplexInfinity(0x1000000000000000)^0 ≡ 1.0+0.0im + @test Base.literal_pow(^, ComplexInfinity(0x1000000000000000), Val(1)) ≡ (ComplexInfinity(0x1000000000000000))^1 ≡ ComplexInfinity(0x1000000000000000) + @test Base.literal_pow(^, ComplexInfinity(0x1000000000000000), Val(-1)) ≡ (ComplexInfinity(0x1000000000000000))^(-1) ≡ 0.0+0.0im end @testset "one/zero/oneunit" begin @@ -433,9 +546,9 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], end @testset "isinteger/round" begin - infinities = (∞, +∞, -∞, ℵ₀, ComplexInfinity(), ComplexInfinity(1/4)) + infinities = (∞, +∞, -∞, ℵ₀, ComplexInfinity(), (1+im)*∞) @test !isinteger(∞) && !isinteger(+∞) && !isinteger(-∞) - @test !isinteger(ComplexInfinity()) && !isinteger(ComplexInfinity(1/4)) + @test !isinteger(ComplexInfinity()) && !isinteger((1+im)*∞) @test isinteger(ℵ₀) # an `InfiniteCardinal` is an `Integer` @test ∞ ∉ 1:5 # `in` asks a range for `isinteger` before comparing for f in (round, floor, ceil, trunc), x in infinities @@ -452,16 +565,22 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], # a zero divisor keeps the direction, and its own sign is the one that counts @test ∞ / 0 ≡ ∞ / 0.0 ≡ (-∞) / (-0.0) ≡ +∞ @test (-∞) / 0 ≡ (-∞) / 0.0 ≡ ∞ / (-0.0) ≡ -∞ - @test ComplexInfinity(0.5) / 2 ≡ ComplexInfinity(0.5) + @test im*∞ / 2 ≡ im*∞ # dividing by a complex turns the direction by its angle @test (+∞) / (1+im) ≡ (1-im)*∞ @test 2 / -∞ ≡ -0.0 - @test 2 / ∞ == ∞ \ 2 == 2 / ℵ₀ == 0 # the type follows `inv`, which returns an `Int` for `∞` + @test 2 / ∞ ≡ ∞ \ 2 ≡ 0 # the type follows `inv`, which returns an `Int` for `∞` + # `ℵ₀` is an `Integer`, and `Base` divides two of them in floating point + @test 2 / ℵ₀ ≡ ℵ₀ \ 2 ≡ true / ℵ₀ ≡ 0.0 && -2 / ℵ₀ ≡ -0.0 + @test big(2) / ℵ₀ isa BigFloat && iszero(big(2) / ℵ₀) # `∞` is positive, so the quotient keeps the dividend exact; a signed infinity # needs a float to carry `-0.0` @test (2//3) / ∞ ≡ (2//3) / ℵ₀ ≡ 0//1 @test (2//3) / (+∞) ≡ 0.0 @test ∞ / ∞ isa NotANumber + # a complex operand on either side makes the undefined quotient complex + @test ComplexInfinity() / ∞ ≡ ∞ / ComplexInfinity() ≡ + ComplexInfinity() / ComplexInfinity() ≡ complex(NotANumber(), NotANumber()) @test isnan(NaN / ∞) && isnan(∞ / NaN) end @@ -521,7 +640,7 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], # ℵ₁ points in the same direction as ∞, even though `ℵ₁ == ∞` is false positive = (∞, +∞, ℵ₀, ℵ₁, ComplexInfinity(), Inf, Inf32, Inf16, big(Inf)) negative = (-∞, -ComplexInfinity(), -Inf, -Inf32, -Inf16, -big(Inf)) - imaginary = (ComplexInfinity(0.5), complex(0.0, Inf)) + imaginary = (im*∞, complex(0.0, Inf)) others = (0, 1.5, -2, -1.5, 0.0, -0.0, NaN, NaN32, prevfloat(Inf), nextfloat(-Inf), nextfloat(0.0), prevfloat(-0.0), "∞", "-∞") @@ -560,15 +679,23 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @testset "NaN arithmetic" begin # the result is the package's own undefined value, as `Inf + ∞` is its own infinity - for nan in (NaN, NaN32, NaN16, big(NaN)), - inf in (∞, +∞, -∞, ℵ₀, ComplexInfinity(), -ComplexInfinity()) - + for nan in (NaN, NaN32, NaN16, big(NaN)), inf in (∞, +∞, -∞, ℵ₀) for op in (+, -, *, div, fld, cld) @test op(nan, inf) ≡ op(inf, nan) ≡ NotANumber() end # `mod(inf, x)` discards `x`, so only one order is needed @test mod(nan, inf) ≡ NotANumber() end + for nan in (NaN, NaN32, NaN16, big(NaN)), inf in (ComplexInfinity(), -ComplexInfinity()) + # arithmetic is defined for a complex operand, so the undefined result is complex + for op in (+, -, *) + @test op(nan, inf) ≡ op(inf, nan) ≡ complex(NotANumber(), NotANumber()) + end + # `Base` defines no integer operation for a `Complex`, and neither do we + for op in (div, fld, cld, mod, rem) + @test_throws MethodError op(nan, inf) + end + end for nan in (NaN, NaN32, NaN16, big(NaN)), inf in (+∞, -∞) @test inf^nan ≡ NotANumber() end @@ -577,7 +704,9 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @testset "NotANumber" begin nan = NotANumber() # every operand it can meet, itself included - operands = (nan, 0, 1.5, ∞, +∞, -∞, ℵ₀, ComplexInfinity(), NaN, NaN32) + reals = (nan, 0, 1.5, ∞, +∞, -∞, ℵ₀, NaN, NaN32) + complexes = (ComplexInfinity(), (1+im)*∞, complex(1.0, 2.0), complex(true, false)) + operands = (reals..., complexes...) @test isnan(nan) && !isinf(nan) && !isfinite(nan) && !iszero(nan) && !isone(nan) && !signbit(nan) @test !isinteger(nan) # a `NaN` of any real type is real, and this is the type-independent one @@ -610,16 +739,19 @@ Base.iterate(s::CharString, i::Integer=1) = i ≤ length(s.chars) ? (s.chars[i], @test isnan(BigFloat(nan)) # anything computed from it is undefined again - for op in (+, -, *, /, ^, div, fld, cld, mod, rem, min, max), x in operands + for op in (+, -, *, /, ^, div, fld, cld, mod, rem, min, max), x in reals @test op(nan, x) ≡ op(x, nan) ≡ nan end # a complex operand makes the undefined result complex, as it does over the floats - for op in (+, -, *, /, ^, div, fld, cld, mod, rem, min, max), x in (complex(1.0, 2.0), complex(true, false)) + for op in (+, -, *, /, ^, div, fld, cld, mod, rem, min, max), x in complexes @test op(nan, x) ≡ op(x, nan) ≡ complex(nan, nan) end - for x in operands + for x in reals @test divrem(nan, x) ≡ divrem(x, nan) ≡ (nan, nan) end + for x in complexes + @test divrem(nan, x) ≡ divrem(x, nan) ≡ (complex(nan, nan), complex(nan, nan)) + end @test nan^(1//2) ≡ (1//2)^nan ≡ ℯ^nan ≡ nan @test -nan ≡ +nan ≡ abs(nan) ≡ inv(nan) ≡ sign(nan) ≡ conj(nan) ≡ nan end @@ -700,8 +832,13 @@ end include("test_cardinality.jl") include("test_ambiguity.jl") +include("test_static.jl") + +include("test_real_infinity.jl") @testset "Project quality" begin Aqua.test_all(Infinities) test_package(Infinities) end + +include("test_forwarddiff.jl") diff --git a/test/test_ambiguity.jl b/test/test_ambiguity.jl index 03e2119..bda419b 100644 --- a/test/test_ambiguity.jl +++ b/test/test_ambiguity.jl @@ -5,17 +5,22 @@ @test_throws MethodError RealInfinity(Base.TwicePrecision(1.0)) @test_throws MethodError RealInfinity(im) - @test ComplexInfinity{Float64}(Base.TwicePrecision(1.0)) == ComplexInfinity(1) + @test ComplexInfinity(Base.TwicePrecision(1.0)) ≡ ComplexInfinity(-∞) @test_throws MethodError ComplexInfinity(im) - for inf in (∞,+∞,ℵ₀,ComplexInfinity()) + for inf in (∞,+∞,ℵ₀) @test mod(inf, 1//2) ≡ NotANumber() @test mod(1//2, inf) ≡ 1//2 @test fld(1//2, inf) == 0 - @test cld(1//2, inf) == 1 + @test cld(1//2, inf) == 0 @test div(1//2, inf) == 0 - @test fld(inf, 1//2) ≡ cld(inf, 1//2) ≡ div(inf, 1//2) ≡ inf - @test fld(inf, ∞) ≡ fld(inf, +∞) ≡ fld(inf, ℵ₀) ≡ fld(inf, ComplexInfinity()) ≡ NotANumber() + @test fld(inf, 1//2) ≡ cld(inf, 1//2) ≡ div(inf, 1//2) == inf + @test fld(inf, ∞) ≡ fld(inf, +∞) ≡ fld(inf, ℵ₀) ≡ NotANumber() + end + # `Base` defines no integer operation for a `Complex`, so neither does a `ComplexInfinity` take one + for op in (mod, fld, cld, div) + @test_throws MethodError op(ComplexInfinity(), 1//2) + @test_throws MethodError op(1//2, ComplexInfinity()) end @testset "rational power" begin diff --git a/test/test_cardinality.jl b/test/test_cardinality.jl index db7a516..85ce63c 100644 --- a/test/test_cardinality.jl +++ b/test/test_cardinality.jl @@ -12,7 +12,7 @@ Base.getindex(::InfVector, ::InfiniteCardinal{0}) = 42 @test !isone(ℵ₀) @test !iszero(ℵ₀) @test sign(ℵ₀) ≡ 1 && !signbit(ℵ₀) - @test angle(ℵ₀) ≡ 0 + @test angle(ℵ₀) ≡ 0.0 @test Integer(∞) ≡ convert(Integer,∞) ≡ Integer(ℵ₀) ≡ convert(Integer, ℵ₀) ≡ ℵ₀ @test abs(ℵ₀) ≡ ℵ₀ @test zero(ℵ₀) ≡ zero(InfiniteCardinal{0}) ≡ 0 @@ -156,8 +156,8 @@ Base.getindex(::InfVector, ::InfiniteCardinal{0}) = 42 @test ℵ₀ ÷ 5 ≡ ℵ₀ @test ℵ₀ ÷ ℵ₀ ≡ NotANumber() @test 5 ÷ ℵ₀ ≡ 0 - @test fld(-5, ℵ₀) ≡ -1 - @test cld(5, ℵ₀) ≡ 1 + @test fld(-5, ℵ₀) ≡ 0 + @test cld(5, ℵ₀) ≡ 0 @test mod(ℵ₀,ℵ₀) ≡ NotANumber() @test mod(ℵ₀,6) ≡ NotANumber() @test mod(5,ℵ₀) ≡ 5 diff --git a/test/test_forwarddiff.jl b/test/test_forwarddiff.jl new file mode 100644 index 0000000..44e66fb --- /dev/null +++ b/test/test_forwarddiff.jl @@ -0,0 +1,49 @@ +using ForwardDiff: ForwardDiff, Dual, value, partials + +@testset "ForwardDiff extension" begin + @test isempty(Test.detect_ambiguities(Base.get_extension(Infinities, :InfinitiesForwardDiffExt))) + same = (actual, expected) -> typeof(actual) === typeof(expected) && isequal(actual, expected) + for inf in (∞, +∞, -∞, NotANumber()), Scalar in (Float32, Float64, BigFloat) + for input in (0, 2, -2, Inf, -Inf, NaN), reverse in (false, true), + op in (+, -, *, /, mod, rem, min, max) + primal, tangents = Scalar(input), (Scalar(0), Scalar(1), Scalar(-1)) + dual = Dual(primal, tangents...) + args = reverse ? (inf, dual) : (dual, inf) + if op === mod && !reverse && !isnan(primal) && !isnan(inf) && signbit(primal) != signbit(inf) + @test_throws ArgumentError mod(primal, inf) + @test_throws ArgumentError op(args...) + continue + end + expected = op((reverse ? (inf, primal) : (primal, inf))...) + result = op(args...) + @test same(value(result), expected) + expected_partials = map(tangents) do tangent + op === (+) && return tangent + op === (-) && return reverse ? -tangent : tangent + op === (*) && return tangent * inf + op === (/) && return reverse ? -(expected / primal) * tangent : tangent / inf + op in (mod, rem) && return reverse || isnan(expected) ? NotANumber() : tangent + reference = op((reverse ? (Scalar(inf), dual) : (dual, Scalar(inf)))...) + return partials(reference)[findfirst(isequal(tangent), tangents)] + end + @test all(same.(Tuple(partials(result)), expected_partials)) + end + for op in (isless, isequal, ==, <, <=, >, >=), input in (2, Inf, -Inf, NaN), tangent in (0, 1) + dual = Dual(Scalar(input), Scalar(tangent)) + @test op(dual, inf) === op(dual, float(inf)) + @test op(inf, dual) === op(float(inf), dual) + end + for constructor in (Dual, Dual{Nothing}, Dual{Nothing, Real}, Dual{Nothing, Real, 2}) + @test value(constructor(inf)) === inf + end + for constructor in (Dual{Nothing, Scalar}, Dual{Nothing, Scalar, 2}) + @test same(value(constructor(inf)), Scalar(inf)) && iszero(partials(constructor(inf))) + end + end + for (op, expected) in ((*, +∞), (+, 1.0), (/, 0.0), (min, 1.0), (max, 0.0), (mod, 1.0), (rem, 1.0)) + @test ForwardDiff.derivative(input -> op(input, ∞), 2.0) === expected + end + @test value((Dual(2.0, 1.0) + ∞) + 3.0) === ∞ + @test ForwardDiff.derivative(input -> ForwardDiff.derivative(inner -> inner^2 + ∞, input), 2.0) === 2.0 + @test isequal(ForwardDiff.gradient(input -> input[1] * ∞, [2.0, 3.0]), [+∞, NotANumber()]) +end \ No newline at end of file diff --git a/test/test_real_infinity.jl b/test/test_real_infinity.jl new file mode 100644 index 0000000..54056d2 --- /dev/null +++ b/test/test_real_infinity.jl @@ -0,0 +1,79 @@ +using Infinities, Test + +struct SignedInfinity <: RealInfinity + negative::Bool +end +Base.signbit(inf::SignedInfinity) = inf.negative + +struct FixedInfinity{Negative} <: RealInfinity end +Base.signbit(::FixedInfinity{Negative}) where {Negative} = Negative + +struct MissingSignInfinity <: RealInfinity end + +@testset "RealInfinity interface" begin + for negative in (false, true), custom in (SignedInfinity(negative), FixedInfinity{negative}()) + canonical = RealInfinity(negative) + for operation in (+, -, signbit, sign, angle, abs, abs2, inv, float, + Float16, Float32, Float64, BigFloat, ComplexInfinity, + isinf, isfinite, isnan, iszero, isone, isinteger, isreal, + real, imag, conj, zero, one, oneunit, round, floor, ceil, trunc, repr) + @test isequal(operation(custom), operation(canonical)) + end + for InfinityType in (typeof(custom), typeof(canonical), RealInfinity) + @test zero(InfinityType) === 0.0 + @test one(InfinityType) === oneunit(InfinityType) === 1.0 + @test float(InfinityType) === Float64 + end + @test RealInfinity(custom) === custom + for InfinityType in (typeof(custom), RealInfinity, Real) + @test convert(InfinityType, custom) === custom + end + for seed in (UInt(0), UInt(123)) + @test hash(custom, seed) == hash(canonical, seed) + end + @test length(Set((custom, canonical, Float64(canonical)))) == 1 + @test Dict(canonical => :found)[custom] === :found + + for scalar in (-2, -0.0, 0, 2, 2.0, big(2.0), 2//1, Inf, -Inf, NaN, + +∞, -∞, SignedInfinity(!negative), FixedInfinity{!negative}(), NotANumber()), + operation in (+, -, *, /, ==, isequal, isless, <, <=, min, max, copysign, flipsign, + div, fld, cld, rem, divrem) + @test isequal(operation(custom, scalar), operation(canonical, scalar)) + @test isequal(operation(scalar, custom), operation(scalar, canonical)) + end + for scalar in (-2, -0.0, 0, 2, NaN) + @test isequal(mod(custom, scalar), mod(canonical, scalar)) + if !isnan(scalar) && signbit(scalar) != negative + @test_throws ArgumentError mod(scalar, custom) + else + @test isequal(mod(scalar, custom), mod(scalar, canonical)) + end + end + for exponent in (-3, -2, 0, 2, 3, 2.0, 2//1, big(2.0), Inf, -Inf, NaN, NotANumber()) + if negative && isinf(exponent) + @test_throws DomainError custom^exponent + else + @test isequal(custom^exponent, canonical^exponent) + end + end + for exponent in (0.5, 1//2) + if negative + @test_throws DomainError custom^exponent + else + @test custom^exponent === canonical^exponent + end + end + for exponent in (-2, -1, 0, 1, 2, 3) + @test isequal(Base.literal_pow(^, custom, Val(exponent)), + Base.literal_pow(^, canonical, Val(exponent))) + end + for scalar in (1+im, im*∞), operation in (+, -, *, /, ==, isequal) + @test isequal(operation(custom, scalar), operation(canonical, scalar)) + @test isequal(operation(scalar, custom), operation(scalar, canonical)) + end + end + for operation in (signbit, Float64, hash, repr, inf -> inf < 0) + @test_throws ArgumentError operation(MissingSignInfinity()) + end + @test_throws "MissingSignInfinity must implement Base.signbit" signbit(MissingSignInfinity()) +end \ No newline at end of file diff --git a/test/test_static.jl b/test/test_static.jl new file mode 100644 index 0000000..67bab8c --- /dev/null +++ b/test/test_static.jl @@ -0,0 +1,52 @@ +using Static: Static, static, dynamic, is_static, known, eq, lt, True, False + +@testset "Static extension" begin + @test !isnothing(Base.get_extension(Infinities, :InfinitiesStaticExt)) + @test isempty(Test.detect_ambiguities(Infinities, Static, Base.get_extension(Infinities, :InfinitiesStaticExt))) + singletons = (∞, +∞, -∞, ℵ₀, ℵ₁, InfiniteCardinal{2}(), NotANumber()) + for value in (singletons..., (∞, -∞, ℵ₀)) + for op in (static, dynamic) + @test @inferred(op(value)) === value + end + for input in (value, typeof(value)) + @test @inferred(is_static(input)) === True() + @test @inferred(known(input)) === value + end + end + for (op, reference) in ((eq, ==), (lt, <), (+, +), (*, *)), + first in (∞, -∞, NotANumber(), static(0), static(2), static(NaN), static(Inf)), second in singletons + + @test @inferred(op(first, second)) === static(reference(dynamic(first), second)) + end + @test @inferred((-∞)^static(2)) === +∞ + for Typ in (RealInfinity, ComplexInfinity) + @test is_static(Typ) === False() + @test isnothing(known(Typ)) + end + @test_throws ErrorException static(ComplexInfinity()) + + operators = (+, -, *, /, div, fld, cld, mod, rem, divrem, min, max, isless, ==, <, <=, >, >=, isequal) + for inf in (∞, +∞, -∞, ℵ₀, ℵ₁, NotANumber(), ComplexInfinity(), im*∞), + value in (false, true, 0, 2, -2, 0.0, -0.0, 1.5, -1.5, NaN, Inf, -Inf), + args in ((inf, static(value)), (static(value), inf)) + + unsupported = inf isa ComplexInfinity ? + (div, fld, cld, mod, rem, divrem, isless, min, max, <, <=, >, >=) : () + unbounded = args[2] === inf && inf isa Union{Infinities.Infinity, RealInfinity, InfiniteCardinal} && + !isnan(value) && signbit(value) != signbit(inf) ? (mod,) : () + invalid = inf isa InfiniteCardinal && value isa Integer && value < 0 ? (unbounded..., *) : unbounded + for op in setdiff(operators, unsupported, invalid) + @test isequal(op(args...), op(map(dynamic, args)...)) + end + for (exception, ops) in ((MethodError, unsupported), (ArgumentError, invalid)), op in ops + @test_throws exception op(args...) + end + end + for inf in (+∞, -∞), exponent in (false, true, -2, 0, 2, 3, -1.5, 1.5, NaN, Inf) + if inf isa NegativeInfinity && !isnan(exponent) && !isinteger(exponent) + @test_throws DomainError inf^static(exponent) + else + @test isequal(inf^static(exponent), inf^exponent) + end + end +end \ No newline at end of file