[vtk-developers] vtkMath::Ceil() and vtkMath::Floor()

[Philippe Pébay]

Hello David,

Thanks for the clarification. Especially this part:
 > On every floating-point system that I've ever heard of, integers can
 > be represented exactly and uniquely by floats.  Even on the rare
 > non-standard ones that don't require normalization of the mantissa.
which I did not know and was the main reason for my wondering.

So, all in all, I think I am most likely to use your implementation of Floor() and Ceil() in most cases, except when I 
have strong reasons to believe that those doubles will go out of int range in typical application cases.

Thanks again for the discussion,
Philippe



On 06/03/2011 11:11 AM, David Gobbi wrote:
> On Fri, Jun 3, 2011 at 11:49 AM, Philippe Pébay<pppebay at sandia.gov>  wrote:
>> Hello David,
>>
>> Yes, I am sorry, I was very unclear, and ended up confusing the entire
>> question.
>>
>> I entirely agree with the fact that, in general, the composition of the 2
>> casts is not the identity and that it should not, because this is exactly
>> what a truncation does: to project an entire segment of reals onto a single
>> integer.
>>
>> So, if I had been clearer yesterday I would have added: if we have this x, a
>> double, which is representable as an int, then we cast it into
>> static_cast<int>( x ), which is an int indeed. So far so good.
>>
>> Then, what happens when we cast back this int into a double, thus obtaining
>> static_cast<double>( static_cast<int>( x ) )? Are we absolutely certain that
>> we obtain again the same FP representation as we had for x, initially?
>
> Yes, it will be exactly the same FP representation unless the original
> x was "-0", in which case the result will be "0".  But these two
> representations (meaning "0" and "-0") compare as equal.
>
>> Or, could it be that we obtain a certain double, y, close to x, but yet
>> different in its FP representation?
>
> Impossible.  An integer-valued double has a unique representation,
> except for the 0, -0 representations that I mentioned above. Maybe you
> are thinking about the case of decimal fractions, where the value
> "1.1" will have multiple floating-point representations.  But here, we
> are talking about integers, not fractions.  The representation of ints
> by doubles is exact and unique.
>
>> My question in fact revolves the question as to what will occur at run-time,
>> on different platforms, when moving expanding the int representation into a
>> double representation.
>
> On every floating-point system that I've ever heard of, integers can
> be represented exactly and uniquely by floats.  Even on the rare
> non-standard ones that don't require normalization of the mantissa.
>
>   - David
>
>
>
>> So, this is what I meant by the non-identity operation: in the case where we
>> would rightfully expect one, will we always have one?
>>
>> Again, sorry for having missed to be explicit about that, which certainly
>> did not help!
>>
>> Thanks for your time
>> Philippe
>>
>>
>> On 06/02/2011 05:01 PM, David Gobbi wrote:
>>>
>>> HI Philippe,
>>>
>>> Of course static_cast<double>(static_cast<int>(x)) is not identity.
>>> It is a truncation operation. There are over a million different double
>>> values that have the same truncated value.  I'm not sure what your
>>> point is regarding that.
>>>
>>> In this inequality:
>>>
>>> x != static_cast<double>(static_cast<int>(x))
>>>
>>> the purpose is to compare the original x (with all its trailing bits)
>>> against the truncation of x (with all the trailing bits removed).
>>>
>>> If "x" is exactly equal to trunc(x), and "y" is exactly equal to trunc(x),
>>> then "x" is exactly equal to "y".  If they are even one bit different,
>>> then they cannot both be exactly equal to the same integer.
>>> (excluding the case of +0 and -0 which IEEE 754 guarantees will
>>> be equal when compared to each other, or the case where "x"
>>> or "y" are too large to be represented by int).
>>>
>>> I can add comments to the code, but only after this is resolved ;)
>>>
>>>   - David
>>>
>>>
>>>
>>> On Thu, Jun 2, 2011 at 5:09 PM, Philippe Pébay<pppebay at sandia.gov>    wrote:
>>>>
>>>> David
>>>>
>>>> Yes, the former. So, what I meant was that we may have, the following
>>>> situation:
>>>>
>>>> * x and y are doubles, different in FP representation (but very close,
>>>> say
>>>> just above machine epsilon)
>>>> * but yet static_cast<int>( x ) is the same as static_cast<int>( y )
>>>>
>>>> and thus the composition of static_cast<double>    with static_cast<int>
>>>>   is not
>>>> the identity. Which may be entirely fine for most cases. But this becomes
>>>> a
>>>> feature which is worth being documented, since those trailing bits are
>>>> thrown away on purpose as you said.
>>>>
>>>> Now it might very well be that this exactly maps what ceil() and floor()
>>>> do
>>>> on all platforms, but I am not sure.
>>>>
>>>> Thanks
>>>> Philippe
>>>>
>>>> On 06/02/2011 03:50 PM, David Gobbi wrote:
>>>>>
>>>>> Just a quick follow up: when you say trailing bits, do you mean the
>>>>> trailing bits after those that represent an integer value, or do you
>>>>> mean the bits that might exist beyond the precision of double (i.e.
>>>>> extended precision).  I was assuming the former.
>>>>>
>>>>>   - David
>>>>>
>>>>>
>>>>> On Thu, Jun 2, 2011 at 4:29 PM, David Gobbi<david.gobbi at gmail.com>
>>>>>   wrote:
>>>>>>
>>>>>> Hi Philippe,
>>>>>>
>>>>>> The trailing bits are important in this line:
>>>>>> int g = ( x<      0 );
>>>>>> which compares the original "x" against zero.
>>>>>>
>>>>>> But the static_cast<int>(x) throws away all trailing bits on purpose,
>>>>>> the intent of the double cast is to check whether throwing away the
>>>>>> trailing bits changes the value of the double.
>>>>>>
>>>>>>   - David
>>>>>>
>>>>>>
>>>>>> On Thu, Jun 2, 2011 at 4:22 PM, Philippe Pébay<pppebay at sandia.gov>
>>>>>>   wrote:
>>>>>>>
>>>>>>> Hello David,
>>>>>>>
>>>>>>> Thanks for the explanation. Also, sorry for the typo, I meant "double
>>>>>>> cast"
>>>>>>> (not "double truncation").
>>>>>>>
>>>>>>> My question was with the != comparison between x and
>>>>>>> static_cast<double>(
>>>>>>> static_cast<int>( x ) ).
>>>>>>>
>>>>>>> So, as I understand your answer, the intent is effectively to give a
>>>>>>> definition of what it means
>>>>>>> to be an integer, and in this case, it means that static_cast<double>(
>>>>>>> static_cast<int>( x ) ) is the same as x.
>>>>>>>
>>>>>>> What I don't understand here is what happens to the trailing bits,
>>>>>>> which
>>>>>>> are
>>>>>>> still significant with respect to machine epsilon. What will happen
>>>>>>> when
>>>>>>> casting x to int?
>>>>>>>
>>>>>>> Thanks
>>>>>>> Philippe
>>>>>>>
>>>>>>>
>>>>>>> On 06/02/2011 03:01 PM, David Gobbi wrote:
>>>>>>>>
>>>>>>>> Hi Philippe,
>>>>>>>>
>>>>>>>> There is only one truncation, not two.  For Floor() the explanation
>>>>>>>> is:
>>>>>>>>
>>>>>>>> x is input (a double, or maybe extended precision in FPU register)
>>>>>>>>
>>>>>>>> // step 1: truncate to integer
>>>>>>>> int r = static_cast<int>(x);
>>>>>>>>
>>>>>>>> // step 2: check if x can be represented exactly as an int,
>>>>>>>> // put the boolean result "1" or "0" in variable "n"
>>>>>>>> int n = ( x != static_cast<double>(r) );
>>>>>>>>
>>>>>>>> We know that the integer "r" can be represented exactly by
>>>>>>>> a double, and so the above checks if "x" can, in turn, be
>>>>>>>> represented exactly by an int.
>>>>>>>>
>>>>>>>> // step 3: check if x is above/below zero, put the boolean
>>>>>>>> // result "1" or "0" in variable "g"
>>>>>>>> int g = ( x<        0 );
>>>>>>>>
>>>>>>>> // step 4: subtract 1 from the result if steps 2 and 3 are true
>>>>>>>> return r - ( n&        g );
>>>>>>>>
>>>>>>>>
>>>>>>>> The idea behind this is the following: when "x" is>= 0, doing
>>>>>>>> static_cast<int>(x) is exactly the same as floor(x).  But if "x"
>>>>>>>> is less than zero (step 3) we need to subtract one unless "x"
>>>>>>>> can be represented exactly by an int (step 2).
>>>>>>>>
>>>>>>>> The purpose of step 4 is to avoid any branching in the code,
>>>>>>>> as compared to doing the same logic with "if" statements.
>>>>>>>>
>>>>>>>>   - David
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> On Thu, Jun 2, 2011 at 3:25 PM, Philippe Pébay<pppebay at sandia.gov>
>>>>>>>>   wrote:
>>>>>>>>>
>>>>>>>>> Hello David,
>>>>>>>>>
>>>>>>>>> Thanks for the prompt answer.
>>>>>>>>>
>>>>>>>>> I still don't understand the double truncation for Floor() and
>>>>>>>>> Ceil():
>>>>>>>>> why
>>>>>>>>> cast the readings to int first, and then to double again? I am not
>>>>>>>>> sure I
>>>>>>>>> understand what "exact results" mean in this contest, as the cast
>>>>>>>>> itself
>>>>>>>>> changes the nature of the representation, and a double cast does
>>>>>>>>> not,
>>>>>>>>> as
>>>>>>>>> far
>>>>>>>>> as I know at least, amount to an identity operation.
>>>>>>>>>
>>>>>>>>> Or does it? Are you saying that static_cast<double>(
>>>>>>>>> static_cast<int>(
>>>>>>>>> x
>>>>>>>>> ) )
>>>>>>>>> has always the exact same FP representation than x?
>>>>>>>>>
>>>>>>>>> Thanks
>>>>>>>>> Philippe
>>>>>>>>>
>>>>>>>>> On 06/02/2011 02:20 PM, David Gobbi wrote:
>>>>>>>>>>
>>>>>>>>>> Hi Philippe,
>>>>>>>>>>
>>>>>>>>>> I wrote the Floor() and Ceil() functions a few months ago.  They
>>>>>>>>>> are
>>>>>>>>>> written specifically for speed.  For double values that are within
>>>>>>>>>> the
>>>>>>>>>> int range, they give exact results.
>>>>>>>>>>
>>>>>>>>>> The vtkMath::Round() is fast, too, but the operation that it
>>>>>>>>>> performs
>>>>>>>>>> does not conform to any of the IEEE 754 standard rounding modes.
>>>>>>>>>>
>>>>>>>>>>   - David
>>>>>>>>>>
>>>>>>>>>> On Thu, Jun 2, 2011 at 3:09 PM, Philippe Pébay<pppebay at sandia.gov>
>>>>>>>>>>   wrote:
>>>>>>>>>>>
>>>>>>>>>>> ... a clarification of what I asked earlier...
>>>>>>>>>>>
>>>>>>>>>>> By "different intent" I was referring to the fact that vtkMath's
>>>>>>>>>>> Floor()
>>>>>>>>>>> and
>>>>>>>>>>> Ceil() perform something which I find a little intriguing, namely,
>>>>>>>>>>> casting a
>>>>>>>>>>> double to an int, and then re-casting the latter to a double to
>>>>>>>>>>> perform
>>>>>>>>>>> a
>>>>>>>>>>> double-double equality comparison, which, I would venture to say,
>>>>>>>>>>> should
>>>>>>>>>>> be
>>>>>>>>>>> avoided whenever possible.
>>>>>>>>>>>
>>>>>>>>>>> On the other hand, Round() does not do this double casting but
>>>>>>>>>>> directly
>>>>>>>>>>> performs double-double inequality comparison, and a single cast to
>>>>>>>>>>> an
>>>>>>>>>>> int
>>>>>>>>>>> at
>>>>>>>>>>> the end. This seems much better to me.
>>>>>>>>>>>
>>>>>>>>>>> Thanks
>>>>>>>>>>> Philippe
>>>>>>>>>>>
>>>>>>>>>>> On 06/02/2011 02:04 PM, Philippe Pébay wrote:
>>>>>>>>>>>>
>>>>>>>>>>>> Hello,
>>>>>>>>>>>>
>>>>>>>>>>>> I was contemplating using vtkMath::Ceil() and vtkMath::Floor() in
>>>>>>>>>>>> some
>>>>>>>>>>>> of my code, instead of the standard ceil() and floor(). However
>>>>>>>>>>>> they
>>>>>>>>>>>> appear to be written with a slightly different intent than
>>>>>>>>>>>> vtkMath::Round(). Can someone comment on that?
>>>>>>>>>>>>
>>>>>>>>>>>> Thanks
>>>>>>>>>>>> Philippe
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> --
>>>>>>>>>>> Philippe Pébay
>>>>>>>>>>> Sandia National Laboratories
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> Powered by www.kitware.com
>>>>>>>>>>>
>>>>>>>>>>> Visit other Kitware open-source projects at
>>>>>>>>>>> http://www.kitware.com/opensource/opensource.html
>>>>>>>>>>>
>>>>>>>>>>> Follow this link to subscribe/unsubscribe:
>>>>>>>>>>> http://www.vtk.org/mailman/listinfo/vtk-developers
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> --
>>>>>>>>> Philippe Pébay
>>>>>>>>> Sandia National Laboratories
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> Philippe Pébay
>>>>>>> Sandia National Laboratories
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> Philippe Pébay
>>>> Sandia National Laboratories
>>>>
>>>>
>>>>
>>>
>>
>>
>> --
>> Philippe Pébay
>> Sandia National Laboratories
>>
>>
>>
>


-- 
Philippe Pébay
Sandia National Laboratories