[Cuis-dev] Problems in class Number

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


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>>) 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
>>
>>
-------------- next part --------------
'From Cuis 5.0 [latest update: #3866] on 8 October 2019 at 8:47:30 pm'!

!Integer methodsFor: 'testing' stamp: 'sqr 10/8/2019 20:46:35'!
isPrimeFast2d

	self < 3 ifTrue: [ ^self = 2 ].
	self even ifTrue: [^false].
	self < 32 ifTrue: [
		3 to: self - 1 by: 2 do: [:each | 
			self \\ each = 0 ifTrue: [ ^false ].
		].
		^true
	].
	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