l_setox.s

VXWORKS源代码

	.long	0x3FFF0000,0xDE60F482,0x5E0E9124,0x9E8BE175
	.long	0x3FFF0000,0xE0CCDEEC,0x2A94E110,0x20032C4B
	.long	0x3FFF0000,0xE33F8972,0xBE8A5A50,0x2004DFF5
	.long	0x3FFF0000,0xE5B906E7,0x7C8348A8,0x1E72F47A
	.long	0x3FFF0000,0xE8396A50,0x3C4BDC68,0x1F722F22
	.long	0x3FFF0000,0xEAC0C6E7,0xDD243930,0xA017E945
	.long	0x3FFF0000,0xED4F301E,0xD9942B84,0x1F401A5B
	.long	0x3FFF0000,0xEFE4B99B,0xDCDAF5CC,0x9FB9A9E3
	.long	0x3FFF0000,0xF281773C,0x59FFB138,0x20744C05
	.long	0x3FFF0000,0xF5257D15,0x2486CC2C,0x1F773A19
	.long	0x3FFF0000,0xF7D0DF73,0x0AD13BB8,0x1FFE90D5
	.long	0x3FFF0000,0xFA83B2DB,0x722A033C,0xA041ED22
	.long	0x3FFF0000,0xFD3E0C0C,0xF486C174,0x1F853F3A

#define	ADJFLAG	 L_SCR2
#define	SCALE	 FP_SCR1
#define	ADJSCALE  FP_SCR2
#define	SC	 FP_SCR3
#define	ONEBYSC	 FP_SCR4

|	xref	__l_t_frcinx
|	xref	__l_t_extdnrm
|	xref	__l_t_unfl
|	xref	__l_t_ovfl


	.text

	.globl	__l_setoxd
__l_setoxd:
|--entry point for EXP(X), X is denormalized
	movel		a0@,d0
	andil		#0x80000000,d0
	oril		#0x00800000,d0		|...sign(X)*2^(-126)
	movel		d0,a7@-
/*	fmoves	&0x3F800000,fp0 */	 .long 0xf23c4400,0x3f800000
	fmovel		d1,fpcr
	fadds		a7@+,fp0
	jra 		__l_t_frcinx

	.globl	__l_setox
__l_setox:
/* |--entry point for EXP(X), here X is finite, non-zero, and not NaN's */

|--Step 1.
	movel		a0@,d0	 |...load part of input X
	andil		#0x7FFF0000,d0	|...biased expo. of X
	cmpil		#0x3FBE0000,d0	|...2^(-65)
	jge 		EXPC1		|...normal case
	jra 		EXPSM

EXPC1:
|--The case |X| >= 2^(-65)
	movew		a0@(4),d0	|...expo. and partial sig. of |X|
	cmpil		#0x400CB167,d0	|...16380 log2 trunc. 16 bits
	jlt 		EXPMAIN	 |...normal case
	jra 		EXPBIG

EXPMAIN:
|--Step 2.
|--This is the normal branch:	2^(-65) <= |X| < 16380 log2.
	fmovex		a0@,fp0	|...load input from	a0@

	fmovex		fp0,fp1
/*	fmuls	&0x42B8AA3B,fp0 */	 .long 0xf23c4423,0x42b8aa3b
	fmovemx	fp2/fp3,a7@-		|...save fp2
	movel		#0,a6@(ADJFLAG)
	fmovel		fp0,d0		|...N = int( X * 64/log2 )
	lea		EXPTBL,a1
	fmovel		d0,fp0		|...convert to floating-format

	movel		d0,a6@(L_SCR1)	|...save N temporarily
	andil		#0x3F,d0		|...D0 is J = N mod 64
	lsll		#4,d0
	addal		d0,a1		|...address of 2^(J/64)
	movel		a6@(L_SCR1),d0
	asrl		#6,d0		|...D0 is M
	addiw		#0x3FFF,d0	|...biased expo. of 2^(M)
	movew		L2,a6@(L_SCR1)	|...prefetch L2, no need in CB

EXPCONT1:
|--Step 3.
|--fp1,fp2 saved on the stack. fp0 is N, fp1 is X,
|--a0 points to 2^(J/64), D0 is biased expo. of 2^(M)
	fmovex		fp0,fp2
/*	fmuls	&0xBC317218,fp0 */	 .long 0xf23c4423,0xbc317218
	fmulx		L2,fp2		|...N * L2, L1+L2 = -log2/64
	faddx		fp1,fp0	 	|...X + N*L1
	faddx		fp2,fp0		|...fp0 is R, reduced arg.
|	MOVE.w		#0x3FA5,EXPA3	|...load EXPA3 in cache

|--Step 4.
|--WE NOW COMPUTE EXP(R)-1 BY A POLYNOMIAL
|-- R + R*R*(A1 + R*(A2 + R*(A3 + R*(A4 + R*A5))))
|--TO FULLY UTILIZE THE PIPELINE, WE COMPUTE S = R*R
|--[R+R*S*(A2+S*A4)] + [S*(A1+S*(A3+S*A5))]

	fmovex		fp0,fp1
	fmulx		fp1,fp1	 	|...fp1 IS S = R*R

/*	fmoves	&0x3AB60B70,fp2 */	 .long 0xf23c4500,0x3ab60b70
|	MOVE.w		#0,a1@(2)	|...load 2^(J/64) in cache

	fmulx		fp1,fp2	 	|...fp2 IS S*A5
	fmovex		fp1,fp3
/*	fmuls	&0x3C088895,fp3 */	 .long 0xf23c45a3,0x3c088895

	faddd		EXPA3,fp2	|...fp2 IS a3+S*A5
	faddd		EXPA2,fp3	|...fp3 IS a2+S*A4

	fmulx		fp1,fp2	 	|...fp2 IS S*(A3+S*A5)
	movew		d0,a6@(SCALE)	|...SCALE is 2^(M) in extended
	clrw		a6@(SCALE+2)
	movel		#0x80000000,a6@(SCALE+4)
	clrl		a6@(SCALE+8)

	fmulx		fp1,fp3	 	|...fp3 IS S*(A2+S*A4)

/*	fadds	&0x3F000000,fp2 */	 .long 0xf23c4522,0x3f000000
	fmulx		fp0,fp3	 	|...fp3 IS R*S*(A2+S*A4)

	fmulx		fp1,fp2	 	|...fp2 IS S*(A1+S*(A3+S*A5))
	faddx		fp3,fp0	 	|...fp0 IS R+R*S*(A2+S*A4),
|					|...fp3 released

	fmovex		a1@+,fp1	|...fp1 is lead. pt. of 2^(J/64)
	faddx		fp2,fp0	 	|...fp0 is EXP(R) - 1
|					|...fp2 released

|--Step 5
|--final reconstruction process
|--EXP(X) = 2^M * ( 2^(J/64) + 2^(J/64)*(EXP(R)-1) )

	fmulx		fp1,fp0	 	|...2^(J/64)*(Exp(R)-1)
	fmovemx	a7@+,fp2/fp3	|...fp2 restored
	fadds		a1@,fp0	|...accurate 2^(J/64)

	faddx		fp1,fp0	 	|...2^(J/64) + 2^(J/64)*...
	movel		a6@(ADJFLAG),d0

|--Step 6
	tstl		d0
	jeq 		NORMAL
ADJUST:
	fmulx		a6@(ADJSCALE),fp0
NORMAL:
	fmovel		d1,fpcr	 	|...restore user fpcr
	fmulx		a6@(SCALE),fp0	|...multiply 2^(M)
	jra 		__l_t_frcinx

EXPSM:
|--Step 7
	fmovemx	a0@,fp0-fp0	|...in case X is denormalized
	fmovel		d1,fpcr
/*	fadds	&0x3F800000,fp0 */	 .long 0xf23c4422,0x3f800000
	jra 		__l_t_frcinx

EXPBIG:
|--Step 8
	cmpil		#0x400CB27C,d0	|...16480 log2
	jgt 		EXP2BIG
|--Steps 8.2 -- 8.6
	fmovex		a0@,fp0	|...load input from	a0@

	fmovex		fp0,fp1
/*	fmuls	&0x42B8AA3B,fp0 */	 .long 0xf23c4423,0x42b8aa3b
	fmovemx	 fp2/fp3,a7@-		|...save fp2
	movel		#1,a6@(ADJFLAG)
	fmovel		fp0,d0		|...N = int( X * 64/log2 )
	lea		EXPTBL,a1
	fmovel		d0,fp0		|...convert to floating-format
	movel		d0,a6@(L_SCR1)			|...save N temporarily
	andil		#0x3F,d0		 |...D0 is J = N mod 64
	lsll		#4,d0
	addal		d0,a1			|...address of 2^(J/64)
	movel		a6@(L_SCR1),d0
	asrl		#6,d0			|...D0 is K
	movel		d0,a6@(L_SCR1)			|...save K temporarily
	asrl		#1,d0			|...D0 is M1
	subl		d0,a6@(L_SCR1)			|...a1 is M
	addiw		#0x3FFF,d0		|...biased expo. of 2^(M1)
	movew		d0,a6@(ADJSCALE)		|...ADJSCALE := 2^(M1)
	clrw		a6@(ADJSCALE+2)
	movel		#0x80000000,a6@(ADJSCALE+4)
	clrl		a6@(ADJSCALE+8)
	movel		a6@(L_SCR1),d0			|...D0 is M
	addiw		#0x3FFF,d0		|...biased expo. of 2^(M)
	jra 		EXPCONT1		|...go back to Step 3

EXP2BIG:
|--Step 9
	fmovel		d1,fpcr
	movel		a0@,d0
	bclr		#sign_bit,a0@		|...setox always returns positive
	cmpil		#0,d0
	jlt 		__l_t_unfl
	jra 		__l_t_ovfl

	.globl	__l_setoxm1d
__l_setoxm1d:
|--entry point for EXPM1(X), here X is denormalized
|--Step 0.
	jra 		__l_t_extdnrm


	.globl	__l_setoxm1
__l_setoxm1:
|--entry point for EXPM1(X), here X is finite, non-zero, non-NaN

|--Step 1.
|--Step 1.1
	movel		a0@,d0	 |...load part of input X
	andil		#0x7FFF0000,d0	|...biased expo. of X
	cmpil		#0x3FFD0000,d0	|...1/4
	jge 		EM1CON1	 |...|X| >= 1/4
	jra 		EM1SM

EM1CON1:
|--Step 1.3
|--The case |X| >= 1/4
	movew		a0@(4),d0	|...expo. and partial sig. of |X|
	cmpil		#0x4004C215,d0	|...70log2 rounded up to 16 bits
	jle 		EM1MAIN	 |...1/4 <= |X| <= 70log2
	jra 		EM1BIG

EM1MAIN:
|--Step 2.
|--This is the case:	1/4 <= |X| <= 70 log2.
	fmovex		a0@,fp0	|...load input from	a0@

	fmovex		fp0,fp1
/*	fmuls	&0x42B8AA3B,fp0 */	 .long 0xf23c4423,0x42b8aa3b
	fmovemx	fp2/fp3,a7@-		|...save fp2
|	MOVE.w		#0x3F81,EM1A4		|...prefetch in CB mode
	fmovel		fp0,d0		|...N = int( X * 64/log2 )
	lea		EXPTBL,a1
	fmovel		d0,fp0		|...convert to floating-format

	movel		d0,a6@(L_SCR1)			|...save N temporarily
	andil		#0x3F,d0		 |...D0 is J = N mod 64
	lsll		#4,d0
	addal		d0,a1			|...address of 2^(J/64)
	movel		a6@(L_SCR1),d0
	asrl		#6,d0			|...D0 is M
	movel		d0,a6@(L_SCR1)			|...save a copy of M
|	MOVE.w		#0x3FDC,L2		|...prefetch L2 in CB mode

|--Step 3.
|--fp1,fp2 saved on the stack. fp0 is N, fp1 is X,
|--a0 points to 2^(J/64), D0 and a1 both contain M
	fmovex		fp0,fp2
/*	fmuls	&0xBC317218,fp0 */	 .long 0xf23c4423,0xbc317218
	fmulx		L2,fp2		|...N * L2, L1+L2 = -log2/64
	faddx		fp1,fp0	 |...X + N*L1
	faddx		fp2,fp0	 |...fp0 is R, reduced arg.
|	MOVE.w		#0x3FC5,EM1A2		|...load EM1A2 in cache
	addiw		#0x3FFF,d0		|...D0 is biased expo. of 2^M

|--Step 4.
|--WE NOW COMPUTE EXP(R)-1 BY A POLYNOMIAL
|-- R + R*R*(A1 + R*(A2 + R*(A3 + R*(A4 + R*(A5 + R*A6)))))
|--TO FULLY UTILIZE THE PIPELINE, WE COMPUTE S = R*R
|--[R*S*(A2+S*(A4+S*A6))] + [R+S*(A1+S*(A3+S*A5))]

	fmovex		fp0,fp1
	fmulx		fp1,fp1		|...fp1 IS S = R*R

/*	fmoves	&0x3950097B,fp2 */	 .long 0xf23c4500,0x3950097b
|	MOVE.w		#0,a1@(2)	|...load 2^(J/64) in cache

	fmulx		fp1,fp2		|...fp2 IS S*A6
	fmovex		fp1,fp3
/*	fmuls	&0x3AB60B6A,fp3 */	 .long 0xf23c45a3,0x3ab60b6a

	faddd		EM1A4,fp2	|...fp2 IS a4+S*A6
	faddd		EM1A3,fp3	|...fp3 IS a3+S*A5
	movew		d0,a6@(SC)		|...SC is 2^(M) in extended
	clrw		a6@(SC+2)
	movel		#0x80000000,a6@(SC+4)
	clrl		a6@(SC+8)

	fmulx		fp1,fp2		|...fp2 IS S*(A4+S*A6)
	movel		a6@(L_SCR1),d0		|...D0 is	M
	negw		d0		|...D0 is -M
	fmulx		fp1,fp3		|...fp3 IS S*(A3+S*A5)
	addiw		#0x3FFF,d0	|...biased expo. of 2^(-M)
	faddd		EM1A2,fp2	|...fp2 IS a2+S*(A4+S*A6)
/*	fadds	&0x3F000000,fp3 */	 .long 0xf23c45a2,0x3f000000

	fmulx		fp1,fp2		|...fp2 IS S*(A2+S*(A4+S*A6))
	oriw		#0x8000,d0	|...signed/expo. of -2^(-M)
	movew		d0,a6@(ONEBYSC)	|...OnebySc is -2^(-M)
	clrw		a6@(ONEBYSC+2)
	movel		#0x80000000,a6@(ONEBYSC+4)
	clrl		a6@(ONEBYSC+8)
	fmulx		fp3,fp1		|...fp1 IS S*(A1+S*(A3+S*A5))
|					|...fp3 released

	fmulx		fp0,fp2		|...fp2 IS R*S*(A2+S*(A4+S*A6))
	faddx		fp1,fp0		|...fp0 IS R+S*(A1+S*(A3+S*A5))
|					|...fp1 released

	faddx		fp2,fp0		|...fp0 IS EXP(R)-1
|					|...fp2 released
	fmovemx	a7@+,fp2/fp3	|...fp2 restored

|--Step 5
|--Compute 2^(J/64)*p

	fmulx		a1@,fp0	|...2^(J/64)*(Exp(R)-1)

|--Step 6
|--Step 6.1
	movel		a6@(L_SCR1),d0		|...retrieve M
	cmpil		#63,d0
	jle 		MLE63
|--Step 6.2	M >= 64
	fmoves		a1@(12),fp1	|...fp1 is t
	faddx		a6@(ONEBYSC),fp1	|...fp1 is t+OnebySc
	faddx		fp1,fp0		|...p+(t+OnebySc), fp1 released
	faddx		a1@,fp0	|...T+(p+(t+OnebySc))
	jra 		EM1SCALE
MLE63:
|--Step 6.3	M <= 63
	cmpil		#-3,d0
	jge 		MGEN3
MLTN3:
|--Step 6.4	M <= -4
	fadds		a1@(12),fp0	|...p+t
	faddx		a1@,fp0	|...T+(p+t)
	faddx		a6@(ONEBYSC),fp0	|...OnebySc + (T+(p+t))
	jra 		EM1SCALE
MGEN3:
|--Step 6.5	-3 <= M <= 63
	fmovex		a1@+,fp1	|...fp1 is T
	fadds		a1@,fp0	|...fp0 is p+t
	faddx		a6@(ONEBYSC),fp1	|...fp1 is T+OnebySc
	faddx		fp1,fp0		|...(T+OnebySc)+(p+t)

EM1SCALE:
|--Step 6.6
	fmovel		d1,fpcr
	fmulx		a6@(SC),fp0

	jra 		__l_t_frcinx

EM1SM:
|--Step 7	|X| < 1/4.
	cmpil		#0x3FBE0000,d0	|...2^(-65)
	jge 		EM1POLY

EM1TINY:
|--Step 8	|X| < 2^(-65)
	cmpil		#0x00330000,d0	|...2^(-16312)
	jlt 		EM12TINY
|--Step 8.2
	movel		#0x80010000,a6@(SC)	|...SC is -2^(-16382)
	movel		#0x80000000,a6@(SC+4)
	clrl		a6@(SC+8)
	fmovex		a0@,fp0
	fmovel		d1,fpcr
	faddx		a6@(SC),fp0

	jra 		__l_t_frcinx

EM12TINY:
|--Step 8.3
	fmovex		a0@,fp0
	fmuld		TWO140,fp0
	movel		#0x80010000,a6@(SC)
	movel		#0x80000000,a6@(SC+4)
	clrl		a6@(SC+8)
	faddx		a6@(SC),fp0
	fmovel		d1,fpcr
	fmuld		TWON140,fp0

	jra 		__l_t_frcinx

EM1POLY:
|--Step 9	exp(X)-1 by a simple polynomial
	fmovex		a0@,fp0	|...fp0 is X
	fmulx		fp0,fp0		|...fp0 is S := X*X
	fmovemx	fp2/fp3,a7@-	|...save fp2
/*	fmoves	&0x2F30CAA8,fp1 */	 .long 0xf23c4480,0x2f30caa8
	fmulx		fp0,fp1		|...fp1 is S*B12
/*	fmoves	&0x310F8290,fp2 */	 .long 0xf23c4500,0x310f8290
/*	fadds	&0x32D73220,fp1 */	 .long 0xf23c44a2,0x32d73220

	fmulx		fp0,fp2		|...fp2 is S*B11
	fmulx		fp0,fp1		|...fp1 is S*(B10 + ...

/*	fadds	&0x3493F281,fp2 */	 .long 0xf23c4522,0x3493f281
	faddd		EM1B8,fp1	|...fp1 is B8+S*...

	fmulx		fp0,fp2		|...fp2 is S*(B9+...
	fmulx		fp0,fp1		|...fp1 is S*(B8+...

	faddd		EM1B7,fp2	|...fp2 is B7+S*...
	faddd		EM1B6,fp1	|...fp1 is B6+S*...

	fmulx		fp0,fp2		|...fp2 is S*(B7+...
	fmulx		fp0,fp1		|...fp1 is S*(B6+...

	faddd		EM1B5,fp2	|...fp2 is B5+S*...
	faddd		EM1B4,fp1	|...fp1 is B4+S*...

	fmulx		fp0,fp2		|...fp2 is S*(B5+...
	fmulx		fp0,fp1		|...fp1 is S*(B4+...

	faddd		EM1B3,fp2	|...fp2 is B3+S*...
	faddx		EM1B2,fp1	|...fp1 is B2+S*...

	fmulx		fp0,fp2		|...fp2 is S*(B3+...
	fmulx		fp0,fp1		|...fp1 is S*(B2+...

	fmulx		fp0,fp2		|...fp2 is S*S*(B3+...)
	fmulx		a0@,fp1	|...fp1 is X*S*(B2...

/*	fmuls	&0x3F000000,fp0 */	 .long 0xf23c4423,0x3f000000
	faddx		fp2,fp1		|...fp1 is Q
|					|...fp2 released

	fmovemx	a7@+,fp2/fp3	|...fp2 restored

	faddx		fp1,fp0		|...fp0 is S*B1+Q
|					|...fp1 released

	fmovel		d1,fpcr
	faddx		a0@,fp0

	jra 		__l_t_frcinx

EM1BIG:
|--Step 10	|X| > 70 log2
	movel		a0@,d0
	cmpil		#0,d0
	jgt 		EXPC1
|--Step 10.2
/*	fmoves	&0xBF800000,fp0 */	 .long 0xf23c4400,0xbf800000
	fmovel		d1,fpcr
/*	fadds	&0x00800000,fp0 */	 .long 0xf23c4422,0x00800000

	jra 		__l_t_frcinx

|	end