Instead of defining a new type, it should be faster to "normalize" complex numbers after every operation. The current, IEEE-based complex numbers not only have various representations of infinity, but also for zero, since they distinguish between +0.0 and -0.0. Thus, the normalization would map all of the four complex zeros to zero (e.g. +0.0+0.0im), and would also choose a single normalized representation of infinity (maybe inf+inf*im, where inf=1.0/0.0).
I would probably be necessary to redefine most complex number arithmetic to take these special cases into account. If you do this for complex numbers, then you may also want to do it for real numbers, just for good measure. Currently, there is +0.0 and -0.0 (and +inf and -inf), and a real number corresponding to an "extended complex number" should probably not distinguish between these either. -erik > On Jan 10, 2015, at 14:37 , Jiahao Chen <jia...@mit.edu> wrote: > > You might be interested in Issue #5234, where I tried to catalogue issues > arising from nonfinite floating-point computations in the ordinary complex > plane (without closure on the Riemann sphere). > > Being able to switch between complex and extended complex (closed at > ComplexInf) would be interesting. I vaguely recall discussing complex > infinity with Mike Innes at some point, but I don't recall if we got anywhere. > > One solution is to wrap ordinary complex numbers in a new type representing > finite extended complex numbers, and create yet another new type representing > complex infinity. Define each function over the exted For each function, and > redefine each function you want to use by either wrapping a finite answer in > the new type, or otherwise checking if the answer is infinite and returning > another new type instead which represents complex infinity. > > #### > > abstract ExtendedComplex{T<:Real} > > immutable ComplexInf{T<:Real} <: ExtendedComplex{T} end > > immutable FiniteExtendedComplex{T<:Real} <: ExtendedComplex{T} > z :: Complex{T} > end > > function ExtendedComplex{T<:Real}(r::T, i::T) > if isfinite(r) && isfinite(i)) > return ExtendedComplex(Complex(r, i)) > else > return ComplexInf{T}() > end > end > > import Base./ > > function /{T}(z:: ComplexInf{T}, w:: FiniteExtendedComplex{T}) > return z > end > > function /{T}(z:: FiniteExtendedComplex{T}, w:: ComplexInf{T}) > return zero(z) > end > > function /{T}(z:: FiniteExtendedComplex{T}, w:: FiniteExtendedComplex{T}) > a = z.z/w.z > return isfinite(a) ? ExtendedComplex(a) : ComplexInf{T} > end > > #### > > With these definitions you can compute as follows: > > julia> FiniteExtendedComplex(2.0, 0.5)/FiniteExtendedComplex(0.5, 2.0) > FiniteExtendedComplex{Float64}(0.47058823529411764 - 0.8823529411764706im) > > julia> FiniteExtendedComplex(2.0, 0.3)/FiniteExtendedComplex(0.0, 0.0) > ComplexInf{Float64}() > > julia> ComplexInf{Float64}()/FiniteExtendedComplex(0.0, 0.0) > ComplexInf{Float64}() > > julia> ComplexInf{Float64}()/ComplexInf{Float64}() #Don't know how to do ∞/∞ > ERROR: `/` has no method matching /(::ComplexInf{Float64}, > ::ComplexInf{Float64}) > > > Thanks, > > Jiahao Chen > Staff Research Scientist > MIT Computer Science and Artificial Intelligence Laboratory > > On Sat, Jan 10, 2015 at 8:55 AM, Ed Scheinerman > <edward.scheiner...@gmail.com> wrote: > Is there a way to have a single complex infinity? This may come at the cost > of computational efficiency I suppose, but I can think of situations where > all of the following give the same result: > > julia> (1+1im)/0 > Inf + Inf*im > > julia> 1im/0 > NaN + Inf*im > > julia> 1/0 + im > Inf + 1.0im > > It would be nice (sometimes) if these were all the same ComplexInf, say. > Perhaps there's an "extended complex numbers" module for this sort of work? > -- Erik Schnetter <schnet...@gmail.com> http://www.perimeterinstitute.ca/personal/eschnetter/ My email is as private as my paper mail. I therefore support encrypting and signing email messages. Get my PGP key from http://pgp.mit.edu/.
signature.asc
Description: Message signed with OpenPGP using GPGMail