zero(a::Real)?

Why is zero() of a Real an Int(0) ?

The example below looks weird, because it’s as if zero(0.3333333333) == 0::Int.

It took me some time to locate the cause of my bug. The bug was that I was getting Vector{Int64} from zero(a_vector_of_floating_point_numbers). I didn’t realize that my original vector was one of Reals and I incorrectly assumed that I were getting Vector{Float64} from zero().

Also, what’s the solution? Normally I avoid “0.0”, not to unnecessarily narrow the range of types that work with my functions. Perhaps zeros(Float64, size(a)) is best, instead of zero(a) ?

func(r) = (r == 0) ? 0 : r/3
a = func.(0.0:1.0:5.0)
display(a)
display(zero(a))

How about this.

julia> func(r) = r == zero(r) ? zero(r) : r/3
func (generic function with 1 method)

julia> a = func.(0.0:1.0:5.0)
6-element Vector{Float64}:
 0.0
 0.3333333333333333
 0.6666666666666666
 1.0
 1.3333333333333333
 1.6666666666666667

julia> zero(a)
6-element Vector{Float64}:
 0.0
 0.0
 0.0
 0.0
 0.0
 0.0
4 Likes

Note that a in your example there is not “a vector of floating point numbers” — it’s a vector that can hold any real-valued numeric type, including Float64 but also Int64 and Rational and Bool and more.

This is happening because your func(r) sometimes returns an integer — 0 — and sometimes returns a Float64, so the broadcasted result must be a general array that can handle both. That’s where I’d try to place the solution: fix that type instability like @greatpet suggests and then you’ll end up with a concretely typed a::Vector{Float64}.

But I know that wasn’t your question. I’d guess the rationale behind zero(Real) === 0 would be because integers promote relatively seamlessly with other numeric types.

7 Likes

Regarding the first question (“Why is zero() of a Real an Int(0)?”) I have nothing to add to @mbauman’s answer. Regarding the second:

Besides @greatpet’s proposal to avoid type instability, I’d also ask: why obtaining a Vector{Int64} instead of Vector{Float64} creates a bug for you?

Maybe you could make your code more generic (e.g. if you had a function expecting a Vector{Float64} input, the type specification could be relaxed to Vector{<:Real}, or maybe even to AbstractVector{<:Real}); and then get rid of such bugs?

Smart!

A slight amendment to the above suggestion:

func(r) = iszero(r) ? zero(r)/1 : r/3
func(r) = iszero(r) ? zero(r/1) : r/3 # alternative, but may behave slightly worse in some settings

The issue with the return values of zero(r) and r/3 is that zero(::Int) isa Int but r::Int/3 isa Float64. For example, when handling an array of Int it would have returned a Real vector again. zero(r)/1 ensures they get the same treatment and return the same type.

Although (unless this is an over-reduced example) one could have just used func(r) = r/3 since the zero case appears to be handled adequately by this.

4 Likes

Makes sense!

You are right that the reason for getting Real wasn’t my question. Initially I was writing a long message including that question, but during the process of writing, I realized the reason for getting Real. So, I deleted that part of my message and posted the message that I did.

I’m adding floating point numbers into the array which is initialized with zeros. Vector{Int64} gives an “inexact” error:

function(a, b) 
  res = zero(a) # initialize with zeros.
  for i = eachindex(a), j = eachindex(b)
      res[i] += function_of(b[j]) # inexact error: Int(1.123456e-16)
  end
  return res
end

The runtime of course gave me the error, but it took me a lot of time to find out why I was getting Vector{Int64} from the zero() funciton.

1 Like

That is an over-reduced example. In my original program

func(r) = (r==0) ? 0 : someotherfunction(r) / r

I just wanted to avoid division by zero (because I knew the function tends to zero when r → +0 and I didn’t bother to calculate the asymptotic form of the complicated function because I didn’t need such an accurate behavior near r = 0).

There is another interesting question here — perhaps zero(v::Vector{Real}) should be defined as zero.(v) instead of playing with the element type? I’d bet this definition dates way back into Julia’s prehistory (edit: prehistory indeed! .j files and all: julia@efb3e23)… and we try to limit how much behavior we embed into element types these days. Sure, it would no longer “save” you from a performance issue, but those are easy enough to find if they’re a trouble.

6 Likes

I have absolutely no type specifications! No types in sight except for typed literals like “0”, “1.3”, etc.

That was my initial thought when I found out I was getting Vector{Real}.

As I said in my initial message, it appears as if zero(0.333333333) were Int(0). (In reality, the zero() function doesn’t look at the values contained in the array. It just uses the eltype of the array, I guess.)

In this case, consider one of

func(r) = someotherfunction(r) / r * !iszero(r)
func(r) = iszero(r) ? zero(r/r) : someotherfunction(r)/r
func(r) = iszero(r) ? someotherfunction(zero(r))/one(r) : someotherfunction(r)/r

The first one will call someotherfunction regardless (unless the compiler manages to drop it, which is possible but uncertain), but multiplying by a Bool is either a no-op (true) or hard-zero (false) that can override even NaN or Inf so the result will be right in the end.
The second assumes that someotherfunction(r) returns something of the same type as r, so that zero(r/r) is a suitable zero. If someotherfunction does not have this property, you’ll need to be a bit more careful. The benefit is that someotherfunction is not called when iszero(r).
The third is much like the first except that it gives constant-propagation a better chance to resolve the iszero case at compile time unless someotherfunction is non-deterministic. It’s probably the “best” version in terms of performance and robustness.

1 Like

It was a “bug” in your code, i.e. it not type stable (that’s not illegal so bug may be the wrong word, at least surprising), but with it giving 0.0, i.e. Float64, it seemingly wouldn’t have been. So is it a bug in Julia? I think 0 is the intended behavior:

[Though I’m not sure, why not then the largest Int, Int128, or smallest Int8, or Bool, maybe this wasn’t given a lot of thought, and Float64 would be perfectly ok, even if Rationals in your Real array. With just Rationals, you would/will still get them for zero.]

I can see when you use concrete types, Int, Int32, Float64 etc, that you would want to get out, 0, 0, 0.0, but for abstract types, I’m not really sure, and Real is one of those (most are prefixed with Abstract, with few exceptions such as Real, because that’s what that general math concept is named).

All integers are Reals, but in general real numbers aren’t integers, they are flaot-like, but even more generally Rationals. What would be best in your view, 0.0, or even 0//0?

If we change your function to be type unstable in a different way, combining integers and a rational:

func(r) = (r == 0) ? 0 : r//3

julia> a = func.(0.0:1.0:5.0)
ERROR: MethodError: no method matching //(::Float64, ::Int64)

you actually get an ERROR (sooner), making a problem obvious, but in some other context you would have gotten away with and then still gotten a Int(0). Would you then have expected 0//0?

It would not be expected, too much work, to get a 0//0 in your case, if e.g. only one rational in a huge array of Reals. So I think you would still be surprised to see it produced if using only integers and floats (or generating, as you did).

What I would like, is that the REPL would warn when you make an type-unstable function, and you could ignore when you don’t care (it’s sometimes powerful to allow, not always a problem, though always slower).

That’s very good practice to avoid non-general code, but then you must also avoid 0. You get away with it though in the check. But 0.0 would also have worked there, and elsewhere then working for Float64:

Note, the is type-unstable for e.g. Float32 passed in, because then the division doesn’t get you Float64, but 0.0 will, that could be fixed with r/3.0. I don’t believe type-stable means getting the same type out as you put in, because you can pass many in an argument list, but some might consider this not, or at least not as generic, wanting Float32 out.

Well all constants are hidden type specification (and should be avoided if possible, also 1 or 2 etc. in indexing arrays, though not for type-stability reasons), and should be avoided. E.g. E=mc^2 was a problem in a thread here some time back. It lead to integer overflow… when and only when c defined as an integer (as in SI standard since 2019). E=mc^2.0 looks a bit odd but avoided that problem…

This /1 is rather obscure Julia knowledge I think, to get Float64 when you get integers in. Is there any other way? It looks odd and like extra (slow) division (but it will be optimized away(?)).

Indeed this is quite esoteric. My point was that if the operations “mirror” each other then you’ll probably be okay.

But since you ask, this should be compiled to nothing. Both because x/1 is a no-op for integers and IEEEFloats and because zero(r) is a constant so zero(r)/c can be resolved at compile time for any constant c (for most Number types, at least).

julia> code_llvm(x -> x/1, (Float64,)) # no-op
;  @ REPL[122]:1 within `#5`
define double @"julia_#5_2200"(double %0) #0 {
top:
  ret double %0
}

julia> code_llvm(x -> x/1, (Int64,)) # must convert to Float64, but still no division
;  @ REPL[125]:1 within `#7`
define double @"julia_#7_2202"(i64 signext %0) #0 {
top:
; ┌ @ int.jl:97 within `/`
; │┌ @ float.jl:294 within `float`
; ││┌ @ float.jl:268 within `AbstractFloat`
; │││┌ @ float.jl:159 within `Float64`
      %1 = sitofp i64 %0 to double
; └└└└
  ret double %1
}

More generally, x*1.0 and x-0.0 (but not x+0.0 except with certain fastmath flags because +(-0.0,+0.0) === +0.0) are no-ops for IEEEFloat numbers. And when inv(c) is exactly representable (such as when c is a power of 2 within range) then x/c can be replaced with x*inv(c).

I’d agree, because if we’re going to play with element type we should go to the narrowest which would mean Booleans.

The widest is not well defined, so we gotta go to the narrowest. Int for Real feels very arbitrary

2 Likes

I wasn’t really asking, stating/assyming, I thought I recalled, but since you confirmed I decided to check with julia -O0 and then predictably I got an fdiv, not a no-op. It’s rare to use this non-default, I see -O1 in effect often used (which DOES get a no-op), e.g. with Plots.jl, and I’ve suggested such to reduce startup cost myself. It’s probably outdated to lower optimization, since 1.9 when precompiling to native code became the default. Not all code is precompiled, not scipts, and I’ve suggested --compile=min for the interpreter in some cases, and it’s slow (but startup is quick) on top of division being slowest operation in machine code (though not much of a worry relite to interpreter slowness.

The actual question still stands, I’m curious is there another way, also for readability…

x*1.0 and x-0.0 are neither good since unlike x/1, they are not type-stable, they convert to Float64 from Int (or even from e.g. Float32).

If you’re commenting on my post, then that’s true at least sometimes. “Int for Real feels very arbitrary” right, and I’m thinking would a non-integer have any ill effect? Float64 simplest, would even have “solved” the original problem. Getting “Real” was non-ideal, because of type-instability, and with Float64 would have hidden that problem, or should I say have gotten away with it, changing it only to a performance-problem. I trying to think of were it would be a problem, for Rationals? Or vise versa getting a Rational a problem for floats?

[I’ll correct these words of mine in a separate message. Sorry for the confusion]
No! But I don’t know what you are getting at. I said it’s “my bug”.

It took me sometime to locate the source of, or the reason for, my bug. So, I wondered what’s the best way to avoid problems of this kind in the future. I also wondered why zero(Real) must return Int(0).

The rest of your discussion (thanks!) indicates that this isn’t a simple problem.

1 Like

You missed my second sentence:

No types [are] in sight except for typed literals like “0”, “1.3”, etc.

So, how would one avoid the concretely-typed 0 ? I wonder if a “generic zero” (and “generic one”) would make sense. . . .

π is generic:

julia> typeof(π)
Irrational{:π}

Would it be possible to have a literal zero which is generic, perhaps like

julia> typeof(00)
Number{:zero}

Also, would it be possible to have a mechanism to tell the compiler that “cond ? trueexpression : falsexpression” should return a value of the same concrete type, when possible to do so. (That’s what you said about a potential REPL warning.)