[Cuis-dev] Problems in class Number
Andres Valloud
ten at smallinteger.com
Tue Oct 8 20:34:12 PDT 2019
I played a bit with the guard clauses and found the attached one is
simpler yet just as fast.
On 10/8/19 20:11, Andres Valloud via Cuis-dev wrote:
> Regarding each+31, sure, 30*k+1 comes first, except when k = 0 because
> why would anyone try dividing by 1. So this is why that case is shifted
> by 30. However, when written this way, the actual divisor evaluation
> order is 31, 7, 11, and so on. It's more likely that a random integer
> is 0 mod 7 than 0 mod 31, and the sooner one detects exact division, the
> sooner the computation can stop. Because of that, I think the each+31
> case should be the last one in the division loop.
>
> On 10/8/19 19:17, Agustín Sansone via Cuis-dev wrote:
>> Hello!
>>
>> I agree with you. I don't think isPrime should send isProbablyPrime
>> because it could fail in the future.
>> I leave you here the implementation with this taken care of.
>> I wrote the (each+31) case first in the trial division loop, because
>> it is testing the 30*k+1 case, wich I also wrote first in the comment.
>>
>> Thanks,
>> Agustín
>>
>> El mar., 8 oct. 2019 a las 8:11, Juan Vuletich via Cuis-dev
>> (<cuis-dev at lists.cuis.st <mailto:cuis-dev at lists.cuis.st>>) escribió:
>>
>> Hi Folks,
>>
>> I agree with Andrés comments, and will just focusing on the proposed
>> changes.
>> (snip)
>>
>> On 10/8/2019 2:20 AM, Andres Valloud via Cuis-dev wrote:
>> > Agustin, nice to see someone looking into these kinds of things
>> :).
>> > ...
>> >> * The *raisedToInteger: exp modulo: m *method**in Integer has
>> a very
>> >> big problem. If we compute, for example, /"5 raisedTo: 0
>> modulo:
>> >> 0"/, this returns 1. This means, that according to
>> Smalltalk, the
>> >> rest of the division by 0 of 1(=5^0) is equal to 1 (Yes,
>> division by
>> >> zero!!). I think you can see the problem. This is due the
>> first line
>> >> of the method, that says /"(exp = 0) ifTrue: [^ 1].", /does
>> >> not check anything else. This problem can be easily fixed by
>> >> checking if m=0 just before.
>> >
>> > I agree, the current code appears to be wrong. The initials on
>> the
>> > code belong to Juan Vuletich and Nicolas Cellier. Guys, is there
>> > reason why e.g. 5 raisedTo: 0 modulo: 0 should answer 1 rather
>> than
>> > fail? I don't see any, but...
>> >
>> > Assuming the code is broken and needs to be fixed,
>> alternatively one
>> > could also write the initial guard clause like this:
>> >
>> > n = 0 ifTrue: [^1 \\ m].
>> >
>> > because the case m = 0 will fail.
>> > ...
>>
>> Just added this suggestion as an update to GitHub. Andrés, I did it
>> with
>> your author initials, it's your code!
>>
>> > ...
>> >> * The *isPrime *method in Integer makes some optimization in
>> order to
>> >> run the algorithm in O(sqrt(self)) instead of the naive
>> way in
>> >> O(self). This is very intelligent, but the constant factor
>> of this
>> >> method can be still improved significantly. I share with
>> you my
>> >> implementation of *isPrimeFast *with a small explanation.
>> This
>> >> implementation runs in general more than 3 times faster
>> than the
>> >> actual one. I leave you a test that checks the correctness
>> of it as
>> >> well, and some other tests that check this complexity I
>> mentioned.
>> >
>> > I see what you did there, but I do not know how to reproduce the
>> time
>> > tests you mention. I built a sample of integers between 1 and
>> 2^32 (I
>> > didn't go up to 2^64 because that would require O(2^32) operations
>> > each, and I want that to finish in reasonable time), and I get
>> > something like a 2x performance improvement rather than 3x. This
>> > seems to make sense because the approach you propose halves the \\
>> > operations (8 remain out of the 16 the current code is doing, for
>> > every batch of 30 potential divisors).
>> >
>> > slicer := 1024.
>> > thickness := 255.
>> > maxK := 1 bitShift: 32.
>> > integers := 1 to: maxK by: maxK // slicer
>> > :: inject: OrderedCollection new
>> > into: [:t :x |
>> > t add: x.
>> > thickness timesRepeat: [t add: t last + 1].
>> > t yourself]
>> > :: asArray.
>> > Time millisecondsToRun:
>> > [1 to: integers size do:
>> > [:x | (integers at: x) isPrime]]
>> >
>> > Using the above code (which I could not format more nicely in this
>> > email), I get about 4.8s for isPrime, and about 2.4s for
>> isPrimeFast.
>> >
>> > Generally, isPrime shouldn't send isProbablyPrime because
>> isPrime is
>> > meant to be deterministic, and one shouldn't assume that the
>> > probabilistic algorithm today will happen to provide the correct
>> > deterministic answer tomorrow.
>> >
>> > Why is the (each+31) case first in the trial division loop?
>> >
>> > Andres.
>>
>> I'll wait for your consensus on what to do here. Making isPrime not
>> send
>> isProbablyPrime sounds reasonable to me, but folks, I prefer to wait
>> for
>> your suggestion.
>>
>> Thanks,
>>
>> -- Juan Vuletich
>> www.cuis-smalltalk.org <http://www.cuis-smalltalk.org>
>> https://github.com/Cuis-Smalltalk/Cuis-Smalltalk-Dev
>> https://github.com/jvuletich
>> https://www.linkedin.com/in/juan-vuletich-75611b3
>> @JuanVuletich
>>
>> -- Cuis-dev mailing list
>> Cuis-dev at lists.cuis.st <mailto:Cuis-dev at lists.cuis.st>
>> https://lists.cuis.st/mailman/listinfo/cuis-dev
>>
>>
-------------- next part --------------
'From Cuis 5.0 [latest update: #3866] on 8 October 2019 at 8:30:48 pm'!
!Integer methodsFor: 'testing' stamp: 'sqr 10/8/2019 20:30:43'!
isPrimeFast2c
self < 3 ifTrue: [^self = 2].
self even ifTrue: [^false].
self \\ 3 = 0 ifTrue: [^false].
self \\ 5 = 0 ifTrue: [^false].
"Now 2, 3 and 5 do not divide self. So, self is of the form
30*k + {1, 7, 11, 13, 17, 19, 23, 29} for integer k >= 0.
The 31 case below is the 30k+1 case, excluding k = 0"
0 to: self sqrtFloor by: 30 do: [:each |
self \\ (each+7) = 0 ifTrue: [ ^false ].
self \\ (each+11) = 0 ifTrue: [ ^false ].
self \\ (each+13) = 0 ifTrue: [ ^false ].
self \\ (each+17) = 0 ifTrue: [ ^false ].
self \\ (each+19) = 0 ifTrue: [ ^false ].
self \\ (each+23) = 0 ifTrue: [ ^false ].
self \\ (each+29) = 0 ifTrue: [ ^false ].
self \\ (each+31) = 0 ifTrue: [ ^false ].
].
^true! !
More information about the Cuis-dev
mailing list