[Cuis-dev] Problems in class Number

Andres Valloud ten at smallinteger.com
Tue Oct 8 20:44:10 PDT 2019


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 that 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>>) 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>>>)
>     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>
>      >> 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>>
>      >> 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