[Cuis-dev] Problems in class Number

Andres Valloud ten at smallinteger.com
Tue Oct 8 21:22:04 PDT 2019


In the larger slicer test,

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) 
isPrimeFast2e]]

I get 2627 vs 2430, or about 7.5% faster.

On 10/8/19 21:19, Andres Valloud via Cuis-dev wrote:
> "The latest code you sent"
> Time millisecondsToRun:
>      [10000 timesRepeat: [1 to: 1000 do: [:x | x isPrimeFast1b]]] 767
> 
> "The code from my last email"
> Time millisecondsToRun:
>      [10000 timesRepeat: [1 to: 1000 do: [:x | x isPrimeFast2e]]] 704
> 
> The observation is that the boundary of 31 is arbitrary, so we might as 
> well tune it according to the break even point.
> 
> On 10/8/19 21:11, Agustín Sansone via Cuis-dev wrote:
>> I don't think there will be any difference by making optimizations for 
>> small numbers. This runs just as fast as the original approach.
>>
>> El mié., 9 oct. 2019 a las 1:01, Andres Valloud via Cuis-dev 
>> (<cuis-dev at lists.cuis.st <mailto:cuis-dev at lists.cuis.st>>) escribió:
>>
>>     Expanding on the idea to treat tiny integers as special cases,
>>     approximating sqrtFloor for tiny integers wins.
>>
>>     On 10/8/19 20:49, Andres Valloud via Cuis-dev wrote:
>>      > See attached hybrid.
>>      >
>>      > On 10/8/19 20:44, Andres Valloud via Cuis-dev wrote:
>>      >> Right, that won't work.  I had tried to avoid doing something
>>     like this,
>>      >>
>>      >>      | mod30Index |
>>      >>      self < 3 ifTrue: [^self = 2].
>>      >>      self < 32 ifTrue: [
>>      >>          ^#(false true true false true false true false false 
>> false
>>      >>              true false true false false false true false true 
>> false
>>      >>              false false true false false false false false true
>>     false
>>      >>              true) at: self].
>>      >>      mod30Index := self \\ 30 + 1.
>>      >>      #(false true false false false false false true false false
>>      >>          false true false true false false false true false true
>>      >>          false false false true false false false false false 
>> true)
>>      >>              at: mod30Index :: ifFalse: [^false].
>>      >>
>>      >>
>>      >> but alas it's not as simple as I thought.
>>      >>
>>      >> Andres.
>>      >>
>>      >> On 10/8/19 20:40, Agustín Sansone via Cuis-dev wrote:
>>      >>> Sorry, I think this does not work for the numbers 3, 5, 7, 
>> 11, 13,
>>      >>> 17, 19, 23, 29 and 31.
>>      >>>
>>      >>> El mié., 9 oct. 2019 a las 0:34, Andres Valloud via Cuis-dev
>>      >>> (<cuis-dev at lists.cuis.st <mailto:cuis-dev at lists.cuis.st>
>>     <mailto:cuis-dev at lists.cuis.st <mailto:cuis-dev at lists.cuis.st>>>)
>>     escribió:
>>      >>>
>>      >>>     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> <mailto:cuis-dev at lists.cuis.st
>>     <mailto:cuis-dev at lists.cuis.st>>
>>      >>>     <mailto:cuis-dev at lists.cuis.st
>>     <mailto:cuis-dev at lists.cuis.st> <mailto: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>
>>     <http://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>
>>     <mailto:Cuis-dev at lists.cuis.st <mailto:Cuis-dev at lists.cuis.st>>
>>      >>>     <mailto:Cuis-dev at lists.cuis.st
>>     <mailto:Cuis-dev at lists.cuis.st> <mailto:Cuis-dev at lists.cuis.st
>>     <mailto:Cuis-dev at lists.cuis.st>>>
>>      >>>      >> https://lists.cuis.st/mailman/listinfo/cuis-dev
>>      >>>      >>
>>      >>>      >>
>>      >>>     --     Cuis-dev mailing list
>>      >>> Cuis-dev at lists.cuis.st <mailto:Cuis-dev at lists.cuis.st>
>>     <mailto:Cuis-dev at lists.cuis.st <mailto:Cuis-dev at lists.cuis.st>>
>>      >>> https://lists.cuis.st/mailman/listinfo/cuis-dev
>>      >>>
>>      >>>
>>      >
>>     --     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
>>
>>


More information about the Cuis-dev mailing list