Whole Number Division

[czardas]

I have hardly any time to write any code recently. Even so, one of the most frustrating problems I, and others like myself, encounter is corruption from floating point innacuracies and integer overflow. I haven't had much time to test this, but the given example works. The method to fix division with whole numbers does not appear to be as complicated as I first anticipated. Divisible integers only!

expandcollapse popup

Func _WholeNumberDivision($iDividend, $iDivisor) ; Input ranges -9223372036854775807 To 9223372036854775807
    If Not (IsInt($iDividend) And IsInt($iDivisor)) Then Return SetError(1, 0, $iDividend / $iDivisor) ; integers only
    If $iDivisor = 0 Then Return SetError(2, 0, $iDividend / $iDivisor) ; division by zero

    Local $aDiv = [$iDividend, $iDivisor], _
    $iSign = 1

    For $i = 0 To 1
        If $aDiv[$i] > 0x7FFFFFFFFFFFFFFF Or $aDiv[$i] < 0x8000000000000001 Then Return SetError(3, 0, $iDividend / $iDivisor) ; input range exceeded
        If VarGetType($aDiv[$i]) = "Double" Then $aDiv[$i] = Number($aDiv[$i], 2) ; convert to Int-64

        If $aDiv[$i] < 0 Then ; force positive integers
            $aDiv[$i] *= -1
            $iSign *= -1 ; to add back later
        EndIf
    Next

    If Mod($aDiv[0], $aDiv[1]) Then Return SetError(4, 0, $iDividend / $iDivisor) ; not divisible
    If $aDiv[0] = 0 Then Return 0
    If $aDiv[1] = 1 Then Return $aDiv[0] * $iSign

    Local $iDivision = Floor($aDiv[0] / $aDiv[1]), $iDifference, $iIntegral

    While $iDivision * $aDiv[1] > $aDiv[0] ; division is overstated
        $iDifference = ($aDiv[1] * $iDivision) - $aDiv[0]
        $iIntegral = Floor($iDifference / $aDiv[1]) ; avoid shooting beyond the target
        If $iIntegral = 0 Then $iIntegral = 1 ; prevents hanging in an infinite loop
        $iDivision -= $iIntegral
    WEnd

    While $iDivision * $aDiv[1] < $aDiv[0] ; division is understated
        $iDifference = $aDiv[0] - ($aDiv[1] * $iDivision)
        $iIntegral = Floor($iDifference / $aDiv[1])
        If $iIntegral = 0 Then $iIntegral = 1 ; prevents hanging
        $iDivision += $iIntegral
    WEnd

    Return $iDivision * $iSign
EndFunc

This function currently works with all int-64 values with one exception - the lowest value 0x8000000000000000, and that's only divisible by powers of 2 anyway.

Edited August 26, 2015 by czardas

[czardas]

Perhaps I ought to have waited a little while before posting the above function, but I saw fit to post a solution to a rather messy problem as proof of concept. The function may be limited to divisible integers right now, and that's mostly how I intend to use it - particularly with the next implementation of my Fraction UDF. The absence of error return values is both temporary and superfluous. There is still a lot I don't understand about floating point maths, however I don't feel so bad about the fact because I know I'm not alone in this. It seems to me that even 'experts' make apparently contradictory statements regarding some of the finer details of double precision floating point arithmetic. One thing you can say for certain is that it is pretty darn fast! 

It's probably going to take me a while to progress the way things are right now: I simply have too many commitments. For this reason I intend to post related snippets as and when they are created - but only if I think they might be useful to someone else. The following function grabs the first 16 (rounded) digits from a float. The 17th digit is unreliable. The sixteenth digit may become unreliable after a mathematical operation - at least I believe this is the case. The extended value must be accessed (immediately): to preserve information about the numeric order of magnitude.
 

Local $fValue = (1.025) ; Try inputting some different floats.
Local $sMachineValue = StringFormat('%.16e', $fValue)
Local $iDigits = _FloatToDigits($fValue)
Local $iExponent = @extended
MsgBox(0, 'Value => ' & $fValue, _
    'Machine Value' & @TAB & $sMachineValue & @LF & _
    'Extracted Digits' & @TAB & $iDigits & " * 10^" & $iExponent)

; Returns a 32-bit or 64-bit signed integer. Sets @extended to the decimal exponent (float = int * 10 ^ exponent).
Func _FloatToDigits($fFloat)
    If VarGetType($fFloat) <> 'Double' Or StringInStr($fFloat, '#') Then Return SetError(1)
    Local $iSign = ($fFloat < 0) ? -1 : 1, $iDigits = 15 ; machine epsilon = 5 × 10^-15
    $fFloat = StringFormat('%.' & $iDigits & 'e', $fFloat) ; rounds to 15 decimal places

    Local $aFloat = StringSplit($fFloat, "e", 2) ; zero-based array
    If $iSign < 0 Then $aFloat[0] = StringTrimLeft($aFloat[0], 1) ; Remove the minus sign.

    $aFloat[0] = StringLeft($aFloat[0], 1) & StringRegExpReplace(StringRight($aFloat[0], $iDigits), '(0+\z)', '') ; Remove the decimal point and trailing zeros.
    $aFloat[1] += 1 - StringLen($aFloat[0]) ; Adjust the exponent to accommodate changes.

    Return SetExtended($aFloat[1], Int($aFloat[0]) * $iSign) ; Add back the minus sign.
EndFunc

 

Edited June 16, 2015 by czardas

[czardas]

A bug has been discovered, affecting the function _WholeNumberDivision(), which caused division by 1 to return incorrect results. A patch has been added. After looking again at this, for a moment, I wondered why I had excluded the input value 0x8000000000000000. If anyone else wondered about this, here is the answer: you cannot change the sign of that number without invoking integer overflow. Error return values may also be added to this function at some point in the future.

Testing reveals further serious bugs. Pity I had no time to test this until now.

Edited August 23, 2015 by czardas

[czardas]

This was not working at all: apart from the few tests I did at the time I wrote this function. It seems to be working now, after several failed attempts to fix it. I've only run a few tests so far and I'm not sure it is entirely reliable, although I know why it was failing before making these changes.

expandcollapse popup

Global $g_Dividend = 0, $g_Divisor = 1

Local $aExamples = [ _
    [9223372036854775807, 7], _
    [9223372036854775552, 2], _
    [9999997800000121, 99999989], _
    [99999999999999990, 5], _
    [9223372036854775552, 4611686018427387776], _
    [99999999999999999, 1], _
    [6, 2] _
    ]

For $i = 0 To UBound($aExamples) -1
    ConsoleWrite( _ ; Internal AutoIt division
        $aExamples[$i][$g_Dividend] & " / " & $aExamples[$i][$g_Divisor] & " = " & _
        Number($aExamples[$i][$g_Dividend] / $aExamples[$i][$g_Divisor], 2) & _
    @LF & _ ; _WholeNumberDivision()
        $aExamples[$i][$g_Dividend] & " / " & $aExamples[$i][$g_Divisor] & " = " & _
        _WholeNumberDivision($aExamples[$i][$g_Dividend], $aExamples[$i][$g_Divisor]) & _
    @LF & @LF)
Next

#cs
9223372036854775807 / 7 = 1317624576693539329
9223372036854775807 / 7 = 1317624576693539401

9223372036854775552 / 2 = 4611686018427387905
9223372036854775552 / 2 = 4611686018427387776

9999997800000121 / 99999989 = 99999989
9999997800000121 / 99999989 = 99999989

99999999999999990 / 5 = 19999999999999997
99999999999999990 / 5 = 19999999999999998

9223372036854775552 / 4611686018427387776 = 2
9223372036854775552 / 4611686018427387776 = 2

99999999999999999 / 1 = 100000000000000001
99999999999999999 / 1 = 99999999999999999

6 / 2 = 3
6 / 2 = 3
#ce

Edited August 24, 2015 by czardas

[czardas]

I have been going about this the wrong way: almost randomly testing methods every time I found a new bug. All I needed was a couple of hours sleep, after which the solution turns out to be really very simple - I ought to kick myself. I decided to return the floats generated by division when errors occur. I contemplated forcing an integer error return value in such cases but couldn't see any advantage in doing so.

The range test below has not revealed any bugs and should not throw any errors. The start (and range) values can be modifed, and timers added if you were curious enough.

expandcollapse popup

Local $iStart1 = 1317624639329, _
      $iRange1 = 500, _
      $iStart2 = 101, _
      $iRange2 = 500

Local $iProduct, $iDiv, $iErrors = 0, $iBugs = 0

For $i = $iStart1 To $iStart1 + $iRange1
    For $j = $iStart2 To $iStart2 + $iRange2
        $iProduct = $i * $j

        ; Test 1
        $iDiv = _WholeNumberDivision($iProduct, $i)
        If @error Then
            ConsoleWrite($iProduct & " / " & $i & " ERROR " & @error & @LF)
            $iErrors += 1
        ElseIf $iDiv * $i <> $iProduct Then
            ConsoleWrite($iProduct & " / " & $i & " BUG FOUND" & @LF)
            $iBugs += 1
        EndIf

        ; Test 2
        $iDiv = _WholeNumberDivision($iProduct, $j)
        If @error Then
            ConsoleWrite($iProduct & " / " & $j & " ERROR " & @error & @LF)
            $iErrors += 1
        ElseIf $iDiv * $j <> $iProduct Then
            ConsoleWrite($iProduct & " / " & $j & " BUG FOUND" & @LF)
            $iBugs += 1
        EndIf
    Next
Next

MsgBox(0, "Range Analysis", "Encountered errors : " & $iErrors & @LF & _
                            "Bugs found : " & $iBugs)


Func _WholeNumberDivision($iDividend, $iDivisor) ; Input ranges -9223372036854775807 To 9223372036854775807
    If Not (IsInt($iDividend) And IsInt($iDivisor)) Then Return SetError(1, 0, $iDividend / $iDivisor) ; integers only
    If $iDivisor = 0 Then Return SetError(2, 0, $iDividend / $iDivisor) ; division by zero

    Local $aDiv = [$iDividend, $iDivisor], _
    $iSign = 1

    For $i = 0 To 1
        If $aDiv[$i] > 0x7FFFFFFFFFFFFFFF Or $aDiv[$i] < 0x8000000000000001 Then Return SetError(3, 0, $iDividend / $iDivisor) ; input range exceeded
        If VarGetType($aDiv[$i]) = "Double" Then $aDiv[$i] = Number($aDiv[$i], 2) ; convert to Int-64

        If $aDiv[$i] < 0 Then ; force positive integers
            $aDiv[$i] *= -1
            $iSign *= -1 ; to add back later
        EndIf
    Next

    If Mod($aDiv[0], $aDiv[1]) Then Return SetError(4, 0, $iDividend / $iDivisor) ; not divisible
    If $aDiv[0] = 0 Then Return 0
    If $aDiv[1] = 1 Then Return $aDiv[0] * $iSign

    Local $iDivision = Floor($aDiv[0] / $aDiv[1]), $iDifference, $iIntegral

    While $iDivision * $aDiv[1] > $aDiv[0] ; division is overstated
        $iDifference = ($aDiv[1] * $iDivision) - $aDiv[0]
        $iIntegral = Floor($iDifference / $aDiv[1]) ; avoid shooting beyond the target
        If $iIntegral = 0 Then $iIntegral = 1 ; prevents hanging in an infinite loop
        $iDivision -= $iIntegral
    WEnd

    While $iDivision * $aDiv[1] < $aDiv[0] ; division is understated
        $iDifference = $aDiv[0] - ($aDiv[1] * $iDivision)
        $iIntegral = Floor($iDifference / $aDiv[1])
        If $iIntegral = 0 Then $iIntegral = 1 ; prevents hanging
        $iDivision += $iIntegral
    WEnd

    Return $iDivision * $iSign
EndFunc

Now the real work can begin.

Edit
Now that this is actually working, further tests reveal the extent of the problems presented by floating point division in extreme fringe cases. While the inaccuracies do not make a great deal of difference in terms of comparative ratios, they make a very real and significant difference if you want to calculate the trajectory of a deep space probe or, in my case, simplify fractions correctly. Accumilative errors are totally unacceptable in these situations. As has often been said; AutoIt is not necessarily the best choice to use for every project, but the accuracy of a fraction comprising two signed 64-bit integers will leave double precision floats flailing around in the mud. I'm not refering to speed: but rather to extremely accurate calibration within a significant numeric range. There can clearly be no comparison, by any stretch of the imagination.
 

Local $iStart1 = 4500000000000000001, _
      $iRange1 = 1000, _
      $iStart2 = 2, _
      $iRange2 = 0

Local $iProduct, $iDiv, $iErrors = 0, $iBugs = 0

For $i = $iStart1 To $iStart1 + $iRange1 Step 10
    For $j = $iStart2 To $iStart2 + $iRange2
        $iProduct = $i * $j

        $iDiv = _WholeNumberDivision($iProduct, $j)
        If @error Then
            ConsoleWrite($iProduct & " / " & $j & " ERROR " & @error & @LF)
            $iErrors += 1
        ElseIf $iDiv * $j <> $iProduct Then
            ConsoleWrite($iProduct & " / " & $j & " BUG FOUND" & @LF)
            $iBugs += 1
        ElseIf $iDiv <> Int($iProduct / $j) Then
            ConsoleWrite($iProduct & " / " & $j & @LF & Int($iProduct / $j) & " wrong" & @LF & $iDiv & " correct" & @LF & @LF)
        EndIf
    Next
Next

MsgBox(0, "Range Analysis", "Encountered errors : " & $iErrors & @LF & _
                            "Bugs found : " & $iBugs)

Sample Results:

expandcollapse popup

9000000000000000502 / 2
4500000000000000001 wrong
4500000000000000251 correct

9000000000000001382 / 2
4500000000000000513 wrong
4500000000000000691 correct

9000000000000001402 / 2
4500000000000000513 wrong
4500000000000000701 correct

9000000000000001422 / 2
4500000000000000513 wrong
4500000000000000711 correct

9000000000000001442 / 2
4500000000000000513 wrong
4500000000000000721 correct

9000000000000001462 / 2
4500000000000000513 wrong
4500000000000000731 correct

9000000000000001482 / 2
4500000000000000513 wrong
4500000000000000741 correct

9000000000000001502 / 2
4500000000000000513 wrong
4500000000000000751 correct

9000000000000001522 / 2
4500000000000000513 wrong
4500000000000000761 correct

9000000000000001542 / 2
4500000000000001025 wrong
4500000000000000771 correct

9000000000000001562 / 2
4500000000000001025 wrong
4500000000000000781 correct

9000000000000001582 / 2
4500000000000001025 wrong
4500000000000000791 correct

It's good to be back.

Edited August 26, 2015 by czardas

[czardas]

Another related function - the next on my list - is detecting integer overflow occurring with addition, subtraction or multiplication. I wasn't quite sure how to approach this problem, but came up with the idea of using doubles to test for proximity: doubles are not accurate enough to make an exact comparison but that should not be necessary.

; test for overflow with the expression ==> 0x7FFFFFFFFFFFFFFF + 1
MsgBox(0, "", _OverflowDetect("+", 0x7FFFFFFFFFFFFFFF, 1))

Func _OverflowDetect($sOperator, $iOperand_1, $iOperand_2)
    If Not StringRegExp($sOperator, '\A\s*[\+\-\*]\s*\z') Then Return SetError(1) ; operator not recognized
    If Not StringInStr(VarGetType($iOperand_1), 'Int') Then Return SetError(2) ; meaningless request
    If Not StringInStr(VarGetType($iOperand_2), 'Int') Then Return SetError(3) ; ditto

    ; execute the expression
    Local $iExecute = Execute($iOperand_1 & $sOperator & $iOperand_2)

    ; execute the expression with the operands converted to doubles
    Local $fCompare = Execute('Number(' & $iOperand_1 & ', 3)' & $sOperator & 'Number(' & $iOperand_2 & ', 3)')

    ; the results should be approximately equal
    Return StringFormat('%.15e', $iExecute) <> StringFormat('%.15e', $fCompare)
EndFunc

This function is not as reliable as the new version posted in the topic: https://www.autoitscript.com/forum/topic/176620-operator64/
The following topic may also be of interest: https://www.autoitscript.com/forum/topic/176690-number-puzzle/

Edited September 1, 2015 by czardas