Math (exponent) Bug

38 views
Skip to first unread message

Sphaerica

unread,
Dec 23, 2011, 4:57:33 PM12/23/11
to Skulpt
Okay, my very first sample program ran into a math issue. I'm hoping
this is an easy fix. I'll play with it, to see if I can fix it
myslef, but given the following program (just does a scientific
computation using the Stefan-Boltzmann Law):

T=255
k=5.6704*(10**-8)
j=k*(T**4)
C=T - 273
F=32+C*9/5
print T,"Kelvin"
print C,"Celsius"
print F,"Fahrenheit"
print j,"W/m<sup>2</sup>"

It gives the following output:

255 Kelvin
-18 Celsius
-0.3999999999999986 Fahrenheit
0 W/m2

Note that that last value (zero) is wrong.

If I change the calculation of j in the program to:

j=k*T*T*T*T

or

import math
j=k*math.pow(T,4)

I get the right answer:

239.75872343999998 W/m2

Note that the following also gives zero, so it appears to be an issue
related to a very small value and a very large value:

j=k*(T*T*T*T)

Kees Bos

unread,
Dec 24, 2011, 3:50:20 AM12/24/11
to sku...@googlegroups.com

Strange. I tried this in javascript:
var T=255
var k=5.6704E-8
var j=k*(T*T*T*T)

and j is 239.75872344

So, its not a underlying javascript limitation.


Sphaerica

unread,
Dec 24, 2011, 9:37:23 AM12/24/11
to Skulpt
Kees,

Yes, I tried a lot of the variations without issue in Javascript, so
clearly it has something to do with how the expression is being
assembled, or how abstract.js is executing the assembled instructions.

Today I dumped some of the assembled javascript code. Some time this
week I'll put some console.log() statements into abstract.js to try to
see where it's going wrong (unless I can figure out a way to run the
raw code instead of the minimized code, in which case then I can use
the Firebug debugger).

Sphaerica

unread,
Dec 25, 2011, 1:16:36 PM12/25/11
to Skulpt
Bug is in lng.prototype.nb$power

Sphaerica

unread,
Dec 25, 2011, 11:39:56 PM12/25/11
to Skulpt
Okay, I could be wrong. Looking at long.js, it's funky, so the "bad"
result of 31745,30731,3,0 may not be bad at all...

It may be further down the line, where the long is multiplied by that
result, or maybe when 5.6704e-8 is converted to a long to do the
multiplication.

Still looking...

Sphaerica

unread,
Dec 25, 2011, 11:58:01 PM12/25/11
to Skulpt
It is in the conversion of 5.6704e-8 to a long. Any number less than
1 converts to 0, while numbers 1+ convert properly.

I'll look at it more closely tomorrow.

Sphaerica

unread,
Dec 26, 2011, 2:28:17 PM12/26/11
to Skulpt
Sigh. Just an update in case Scott happens to look at this today.

The difficulty is that the conversion to long function (lng.fromInt$)
uses bitwise operations to parse the number into the new bigint
format, but Javascript represents numbers internally in 64 bit
floating point form, while the bitwise operations only work on a 32-
bit integer version of a number. To do this Javascript basically
truncates the fractional portion of any number, returns the integer
portion as a 32 bit number (losing precision), and then this is what
lng.fromInt$ converts.

To solve the problem I will probably do something kludgey like convert
the number to binary using toString(2), detect the decimal point, and
if one is found, work on the string rather than doing actual bit ops.
This will be as slow as whatever.

Alternately, at least for my purposes, lng could just be altered to
use the internal 64-bit representation. I'm not sure what this will
mean in terms of compatibility with other Python implementations and
expected expression results, but for my purposes it's better (gives a
faster and guaranteed to be bug-free result). I'd just modify lng to
do nothing more than return the simplest Javascript version of an
operation.

I'd like some input back from Scott as to why lng was added instead of
just letting Javascript do numbers its way, so I can decide how to
approach this. And if I've finished by the time he sees this, so be
it.

Sphaerica

unread,
Dec 26, 2011, 2:54:23 PM12/26/11
to Skulpt
Looking further, it seems lng does not implement division. The
following program throws

TypeError: unsupported operand type(s) for Div: 'long' and 'long'

T=255
k=1.5777
j=k*(T**4)
print T
print k
print j
k=10
for i in range(1,10):
j=j/k
print j


Maybe I should instead be asking why these simple values (255**4,
5.6704e-8) are being converted to "long" in the first place? Perhaps
the problem, that lng is being invoked when it's not necessary.


Sk.builtin.lng.divremFull$ = function(v1, w1)
{
throw "todo;";
/*
var size_v = Math.abs(v1.size$);
var size_w = Math.abs(w1.size$);
var d = Sk.builtin.lng.BASE$ / (w1.digit[size_w - 1] + 1);
var v = Sk.builtin.lng.mulInt$(v1, d);
var w = Sk.builtin.lng.mulInt$(w1, d);
*/
};

Sphaerica

unread,
Dec 26, 2011, 3:42:48 PM12/26/11
to Skulpt
Okay, so the simplest solution for now was to prevent the conversion
to long by setting:

Sk.builtin.lng.threshold$ = Math.pow(2, 46);

Values above 46 cause py test failures and hang the m script, I'm
unsure why. I'd like to understand that, but I'm out of energy on
this.

Really, however, I should just alter the tests in abstract.js and
ast.js to skip this, because it's more a question of precision than
size, and lng doesn't improve the precision as far as I can tell, so
it doesn't help.

Using long is of little value if division isn't implemented. I'll
have to look into that. For now, however, this change allows my
original script to run. I may alter long to throw an error if it
attempts to convert a number with decimals, or else find a workaround
for the bit problem, or just force lng to work with Javascript numbers
and forget the bigint implementation. I'll also look into things
like:

http://leemon.com/crypto/BigInt.html
http://www.onicos.com/staff/iz/amuse/javascript/expert/BigInt.txt

or anything else I can find that might be easily swapped in to the
existing lng functions.

Brad Miller

unread,
Dec 26, 2011, 6:02:01 PM12/26/11
to sku...@googlegroups.com
Hi,

This is on my todo list.  I've run across this bug or a very similar one before myself.   Divide and other operators are not implemented for a infinite precision integers. 

Brad

-- 
Brad Miller
Sent with Sparrow

Sphaerica

unread,
Dec 26, 2011, 6:37:52 PM12/26/11
to Skulpt
Am I right in thinking that the switch to lng shouldn't have happened,
though? Isn't Python just IEEE 754 64-bit double precision, like
javascript? Why include lng at all? Or, if it's included, shouldn't
it only be activated when the necessary precision is exceeded?

I think for what I'm doing I'm just going to go into abstract.js and
act.as and turn off those tests that cut over to using lng. I'm
trying to just move beyond this and get back to using Skulpt instead
of digging under the hood.

Thoughts?

-- Bob

On Dec 26, 6:02 pm, Brad Miller <bonel...@gmail.com> wrote:
> Hi,
>
> This is on my todo list.  I've run across this bug or a very similar one before myself.   Divide and other operators are not implemented for a infinite precision integers.
>
> Brad
>
> --
> Brad Miller
> Sent with Sparrow (http://www.sparrowmailapp.com)

Sphaerica

unread,
Dec 27, 2011, 9:07:12 AM12/27/11
to Skulpt
Okay, I think I got the right answer, at least as far as the way the
code was intended.

ast.jsp has:

if ((k > Sk.builtin.lng.threshold$ || k < -Sk.builtin.lng.threshold$)
&& Math.floor(k) === k
&& (s.indexOf('e') === -1 && s.indexOf('E') === -1))

while abstract.js only does the threshold test. The Floor() call is
clearly meant to be certain that there is no decimal component, and so
should be included. The 'e' tests make sure there is no scientific
component, but it's being applied to the conversion of a string to a
number. I'm not certain how one would test for this, other than to
convert k to a string and test for either the 'e' or an unbroken
string of 0's at the end.

For my purposes, though, I think the the floor test is adequate (i.e.
don't lose a fractional portion of a value).

Also... I wonder if the threshold test should actually test for k ==
Number.POSITIVE_INFINITY or k == Number.NEGATIVE_INFINITY, as this is
the best evidence of an overflow. For now, though, I'll leave
threshold at the new Math.pow(2,46)

-- Bob

Sphaerica

unread,
Dec 27, 2011, 9:58:00 AM12/27/11
to Skulpt
Okay, it's just a little trickier than that, because the threshold has
to be left at 30 unless lng is fixed to deal with the fact that shift
operations only work on 32, not 64 bits. Putting it higher starts an
infinite loop somewhere in lng.

But the other changes also don't help me, because I get 255*4, which
needs to convert to a lng, but then the math fails when that is
multiplied by 5.6704e-8 (i.e. a decimal component, yielding a value
that no longer needs to be a long).

So ultimately I need to fix lng, and lng multiplication should be
expanded to yield a "float" when necessary.

For my purposes, the easier solution is to expand the lng conversion
to work past 32 bits, and then to set threshold$ up to 46. I can
probably do the 32 bit thing just by dividing the incoming number by
2**32 and doing the shift-thing twice.

Sphaerica

unread,
Dec 27, 2011, 5:44:09 PM12/27/11
to Skulpt
Heading out to dinner, but this is done (at least so far as fixing
lng.fromInt$, so that threshold$ can safely be raised to 2**46.

I'll post the code later tonight.

-- Bob

Sphaerica

unread,
Dec 27, 2011, 11:20:28 PM12/27/11
to Skulpt
Thinking about it, I don't like the way I wrote it. I have a shorter,
cleaner way I'll implement in the AM, then I'll post that.

Sphaerica

unread,
Dec 28, 2011, 10:34:30 AM12/28/11
to Skulpt
Here's the better code. Note that I didn't count digits first and
then apply them after, because the logic uses division instead of
shifting (although an intelligent implementation should do so behind
the scenes since the divisor is a power of 2, and I assume this is
what happens).

Also, the safety net isn't really needed, but I put it in just in
case... as my earlier thread pointed out, and infinite loop within
eval() doesn't get caught by "script running too long" logic and so
requires a force quit on the browser, so I'm avoiding that at all
costs these days.

Sk.builtin.lng.SHIFT$ = 15;
Sk.builtin.lng.BASE$ = 1 << Sk.builtin.lng.SHIFT$;
Sk.builtin.lng.MASK$ = Sk.builtin.lng.BASE$ - 1;
Sk.builtin.lng.threshold$ = Math.pow(2, 46); // Was Math.pow(2,30)
Sk.builtin.lng.TWOTOTHESHIFT$ = Math.pow(2,Sk.builtin.lng.SHIFT$);

Sk.builtin.lng.fromInt$ = function(ival)
{
var negative = false;
if (ival < 0)
{
ival = -ival;
negative = true;
}

var t = ival
var d = new Array(0)
while (t) {
var b = t & Sk.builtin.lng.MASK$
d.push(b)
t = (t - b) / Sk.builtin.lng.TWOTOTHESHIFT$
if (d.length > 500) // Safety net, to avoid an infinite loop
throw "Internal Error: Too many digits encountered converting to
long."
}

var ret = new Sk.builtin.lng(undefined, d.length)
if (negative) ret.size$ = -ret.size$

for (var i in d) {
ret.digit$[i] = d[i]
}

return ret
};
Reply all
Reply to author
Forward
0 new messages