title 'Eratosthenes Sieve for 80x86 Real Mode'
name sieve
page 50,132
;
; Eratosthenes Sieve for 80x86 Real Mode
; Implemented by Ray Duncan, April 1987
; After Gilbreath, Byte September 1981, and January 1983
;
niter equ 100 ; number of iterations
asize equ 16380 ; size of array "flags"
cr equ 0dh ; ASCII carriage return
lf equ 0ah ; ASCII line feed
stdin equ 0 ; handle for standard input
stdout equ 1 ; handle for standard output
_TEXT segment para public 'CODE'
assume cs:_TEXT,ds:_DATA,es:_DATA
sieve proc near
mov ax,seg _DATA
mov ds,ax
mov es,ax
mov dx,0 ; convert number of iterations
mov ax,niter
mov cx,10
mov si,offset msg1a+3
call binasc
mov dx,offset msg1 ; display message
mov cx,msg1_len ; "Starting N iterations of Sieve"
call pmsg
call getmsec ; get current time in msec
push dx ; and save it ...
push ax
mov counter,niter ; initialize iterations counter
sieve1: ; a sieve iteration starts here...
mov di,offset flags ; initialize flags array
mov cx,asize ; to all bytes = TRUE
mov al,1
cld
rep stosb
xor si,si ; SI= index to flags array
xor di,di ; DI = primes counter
sieve2: ; main loop of sieve
test byte ptr flags[si],1 ; is this a prime?
jnz short sieve4 ; jump if prime
sieve3: inc si ; bump to next slot in "flags"
cmp si,asize ; are we done?
jle sieve2 ; jump to test another
dec word ptr counter ; more iterations?
jnz sieve1 ; jump, another iteration needed
jmp sieve7
sieve4: ; prime found, zap its multiples
mov bx,si ; copy i
mov dx,bx ; DX = prime = i + i + 3
add dx,dx
add dx,3
xor al,al
jmp short sieve6
sieve5: mov byte ptr flags[bx],al ; zero this multiple
sieve6: add bx,dx ; find next multiple of prime
cmp bx,asize ; have we exhausted the array?
jle sieve5 ; not yet, zap it
inc di ; count primes and try next
jmp sieve3
sieve7: ; done with all iterations...
call getmsec ; get current time
push dx ; and save it...
push ax
mov ax,di ; convert number of primes
mov dx,0
mov cx,10
mov si,offset msg2a+4
call binasc
mov dx,offset msg2 ; display "Number of primes: "
mov cx,msg2_len
call pmsg
pop ax ; stop time: low word
pop dx ; high word
pop bx ; start time: low word
pop cx ; high word
sub ax,bx ; DX:AX = stop - start
sbb dx,cx
mov cx,niter ; divide by number of iterations
idiv cx
mov dx,0 ; convert msec to ASCII
mov cx,10
mov si,offset msg3a+4
call binasc
mov dx,offset msg3 ; display "Elapsed time:"
mov cx,msg3_len
call pmsg
mov ax,04C00h ; exit to DOS with
int 21H ; return code = 0
sieve endp
getmsec proc near ; DX:AX := current time in msec.
mov ah,2ch ; read time from MS-DOS
int 21h
push dx ; save seconds, hundredths
mov al,ch ; AX := hours
cbw
mov bx,60 ; hours -> minutes
imul bx
xor ch,ch ; isolate system minutes
add ax,cx ; and find total minutes
mov bx,60 ; minutes -> seconds
imul bx
pop cx ; recover seconds, hundredths
mov bl,ch ; get seconds
xor bh,bh
add ax,bx ; find total seconds
adc dx,0 ; carry if necessary
xor ch,ch ; save centisec.
mov bp,cx
; total seconds * 100 the hard way
mov bx,ax ; double multiply * 10
mov cx,dx
add ax,ax ; * 2
adc dx,dx
add ax,ax ; * 4
adc dx,dx
add ax,bx ; * 5
adc dx,cx
add ax,ax ; * 10
adc dx,dx
mov bx,ax ; double multiply * 10
mov cx,dx
add ax,ax ; * 2
adc dx,dx
add ax,ax ; * 4
adc dx,dx
add ax,bx ; * 5
adc dx,cx
add ax,ax ; * 10
adc dx,dx
add ax,bp ; add in hundredths of seconds
adc dx,0
; now convert total to msec.
mov bx,ax ; double multiply * 10
mov cx,dx
add ax,ax ; * 2
adc dx,dx
add ax,ax ; * 4
adc dx,dx
add ax,bx ; * 5
adc dx,cx
add ax,ax ; * 10
adc dx,dx
ret ; return DX:AX = msec.
getmsec endp
;
; BINASC: Convert 32 bit binary value to ASCII string.
;
; Call with DX:AX = signed 32 bit value
; CX = radix
; SI = last byte of area to store resulting string
; (make sure enough room is available to store
; the string in the radix you have selected.)
;
; Destroys AX, BX, CX, DX, and SI.
;
binasc proc near ; convert DX:AX to ASCII.
mov byte ptr [si],'0' ; force storage of at least 1 digit.
or dx,dx ; test sign of 32 bit value,
pushf ; and save sign on stack.
jns bin1 ; jump if it was positive.
not dx ; negative, take 2's complement
not ax ; of the value.
add ax,1
adc dx,0
bin1: ; divide 32 bit value by radix
; to extract next digit for the
; forming string.
mov bx,ax ; is the value zero yet?
or bx,dx
jz bin3 ; yes, we are done converting.
call divide ; no, divide by radix.
add bl,'0' ; convert remainder to ASCII digit.
cmp bl,'9' ; might be converting to hex ASCII,
jle bin2 ; jump if in range 0-9,
add bl,'A'-'9'-1 ; correct it if in range A-F.
bin2: mov [si],bl ; store this character into string.
dec si ; back up through string,
jmp bin1 ; and do it again.
bin3: ; restore sign flag,
popf ; was original value negative?
jns bin4 ; no, jump
mov byte ptr [si],'-' ; yes,store sign into output.
bin4: ret ; back to caller.
binasc endp
;
; General purpose 32 bit by 16 bit unsigned divide.
; This must be used instead of the plain machine unsigned divide
; for cases where the quotient may overflow 16 bits (for example,
; dividing 100,000 by 2). If called with a zero divisor, this
; routine returns the dividend unchanged and gives no warning.
;
; Call with DX:AX = 32 bit dividend
; CX = divisor
;
; Returns DX:AX = quotient
; BX = remainder
; CX = divisor (unchanged)
;
divide proc near ; Divide DX:AX by CX
jcxz div1 ; exit if divide by zero
push ax ; 0:dividend_upper/divisor
mov ax,dx
xor dx,dx
div cx
mov bx,ax ; BX = quotient1
pop ax ; remainder1:dividend_lower/divisor
div cx
xchg bx,dx ; DX:AX = quotient1:quotient2
div1: ret ; BX = remainder2
divide endp
pmsg proc near ; print a message on std output
; call with DS:EDX = address
; ECX = length
mov ah,40h
mov bx,stdout
int 21h
ret
pmsg endp
_TEXT ends
_DATA segment para public 'DATA'
flags db asize+1 dup (?)
counter dw ?
msg1 db cr,lf,'Starting '
msg1a db ' iterations of Sieve...',cr,lf
msg1_len equ $-msg1
msg2 db cr,lf,'Primes found: '
msg2a db ' ',cr,lf
msg2_len equ $-msg2
msg3 db cr,lf,'Elapsed time: '
msg3a db ' msec. per iteration',cr,lf
msg3_len equ $-msg3
_DATA ends
_STACK segment byte stack 'stack'
db 4096 dup (?)
_STACK ends
end sieve
