'Subject: Sieve of Erathostenes: a variation of Chartres' Algorithm 311.
'         Generates PrimFlgs.bin, an encoded list of prime numbers
'         read by library function Nxtprm(). The prime / not prime
'         information is stored for 90 integers in 4 bytes.
'Author : Sjoerd.J.Schaper
'Date   : 07-01-2005
'Code   : QuickBasic, PDS, VBdos, FreeBasic
REM $LANG: "qb"
DEFLNG A-Z
DIM dword AS STRING * 4
DIM bi(7 TO 95), bo(7 TO 95)
CLS : tim# = TIMER: lin = 7
Slate = FREEFILE

CONST Fx = 2948947 '                 set Quickbasic or
'CONST Fx = 2147483629 '              Freebasic limit

N = 1: dt = 4: t = 5
FOR k = 1 TO 30
   dt = 6 - dt: t = t + dt
   bi(t) = N: bo(t) = NOT N '         store masks
   N = N * 2
NEXT k

begin:
LOCATE 1, 1: f$ = "0"
PRINT "   +---------------------+"
PRINT "   |  input primelimit,  |"
PRINT "   |  0 to exit.         |"
PRINT "   |                     |"
PRINT "   +---------------------+"
LOCATE 4, 7: INPUT ; "l"; il$
IF INSTR(il$, "'") THEN GOTO begin
IF VAL(il$) < 1 THEN GOTO eind
IF VAL(il$) > Fx THEN il$ = ""
lP = VAL(il$)

LOCATE lin, 2: PRINT "Filling slate..."
IF lP = 0 THEN
   lP = Fx: PRINT "limit set to"; lP
END IF
lv = (lP - 7) \ 90
DIM SHARED bv(lv) '                   preprocessed bit vector
bv(0) = &HB76BDBF '                   5 & 7-folds excluded
FOR i = 1 TO lv - 1
   bv(i) = &H1BF6FDBF '               5-folds excluded
NEXT i

PRINT " Erasing multiples..."
li = (SQR(lP) - 7) \ 90: dt = 4
FOR i = 0 TO li
   c = bv(i): t = 5 '                 read 30-bit chunk
   FOR k = 1 TO 30
      dt = 6 - dt: t = t + dt '       skip multiples of 2 and 3
      IF c AND bi(t) THEN
         P = i * 90 + t '             calculate prime,
         dP = dt * P: P6 = 6 * P '     steps
         N = P * P '                  initial N = k * P
         WHILE N < lP
            j = (N - 7) \ 90 '        point to multiple N
            m = N - j * 90
            bv(j) = bv(j) AND bo(m) ' clear bit
            dP = P6 - dP: N = N + dP ' step to next N
         WEND
      END IF
   NEXT k
NEXT i
f$ = STR$(CSNG(TIMER - tim#))

PRINT " Writing slate to disk..."
g$ = LTRIM$(RTRIM$(ENVIRON$("LARGEINT")))
IF g$ = "" THEN g$ = CURDIR$ '        working directory
IF RIGHT$(g$, 1) <> "\" THEN g$ = g$ + "\"
OPEN g$ + "primflgs.bin" FOR BINARY SHARED AS Slate
FOR i = 0 TO lv - 1
   dword = MKL$(bv(i))
   PUT #Slate, , dword
NEXT i
PRINT " Done": PRINT
lin = CSRLIN

eind:
CLOSE Slate
OPEN g$ + "primflgs.log" FOR OUTPUT AS Slate
f$ = "Timer:" + f$ + " s"
LOCATE lin, 1: PRINT f$
PRINT #Slate, f$
CLOSE
SYSTEM

