FB II Compiler

PG PRO

Debugging

Memory

System

Mathematics

Resources

Disk I/O

Windows

Controls

Mouse

Keyboard

Text

Fonts

Drawing

Sound

Clipboard

Printing

Communication

ASM

# MATHEMATICS

## Find prime numbers

See how this does for you:

```long if myNumber% / 2 <> int(myNumber% / 2)   isPrime% = _true 'or that _ztrue thing if you like   for loop% = 3 to myNumber% / 2 step 2 'Skip the even numbers     long if myNumber% / loop% = int (myNumber% / loop%) 'Is it a factor?        isPrime% = _false 'not a prime!        loop% = myNumber% / 2 'Will this dump us out of the loop ok? (optional)     end if   next loop% xelse   isPrime% = _false 'it was an even number   if myNumber% = 1 then isPrime% = _true 'almost forgot this one! 1 is prime, right? end if```

Paul Bruneau

At Last, something I can answer !

```LOCAL FN IsItAPrime% (Num&)   IsPrime% = _True   LONG IF (Num&+1 AND 1)     IsPrime% = _False   XELSE     HalfWay& = Num& / 2     div& = 2     DO      INC(div&)      IF Num& MOD div& = 0 THEN IsPrime% = _False     UNTIL IsPrime% = _False OR div& > HalfWay&   END IF END FN = IsPrime%```

Ian Mann

Try this. I just wrote it up. It works for me.

Tedd

```COMPILE 0,_dimmedVarsOnly END GLOBALS '-------------------------- LOCAL FN buildWind   WINDOW 1,"test",(0,0)-(500,500),_doczoom END FN '-------------------------- ' is a number prime? LOCAL FN isPrime(num)   DIM i   DIM test   test = _true   PRINT "Number is = ";num   FOR i = 2 TO num -1     LONG IF num MOD i = 0       PRINT "Divisors = ";i       test = _false     END IF   NEXT END FN = test '-------------------------- ' find the primes within a number LOCAL FN findPrime(num)   DIM i, k   DIM test   DIM a   PRINT "Number is = ";num   FOR k = 1 TO num     test = _true     FOR i = 2 TO k-1       LONG IF k MOD i = 0         test = _false       END IF     NEXT     LONG IF test = _true       PRINT "Primes are = ";i     END IF   NEXT END FN '-------------------------- "MAIN" DIM i DIM a\$ DIM num DIM test FN buildWind num = 242 'check with a non-prime number test = FN isPrime (num) LONG IF test = _true   PRINT num;" is prime" XELSE   PRINT num;" is not prime" END IF PRINT num = 11 'check with a prime number test = FN isPrime (num) LONG IF test = _true   PRINT num;" is prime" XELSE   PRINT num;" is not prime" END IF PRINT FN findPrime(100) 'find all primes wihtin this number PRINT INPUT"Enter anything to end";a\$ END```

Why don't you, once you have the prime numbers, write them out as a binary string, with the bit set if the corresponding number is a prime -- i.e., the first bit is 0, because "0" is not a prime, the second is 0 because "1" is not a prime, the third and fourth are set to 1, because "2" and "3" are prime, etc.
This takes only 832 bits, is 104 bytes, plus an FN BITTEST, surely less than a complete function, and very fast...

Hans van Maanen

Cool suggestion. I, however, took it as a challenge. I couldn't get my code down to 104 bytes, but at 146 byte it's not _that_ far off. This works for numbers from 2 to 832:

```LOCAL FN IsItAPrime (Num&)   FOR factor& = Num&/2 TO 2 STEP -1     IF Num& MOD factor& = 0 THEN factor& = - factor&   NEXT END FN = factor& > -3```

Of course your suggestion blows this away speed-wise. :-)

Jay

I noticed most people's algorithms were checking divisors all the way up to num&/2. It's only necessary to check up to SQR(num&):

```LOCAL FN IsItPrime(num&)   SELECT CASE     CASE num& = 2 OR num& = 3       prime = _zTrue     CASE num& MOD 2 = 0       prime = _false     CASE ELSE       prime = _zTrue '(initially)       root = SQR(num&)       FOR i = 3 TO root STEP 2         LONG IF num& MOD i = 0           prime = _false           i = root 'force early exit from loop         END IF       NEXT   END SELECT END FN = prime```

In fact, it's only necessary to check _prime_ potential divisors up to SQR(num&), which further shortens the list. So to check whether 832 is prime, you only need to see if it's divisible by any of these numbers: 2,3,5,7,11,13,17,19,23. My short function above doesn't have that optimization, but you could potentially use an array of prime numbers to use as test divisors.

Rick

There's no contest. Here's a short program comparing my 5-line FN (which probably is not the fastest) with Hans's BITTST suggestion. Evaluating 10,000 random numbers in the range 1-832 with my code takes my machine at least 155 ticks. Looking up the same 10,000 numbers in Hans's list takes less than 3 ticks!. (BTW, be prepared for a time increase of 10-20x if you run this in the debugger.) :-)

Jay

```GLOBALS _trials = 10000 'Set number of trials here DIM gPrimeBits.104 DIM randomList(_trials) DIM gTime& DIM r&,isPrime,d\$ END GLOBALS LOCAL FN IsItAPrime (Num&) 'Jay's FN   'Only works for integers > 1   FOR factor& = Num&/2 TO 2 STEP -1     IF Num& MOD factor& = 0 THEN factor& = - factor&   NEXT END FN = factor& > -3 LOCAL FN makeLists   FOR r = 2 TO 832 'Set up Hans' list using Jay's FN     IF FN IsItAPrime(r) THEN CALL BITSET(#@gPrimeBits,r)   NEXT   FOR r = 1 TO _trials 'Choose #'s to test     randomList(r) = RND(831)+1   NEXT END FN '------------Main--------------- FN makeLists CLS PRINT "For";_trials;"trials..." gTime& = FN TICKCOUNT+1 WHILE gTime& < FN TICKCOUNT WEND FOR r& = 1 TO _trials   isPrime = FN IsItAPrime(randomList(r&)) NEXT PRINT "The code took"; FN TICKCOUNT-gTime&;"ticks." gTime& = FN TICKCOUNT+1 WHILE gTime& < FN TICKCOUNT:WEND FOR r& = 1 TO _trials   isPrime = -FN BITTST(#@gPrimeBits,randomList(r&)) NEXT PRINT "Bit test took";FN TICKCOUNT-gTime&;"ticks." PRINT CURSOR _arrowcursor COLOR _zred INPUT "Press Return to end"; d\$ PRINT```

I claim (at least temporarily) the prize. The program below determines primeness and factors of any number from 1 to 4294967295 which is the maximum value of an unsigned LONG. Rick's routine, posted last week, works up to 2147483647. In a head-to-head test (on all numbers from 1 to 1 million), my routine is 11 times faster, partly because I replaced SQR(num&) by a very fast integer square root. Along the way it turned out that the obvious replacement for SQR, USR _sqRoot, has some defects, which are remedied in FN SquareRoot&.

BTW, Tedd, isn't it time you set your date to 1998? Are you are hoping to postpone Y2K to 2001? ;-)

Robert Purves

```LOCAL FN SquareRoot&(a&) ' Replacement for USR _sqRoot (a&) which has three faults:- ' 1. It crashes for a&=_Maxlong/2 ' 2. It gives wrong answer for a&>_Maxlong/2 ' 3. it gves "high-by-one" answers in many cases, e.g. USR _sqRoot(15)=4 ' This routine works for all values of a& unsigned (0-4294967295) ' and is slightly faster than USR _sqRoot. It DOES NOT WORK on 68000 machines. DIM root& ` move.l ^a&,d3 ` move.l d3,d0 ` beq.s SqDone ; skip if 0 ` move.l d3,d1 `ApprLoop ; init root by halving the number of digits in a& ` lsr.l #1,d0 ` lsr.l #2,d1 ` bne.s ApprLoop ` addq.l #1,d0 ; force non-zero ` move.w #4,d7 ; max loop count `SqLoop ; iterate root& = (a&/root&+root& )/2 ` move.l d0,d5 ; copy of root ` clr.l d1 ; hi.l of a&=0 ` move.l d3,d2 ; lo.l of a& ` dc.w \$4C40 ; DIVU.L D0,D1:D2 ` dc.w \$2601 ; d2.l= d1(hi) d2(lo) / d0.l d1.l = remainder ` add.l d2,d0 ` lsr.l #1,d0 ` cmp.l d0,d5 ; equals old? If so, we are done ` dbeq d7, SqLoop ; branch back if <> old and loop count still >0 ' IF root&*root& > a& THEN dec(root&) ` move.l d0,d2 ` dc.w \$4c00 ; MULU.L D0,D1:D2 ` dc.w \$2401 ;d1 (hi) and d2 (lo) = d0.l * d2.l ` cmp.l d2,d3 ` bge.s SqDone ` subq.l #1,d0 `SqDone ` move.l d0,^root& END FN=root& LOCAL FN IsItPrimeRDP(num&,fact1Ptr&,fact2Ptr&) ' Brute force test for primeness ' Works for unsigned num& 0<=num&<= 2^32-1 (4294967295) ' DOES NOT WORK on 68000 machines DIM prime SELECT CASE   CASE num& = 2 OR num& = 3    prime = _zTrue    POKE LONG fact1Ptr&,1: POKE LONG fact2Ptr&,num&   CASE (num& AND 1) = 0' even    prime = _false    POKE LONG fact1Ptr&,2: POKE LONG fact2Ptr&,num&>>1   CASE ELSE    REGISTER(d3)=FN SquareRoot&(num&)    ` move.l ^num&,d4    ` move.w #3,d0 ; i=3 divisor    `loop    ` clr.l d1 ; hi of num&=0    ` move.l d4,d2 ; lo of num&    ` dc.w \$4C40 ; DIVU.L D0,D1:D2    ` dc.w \$2601 ; d2.l= d1(hi) d2(lo) / d0.l d1.l = remainder    ` tst.l d1    ` beq.s notPrime ; branch if remainder=0    ` addq.w #2,d0    ` cmp.w d3,d0    ` bls.s loop    ` move.w #-1,d1 ; flag _ztrue    ` move.l #1, d0    ` move.l ^num&,d2    `notPrime    ` movea.l ^fact1Ptr&,a0    ` ext.l d0    ` move.l d0,(a0)    ` movea.l ^fact2Ptr&,a0    ` move.l d2,(a0)    ` move.w d1,^prime    `done END SELECT END FN = prime LOCAL FN StringToLongUnsigned&(s\$) ' Convert number up to 4294967295 DIM j,a&, 1 char\$ a&=0 FOR j=1 TO LEN(s\$)   char\$=MID\$(s\$,j,1)   LONG IF char\$>="0" AND char\$<="9"'ignore non-numeric    a&=a&*10 + VAL(char\$)   END IF NEXT j END FN=a& WINDOW 1 DEFSTR LONG DIM s\$,a&,f1&,f2& DO INPUT "Enter number (return to end) :"; s\$ IF s\$="" THEN END a&=FN StringToLongUnsigned&(s\$) PRINT "Number is:" VAL(UNS\$(a&)) LONG IF FN IsItPrimeRDP(a&,@f1&,@f2&)   PRINT "P"; XELSE   PRINT "Not p"; END IF PRINT "rime. Factors: " VAL(UNS\$(f1&))"and"VAL(UNS\$(f2&)) UNTIL 0```

Rick's and Hans' suggestions are actually equivalent, that is they lead to the same sqeuence of divisors. Nice call, guys.
And by extension we can remove multiples of 5 from the list of divisors as in the following program. It's in plain vanilla FB for comprehensibility.
Even though the previous version contained screeds of assembly code, this one is actually faster. The down-side of no ASM is that it can't factor numbers that are >= 1073741823

```COMPILE _dimmedVarsOnly DIM gIncr(7) END GLOBALS LOCAL FN InitIncrements gIncr(0)=4: gIncr(1)=2: gIncr(2)=4: gIncr(3)=2 gIncr(4)=4: gIncr(5)=6: gIncr(6)=2:gIncr(7)=6 END FN LOCAL FN RDPIsItPrime2(num&,fact1Ptr&,fact2Ptr&) DIM prime,root&,divisor&,index prime = _false _maxLongBy2=1073741823 SELECT CASE   CASE num&=0    & fact1Ptr&,0: & fact2Ptr&,0   CASE num&=1    & fact1Ptr&,1: & fact2Ptr&,1   CASE num&>=_maxLongBy2 OR num&<0    PRINT "Too big for USR _sqRoot "    & fact1Ptr&,0: & fact2Ptr&,0   CASE num& = 2 OR num& = 3 OR num&=5    prime = _zTrue    & fact1Ptr&,1: & fact2Ptr&,num&   CASE (num& AND 1) = 0' even    & fact1Ptr&,2: & fact2Ptr&,num&>>1   CASE num& MOD 3 = 0' divisible by 3    & fact1Ptr&,3: & fact2Ptr&,num&/3   CASE num& MOD 5 = 0' divisible by 5    & fact1Ptr&,5: & fact2Ptr&,num&/5   CASE ELSE    prime = _zTrue 'initially    & fact1Ptr&,1: & fact2Ptr&,num&    ' Start with divisor 7 and add {4,2,4,2,4,6,2,6} in sequence,    ' thus skipping divisors that are multiples of 2,3 or 5    root& =USR _sqRoot (num&)    divisor&=7    index=0    WHILE divisor&<=root&     LONG IF num& MOD divisor&= 0       prime = _false       & fact1Ptr&,divisor&: & fact2Ptr&,num&/divisor&       divisor=root&'force early exit     END IF     divisor&=divisor&+gIncr(index)     index=(index+1) AND 7    WEND END SELECT END FN = prime WINDOW 1 DEFSTR LONG DIM s\$,a&,f1&,f2& FN InitIncrements DO   INPUT "Enter number (return to end) :"; s\$   IF s\$="" THEN END   a&=VAL(s\$)   PRINT "Number is:" VAL(UNS\$(a&))   LONG IF FN RDPIsItPrime2(a&,@f1&,@f2&)     PRINT "P";   XELSE     PRINT "Not p";   END IF   PRINT "rime. Factors: " VAL(UNS\$(f1&))"and"VAL(UNS\$(f2&)) UNTIL 0```

Robert