Adding the possibility to specify the number of decimal digits instead of using the number of bits

[SergeStinckwich]

From @SergeStinckwich on June 28, 2015 20:0

We should be able to specify the precision in number of digital digits instead of just number of bits :

1 asArbitraryPrecisionFloatNumDecimalDigits: 37

Copied from original issue: SergeStinckwich/SciSmalltalk#22

[SergeStinckwich]

From @nicolas-cellier-aka-nice on June 28, 2015 20:4

Hi Serge, Yes, it would ne nice, but we have to find the meaning and then the corresponding formula... Is it simply the smallest p such that 2^p > 10^n?

2015-06-28 22:00 GMT+02:00 Serge Stinckwich [email protected]:

We should be able to specify the precision in number of digital digits instead of just number of bits :

1 asArbitraryPrecisionFloatNumDecimalDigits: 37

— Reply to this email directly or view it on GitHub https://github.com/SergeStinckwich/SciSmalltalk/issues/22.

[SergeStinckwich]

At the moment this is implemented as:

Number>>asArbitraryPrecisionFloatNumDecimalDigits: n
    ^ self asArbitraryPrecisionFloatNumBits: (n / (2 log: 10)) ceiling

Do you agree with this ?

Thank you Nicolas.

[SergeStinckwich]

From @nicolas-cellier-aka-nice on June 29, 2015 23:12

2015-06-28 22:09 GMT+02:00 Serge Stinckwich [email protected]:

At the moment this is implemented as:

Number>>asArbitraryPrecisionFloatNumDecimalDigits: n ^ self asArbitraryPrecisionFloatNumBits: (n / (2 log: 10)) ceiling

Do you agree with this ?

Thank you Nicolas.

— Reply to this email directly or view it on GitHub https://github.com/SergeStinckwich/SciSmalltalk/issues/22#issuecomment-116322386 .

That sounds fine. The only question I have is about the accuracy of the formulation. I just tried this snippet:

{ (1 to: 10000) select: [:n | (10 raisedTo: n) highBit ~= (n / (2 log: 10)) ceiling]. (1 to: 10000) select: [:n | (10 raisedTo: n) highBit ~= (n * (10 log: 2)) ceiling]. }

Both are empty on Mac which means that above floating point approximations are accurate enough to support at least 10 000 figures.

To be a bit faster: (1 to: 1000000) inject: 10 into: [:pow :n | pow highBit = (n / (2 log: 10)) ceiling ifFalse: [^n halt]. pow*10]. nil is nil on my Mac, so formulation is valid at least up to a million figures.

It also seems that on my Mac, both expressions are true: { (2 log: 10) < ((2 asArbitraryPrecisionFloatNumBits: 200) log: 10). (10 log: 2) > ((2 asArbitraryPrecisionFloatNumBits: 200) log: 10). } Which might be a good property w.r.t. ceiling evaluation... It's just that I don't know whether or not this property will hold with various libm log functions used on various VM. Nor what is the limit?

Eventually, we could have used an ArbitraryPrecisionFloat for evaluating log:. Since we know that (10 log: 2) is irrational, we are sure that (n * (10 log: 2)) is not a whole number. So, if float approximation is near enough to a whole number (say 2 ulp to account for faithfully rounded log) then we just have to extend precision and retry.

So something like this would be safer, but maybe overkill: | precision nbits | precision := n highBit * 2 max: 24. [nbits := n * ((10 asArbitraryPrecisionFloatNumBits: precision) log: 2). (nbits ceiling - nbits) abs > (2 * nbits ulp)] whileTrue: [precision := precision + 10]. ^nbits ceiling

Last thing; Werner already told me that log: needs to be redefined in ArbitraryPrecisionFloat and use extended precision to be at least faithfully rounded...

Nicolas

[SergeStinckwich]

From @nicolas-cellier-aka-nice on June 29, 2015 23:15

2015-06-30 1:12 GMT+02:00 Nicolas Cellier < [email protected]>:

2015-06-28 22:09 GMT+02:00 Serge Stinckwich [email protected]:

At the moment this is implemented as:

Number>>asArbitraryPrecisionFloatNumDecimalDigits: n ^ self asArbitraryPrecisionFloatNumBits: (n / (2 log: 10)) ceiling

Do you agree with this ?

Thank you Nicolas.

— Reply to this email directly or view it on GitHub https://github.com/SergeStinckwich/SciSmalltalk/issues/22#issuecomment-116322386 .

That sounds fine. The only question I have is about the accuracy of the formulation. I just tried this snippet:

{ (1 to: 10000) select: [:n | (10 raisedTo: n) highBit ~= (n / (2 log: 10)) ceiling]. (1 to: 10000) select: [:n | (10 raisedTo: n) highBit ~= (n * (10 log: 2)) ceiling]. }

Both are empty on Mac which means that above floating point approximations are accurate enough to support at least 10 000 figures.

To be a bit faster: (1 to: 1000000) inject: 10 into: [:pow :n | pow highBit = (n / (2 log: 10)) ceiling ifFalse: [^n halt]. pow*10]. nil is nil on my Mac, so formulation is valid at least up to a million figures.

It also seems that on my Mac, both expressions are true: { (2 log: 10) < ((2 asArbitraryPrecisionFloatNumBits: 200) log: 10). (10 log: 2) > ((2 asArbitraryPrecisionFloatNumBits: 200) log: 10). } Which might be a good property w.r.t. ceiling evaluation... It's just that I don't know whether or not this property will hold with various libm log functions used on various VM. Nor what is the limit?

Eventually, we could have used an ArbitraryPrecisionFloat for evaluating log:. Since we know that (10 log: 2) is irrational, we are sure that (n * (10 log: 2)) is not a whole number. So, if float approximation is near enough to a whole number (say 2 ulp to account for faithfully rounded log) then we just have to extend precision and retry.

So something like this would be safer, but maybe overkill: | precision nbits | precision := n highBit * 2 max: 24. [nbits := n * ((10 asArbitraryPrecisionFloatNumBits: precision) log: 2). (nbits ceiling - nbits) abs > (2 * nbits ulp)] whileTrue: [precision := precision + 10]. ^nbits ceiling

Ahem... whileFalse: !!! coding that late in a dead environment like gmail should be forbidden!

Last thing; Werner already told me that log: needs to be redefined in ArbitraryPrecisionFloat and use extended precision to be at least faithfully rounded...

Nicolas