fp.asm

汇编编程艺术

;
StdGrp		group	StdLib, StdData
;
StdData		segment	para public 'sldata'
;
; Floating point package.
;
;
; Released to the public domain
; Created by: Randall Hyde
; Date: 8/13/90
;	8/28/91
;
;
; FP format:
;
; 80 bits:
; bit 79            bit 63                           bit 0
; |                 |                                    |
; seeeeeee eeeeeeee mmmmmmmm m...m m...m m...m m...m m...m
;
; e = bias 16384 exponent
; m = 64 bit mantissa with NO implied bit!
; s = sign (for mantissa)
;
;
; 64 bits:
; bit 63       bit 51                                               bit 0
; |            |                                                        |
; seeeeeee eeeemmmm mmmmmmmm mmmmmmmm mmmmmmmm mmmmmmmm mmmmmmmm mmmmmmmm
;
; e = bias 1023 exponent.
; s = sign bit.
; m = mantissa bits.  Bit 52 is an implied one bit.
;
; 32 bits:
; Bit 31    Bit 22              Bit 0
; |         |                       |
; seeeeeee emmmmmmm mmmmmmmm mmmmmmmm
;
; e = bias 127 exponent
; s = sign bit
; m = mantissa bits, bit 23 is an implied one bit.
;
;
;
; WARNING: Although this package uses IEEE format floating point numbers,
;	   it is by no means IEEE compliant.  In particular, it does not
;	   support denormalized numbers, special rounding options, and
;	   so on.  Why not?  Two reasons:  I'm lazy and I'm ignorant.
;	   I do not know all the little details surround the IEEE
;	   implementation and I'm not willing to spend more of my life
;	   (than I already have) figuring it out.  There are more
;	   important things to do in life.  Yep, numerical analysts can
;	   rip this stuff to shreads and come up with all kinds of degenerate
;	   cases where this package fails and the IEEE algorithms succeed,
;	   however, such cases are very rare.  One should not get the idea
;	   that IEEE is perfect.  It blows up with lots of degenerate cases
;	   too.  They just designed it so that it handles a few additional
;	   cases that mediocre packages (like this one) do not.  For most
;	   normal computations this package works just fine (what it lacks
;	   it good algorithms it more than makes up for by using an 88-bit
;	   internal format during internal computations).
;
;	   Moral of the story: If you need highly accurate routines which
;          produce okay results in the worst of cases, look elsewhere please.
;	   I don't want to be responsible for your blowups.  OTOH, if you need
;	   a fast floating point package which is reasonably accurate and
;	   you're not a statistician, astronomer, or other type for whom
;	   features like denormalized numbers are important, this package
;	   may work out just fine for you.
;
;						Randy Hyde
;						August 1990
;						(Hard to believe I started this
;						 a year ago and I'm just coming
;						 back to it now!)
;
;						UC Riverside &
;						Cal Poly Pomona.
;
; FPACC- Floating point accumuator.
; FPOP-  Floating point operand.
;
; These variables use the following format:
;
; 88 bits:
; sxxxxxxx eeeeeeee eeeeeeee m..m m..m m..m m..m m..m m..m m..m m..m
; Sign          exponent                   mantissa (64 bits)
;
; Only H.O. bit of Sign byte is significant.  The rest is garbage.
; Exponent is bias 32767 exponent.
; Mantissa does NOT have an implied one bit.
;
; This format was picked for convenience (it is easy to work with) and it
; exceeds the 80-bit format used by Intel on the 80x87 chips.
;
fptype		struc
Mantissa	dw	4 dup (?)
Exponent	dw	?
Sign		db	?
		db	?		;Padding
fptype		ends
;
;
;
;
		public	fpacc
fpacc		fptype	<>
;
		public	fpop
fpop		fptype  <>
;
;
; FProd- Holds 144-bit result obtained by multiplying fpacc.mant x fpop.mant
;
Quotient	equ	this word
fprod		dw	9 dup (?)
;
;
; Variables used by the floating point I/O routines:
;
TempExp		dw	?
ExpSign		db	?
DecExponent	dw	?
DecSign		db	0
DecDigits	db	31 dup (?)
;
;
;
StdData		ends
;
;
stdlib		segment	para public 'slcode'
		assume	cs:stdgrp, ds:nothing, es:nothing, ss:nothing
;
;
;
;
;
;
;
;
;
;
;---------------------------------------------------------------------------
;		Floating Point Load/Store Routines
;---------------------------------------------------------------------------
;
; sl_AccOp	Copies the floating point accumulator to the floating point
;		operand.
;
		public	sl_AccOp
sl_AccOp	proc	far
		assume	ds:StdGrp
		push	ax
		push	ds
		mov	ax, StdGrp
		mov	ds, ax
;
		mov	ax, FPacc.Exponent
		mov	FPop.Exponent, ax
		mov	ax, FPacc.Mantissa
		mov	FPop.Mantissa, ax
		mov	ax, FPacc.Mantissa+2
		mov	FPop.Mantissa+2, ax
		mov	ax, FPacc.Mantissa+4
		mov	FPop.Mantissa+4, ax
		mov	ax, FPacc.Mantissa+6
		mov	FPop.Mantissa+6, ax
		mov	al, Fpacc.Sign
		mov	FPop.Sign, al
;
		pop	ds
		pop	ax
		ret
sl_AccOp	endp
		assume	ds:nothing
;
;
; sl_XAccOp-	Exchanges the values in the floating point accumulator
;		and floating point operand.
;
		public	sl_XAccOp
sl_XAccOp	proc	far
		assume	ds:StdGrp
		push	ax
		push	ds
		mov	ax, StdGrp
		mov	ds, ax
;
		mov	ax, FPacc.Exponent
		xchg	ax, FPop.Exponent
		mov	FPacc.Exponent, ax
;
		mov	ax, FPacc.Mantissa
		xchg	ax, FPop.Mantissa
		mov	FPacc.Mantissa, ax
;
		mov	ax, FPacc.Mantissa+2
		xchg	ax, FPop.Mantissa+2
		mov	FPacc.Mantissa+2, ax
;
		mov	ax, FPacc.Mantissa+4
		xchg	ax, FPop.Mantissa+4
		mov	FPacc.Mantissa+4, ax
;
		mov	ax, FPacc.Mantissa+6
		xchg	ax, FPop.Mantissa+6
		mov	FPacc.Mantissa+6, ax
;
		mov	al, FPacc.Sign
		xchg	al, FPop.Sign
		mov	FPacc.Sign, al
;
		pop	ds
		pop	ax
		ret
sl_XAccOp	endp
		assume	ds:nothing
;
;
;
; sl_LSFPA- 	Loads a single precision (32-bit) IEEE format number into
;		the floating point accumulator.  ES:DI points at the # to
;		load into FPACC.
;
		public	sl_LSFPA
sl_LSFPA	proc	far
		push	ax
		push	bx
		mov	ax, es:[di]
		mov	word ptr StdGrp:fpacc.mantissa[5], ax
		mov	ax, es:2[di]
		mov	bx, ax
		shl	ax, 1
		mov	al, ah
		mov	ah, 0
		add	ax, 32767-127		;Adjust exponent bias.
		mov	word ptr StdGrp:fpacc.exponent, ax
		mov	StdGrp:fpacc.sign, bh	;Save sign away.
		mov	al, es:2[di]
		and	al, 7fh			;Strip out L.O. exp bit.
		or	al, 80h			;Add in implied bit.
		mov	byte ptr StdGrp:fpacc.mantissa[7], al ;Save H.O. mant byte.
		xor	ax, ax
		mov	word ptr StdGrp:fpacc.mantissa, ax
		mov	word ptr StdGrp:fpacc.mantissa[2], ax
		mov	byte ptr StdGrp:fpacc.mantissa[4], al
		pop	bx
		pop	ax
		ret
sl_LSFPA	endp
;
;
;
;
; sl_SSFPA-	Stores FPACC into the single precision variable pointed at by
;		ES:DI.  Performs appropriate rounding.  Returns carry clear
;		if the operation is successful, returns carry set if FPACC
;		cannot fit into a single precision variable.
;
		public	sl_SSFPA
sl_SSFPA	proc	far
		assume	ds:stdgrp
		push	ds
		push	ax
		push	bx
		mov	ax, StdGrp
		mov	ds, ax
		push	fpacc.Exponent
		push	fpacc.Mantissa   	;Save the stuff we tweak
		push	fpacc.Mantissa[2]	; so that this operation
		push	fpacc.Mantissa[4]	; will be non-destructive.
		push	fpacc.Mantissa[6]
;
; First, round FPACC:
;
		add	fpacc.Mantissa [4], 80h
		adc	fpacc.Mantissa [6], 0
		jnc	StoreAway
		rcl	fpacc.Mantissa [6], 1
		rcl	fpacc.Mantissa [4], 1
		inc	fpacc.Exponent
		jz	BadSSFPA		;If exp overflows.
;
; Store the value away:
;
StoreAway:	mov	ax, fpacc.Exponent
		sub	ax, 32767-127		;Convert to bias 127
		cmp	ah, 0
		jne	BadSSFPA
		mov	bl, fpacc.Sign
		shl	bl, 1			;Merge in the sign bit.
		rcr	al, 1
		mov	es:[di] + 3, al		;Save away exponent/sign
		pushf				;Save bit shifted out.
		mov	ax, fpacc.Mantissa [6]
		shl	ax, 1			;Get rid of implied bit and
		popf				; shift in the L.O. exponent
		rcr	ax, 1			; bit.
		mov	es:[di] + 1, ax
		mov	al, byte ptr fpacc.Mantissa [5]
		mov	es:[di], al
		clc
		jmp	SSFPADone
;
BadSSFPA:	stc
SSFPADone:	pop	fpacc.Mantissa[6]
		pop	fpacc.Mantissa[4]
		pop	fpacc.Mantissa[2]
		pop	fpacc.Mantissa
		pop	fpacc.Exponent
		pop	bx
		pop	ax
		pop	ds
		ret
		assume	ds:nothing
sl_SSFPA	endp
;
;
; sl_LDFPA-	Loads the double precision (64-bit) IEEE format number pointed
;		at by ES:DI into FPACC.
;
		public	sl_LDFPA
sl_LDFPA	proc	far
		push	ax
		push	bx
		push	cx
		mov	ax, es:6[di]
		mov	StdGrp:fpacc.sign, ah	;Save sign bit.
		mov	cl, 4
		shr	ax, cl			;Align exponent field.
		and	ah, 111b		;Strip the sign bit.
		add	ax, 32767-1023		;Adjust bias
		mov	StdGrp:fpacc.exponent, ax
;
; Get the mantissa bits and left justify them in the FPACC.
;
		mov	ax, es:5[di]
		and	ax, 0fffh		;Strip exponent bits.
		or	ah, 10h			;Add in implied bit.
		mov	cl, 3
		shl	ax, cl
		mov	bx, es:3[di]
		rol	bx, cl
		mov	ch, bl
		and	ch, 7
		or	al, ch
		mov	StdGrp:fpacc.mantissa[6], ax
;
		and	bl, 0f8h
		mov	ax, es:1[di]
		rol	ax, cl
		mov	ch, al
		and	ch, 7
		or	bl, ch
		mov	StdGrp:fpacc.mantissa[4], bx
;
		and	al, 0f8h
		mov	bh, es:[di]
		rol	bh, cl
		mov	ch, bh
		and	ch, 7
		or	al, ch
		mov	StdGrp:fpacc.mantissa[2], ax
		and	bh, 0f8h
		mov	bl, 0
		mov	StdGrp:fpacc.Mantissa[0], bx
;
		pop	cx
		pop	bx
		pop	ax
		ret
sl_LDFPA	endp
;
;
;
;
; sl_SDFPA-	Stores FPACC into the double precision variable pointed
;		at by ES:DI.
;
		public	sl_sdfpa
sl_SDFPA	proc	far
		assume	ds:stdgrp
		push	ds
		push	ax
		push	bx
		push	cx
		push	dx
		push	di
;
		mov	bx, StdGrp
		mov	ds, bx
;
		push	fpacc.Mantissa [0]
		push	fpacc.Mantissa [2]
		push	fpacc.Mantissa [4]
		push	fpacc.Mantissa [6]
		push	fpacc.Exponent
;
; First, round this guy to 52 bits:
;
		add	byte ptr fpacc.Mantissa [1], 8
		jnc	SkipRndUp
		inc	fpacc.Mantissa [2]
		jnz	SkipRndUp
		inc	fpacc.Mantissa [4]
		jnz	SkipRndUp
		inc	fpacc.Mantissa [6]
		jnz	SkipRndUp
;
; Whoops!  Got an overflow, fix that here:
;
		stc
		rcr	fpacc.Mantissa [6], 1
		rcr	fpacc.Mantissa [4], 1
		rcr	fpacc.Mantissa [2], 1
		rcr	byte ptr fpacc.Mantissa [1], 1
		inc	fpacc.Exponent
		jz	BadSDFPA		;In case exp was really big.
;
; Okay, adjust and store the exponent-
;
SkipRndUp:	mov	ax, fpacc.Exponent
		sub	ax, 32767-1023		;Adjust bias
		cmp	ax, 2048		;Make sure the value will still
		jae	BadSDFPA		; fit in an 8-byte real.
		mov	cl, 5
		shl	ax, cl			;Move exponent into place.
		mov	bl, fpacc.Sign
		shl	bl, 1
		rcr	ax, 1			;Merge in sign bit.
;
; Merge in the upper four bits of the Mantissa (don't forget that the H.O.
; Mantissa bit is lost due to the implied one bit).
;
		mov	bl, byte ptr fpacc.Mantissa [7]
		shr	bl, 1
		shr	bl, 1
		shr	bl, 1
		and	bl, 0fh			;Strip away H.O. mant bit.
		or	al, bl
		mov	es:[di]+6, ax		;Store away H.O. word.
;
; Okay, now adjust and store away the rest of the mantissa:
;
		mov	ax, fpacc.Mantissa [0]
		mov	bx, fpacc.Mantissa [2]
		mov	cx, fpacc.Mantissa [4]
		mov	dx, fpacc.Mantissa [6]
;
; Shift the bits to their appropriate places (to the left five bits):
;
		shl	ax, 1
		rcl	bx, 1
		rcl	cx, 1
		rcl	dx, 1
;
		shl	ax, 1
		rcl	bx, 1
		rcl	cx, 1
		rcl	dx, 1
;
		shl	ax, 1
		rcl	bx, 1
		rcl	cx, 1
		rcl	dx, 1
;
		shl	ax, 1
		rcl	bx, 1
		rcl	cx, 1
		rcl	dx, 1
;
		shl	ax, 1
		rcl	bx, 1
		rcl	cx, 1
		rcl	dx, 1
;
; Store away the results:
;
		mov	es:[di], bx
		mov	es:[di] + 2, cx
		mov	es: [di] + 4, dx
;
; Okay, we're done.  Return carry clear to denote success.
;
		clc
		jmp	short QuitSDFPA
;
BadSDFPA:	stc				;If an error occurred.
QuitSDFPA:	pop	fpacc.Exponent
		pop	fpacc.Mantissa [6]
		pop	fpacc.Mantissa [4]
		pop	fpacc.Mantissa [2]
		pop	fpacc.Mantissa [0]
		pop	di
		pop	dx
		pop	cx
		pop	bx
		pop	ax
		pop	ds
		ret
;
		assume	ds:nothing
sl_SDFPA	endp
;
;
;
;
; sl_LEFPA-	Loads an extended precision (80-bit) IEEE format number
;		into the floating point accumulator.  ES:DI points at the
;		number to load into FPACC.
;
		public	sl_LEFPA
sl_LEFPA	proc	far
		push	ax
		mov	ax, es:8[di]
		mov	StdGrp:fpacc.Sign, ah
		and 	ah, 7fh
		add	ax, 4000h
		mov	StdGrp:fpacc.Exponent, ax
		mov	ax, es:[di]
		mov	StdGrp:fpacc.Mantissa, ax
		mov	ax, es:2[di]
		mov	StdGrp:fpacc.Mantissa[2], ax
		mov	ax, es:4[di]
		mov	StdGrp:fpacc.Mantissa[4], ax
		mov	ax, es:6[di]
		mov	StdGrp:fpacc.Mantissa[6], ax
		pop	ax
		ret
sl_LEFPA	endp
;
;
; sl_LEFPAL-	Loads an extended precision (80-bit) IEEE format number
;		into the floating point accumulator.  The number to load
;		into FPACC follows the call in the code stream.
;
		public	sl_LEFPAL
sl_LEFPAL	proc	far
		push	bp
		mov	bp, sp
		push	es
		push	di
		push	ax
		les	di, 2[bp]
;
		mov	ax, es:8[di]
		mov	StdGrp:fpacc.Sign, ah
		and 	ah, 7fh
		add	ax, 4000h
		mov	StdGrp:fpacc.Exponent, ax
		mov	ax, es:[di]
		mov	StdGrp:fpacc.Mantissa, ax
		mov	ax, es:2[di]
		mov	StdGrp:fpacc.Mantissa[2], ax
		mov	ax, es:4[di]
		mov	StdGrp:fpacc.Mantissa[4], ax
		mov	ax, es:6[di]
		mov	StdGrp:fpacc.Mantissa[6], ax
;
; Adjust the return address to point past the floating point number we
; just loaded.
;
		add	word ptr 2[bp], 10
;
		pop	ax
		pop	di
		pop	es
		pop	bp
		ret
sl_LEFPAL	endp
;
;
; sl_SEFPA-	Stores FPACC into in the extended precision variable
;		pointed at by ES:DI.
;
		public	sl_sefpa
sl_SEFPA	proc	far
		assume	ds:stdgrp
		push	ds
		push	ax
		mov	ax, StdGrp
		mov	ds, ax
		push	fpacc.Mantissa [0]
		push	fpacc.Mantissa [2]
		push	fpacc.Mantissa [4]
		push	fpacc.Mantissa [6]
		push	fpacc.Exponent
;
		mov	ax, fpacc.Exponent
		sub	ax, 4000h
		cmp	ax, 8000h