mathhardalib.s

vxworks源码源码解读是学习vxworks的最佳途径

/* mathHardALib.s - C-callable math routines for the i80387 */

/* Copyright 1984-1991 Wind River Systems, Inc. */
	.data
	.globl	_copyright_wind_river
	.long	_copyright_wind_river

/*
modification history 
--------------------
01e,29sep97,hdn  fixed a bug in _mathHardPow.
01d,17jun96,hdn  fixed a bug in _mathHardSincos.
01c,16jun93,hdn  updated to 5.1.
01b,14oct92,hdn  aligned all functions.
01a,16sep92,hdn  written by modifying Tron's mathHardALib.s.
*/

/*
DESCRIPTION
This library provides a C interface to the high-level math functions
on the i80387 floating-point coprocessor.  All angle-related
parameters and return values are expressed in radians.  Functions
capable errors, will set errno upon an error. All functions
included in this library whos names correspond to the ANSI C specification
are, indeed, ANSI-compatable. In the spirit of ANSI, HUGE_VAL is now
supported.

WARNING
This library works only if an i80387 coprocessor is in the system! 
Attempts to use these routines with no coprocessor present 
will result in illegal instruction traps.

SEE ALSO:
fppLib (1), floatLib (1), The C Programming Language - Second Edition

INCLUDE FILE: math.h

INTERNAL
Each routine has the following format:
    o calculate floating-point function using double parameter 
    o store result to %st0 register
*/

#define _ASMLANGUAGE
#include "vxWorks.h"
#include "asm.h"
#include "errno.h"

/* rounding mode in Control Reg. */
#define		FPP_ROUND_NEAREST	0x0000
#define		FPP_ROUND_ZERO		0x0c00
#define		FPP_ROUND_UP		0x0800
#define		FPP_ROUND_DOWN		0x0400
#define		FPP_ROUND_MASK		0xf3ff


	/* internal */

        .globl  _mathHardLog2
        .globl  _mathHardLog10
        .globl  _mathHardLog
        .globl  _mathHardExp
	.globl  _mathHardAsin
	.globl  _mathHardAcos
	.globl  _mathHardAtan
	.globl  _mathHardAtan2
        .globl  _mathHardTan
        .globl  _mathHardCos
        .globl  _mathHardSin
        .globl  _mathHardPow
        .globl  _mathHardSqrt
        .globl  _mathHardFabs
	.globl  _mathHardFmod

	.globl  _mathHardSincos
	.globl  _mathHardFloor
	.globl  _mathHardCeil
	.globl  _mathHardTrunc
	.globl  _mathHardRound
	.globl  _mathHardIround
	.globl  _mathHardIrint
	.globl  _mathHardInfinity

	/* external */

	.globl	_errnoSet

        .text
	.align 4


/*******************************************************************************
*
* mathHardAcos - ANSI-compatable hardware floating-point arc-cosine
*
* RETURNS: The arc-cosine in the range -pi/2 to pi/2 radians.
*
* double mathHardAcos (dblParam)
*     double dblParam;	/* angle in radians *
*
*/

_mathHardAcos:			/* acos(x) = pi/2 - atan(x/sqrt(1-x**2)) */
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fldl	DARG1(%ebp)
	fmul	%st,%st(1)
	fxch	%st(1)
	fld1
	fsubp	%st,%st(1)
	fsqrt
	fpatan
	fld1
	fld1
	faddp	%st,%st(1)
	fldpi
	fdivp	%st,%st(1)
	fsubp	%st,%st(1)

	leave
	ret

/*******************************************************************************
*
* mathHardAsin - ANSI-compatable hardware floating-point arc-sine
*
* RETURNS: The arc-sine in the range 0.0 to pi radians.
* 
* SEE ALSO: 
* floatLib (1), "The C Programming Language - Second Edition"
*
* double mathHardAsin (dblParam)
*     double dblParam;	/* angle in radians *
*
*/

	.align 4,0x90
_mathHardAsin:			/* asin(x) = atan(x/sqrt(1-x**2)) */
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fldl	DARG1(%ebp)
	fmul	%st,%st(1)
	fxch	%st(1)
	fld1
	fsubp	%st,%st(1)
	fsqrt
	fpatan

	leave
	ret

/*******************************************************************************
*
* mathHardAtan - ANSI-compatable hardware floating-point arc-tangent
*
* RETURNS: The arc-tangent of dblParam in the range -pi/2 to pi/2.
*
* SEE ALSO: floatLib (1), acos (2), asin (2)
*
* double mathHardAtan (dblParam)
*     double dblParam;	/* angle in radians *
*
*/

	.align 4,0x90
_mathHardAtan:
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fld1
	fpatan

	leave
	ret

/*******************************************************************************
*
* mathHardAtan2 - hardware floating point function for arctangent of (dblY/dblX)
*
* RETURNS:
*    The arc-tangent of (dblY/dblX) in the range -pi to pi.
*
* SEE ALSO: 
* floatLib (1), "The C Programming Language - Second Edition"
*
* double mathHardAtan2 (dblY, dblX)
*     double dblY;		/* Y *
*     double dblX;		/* X *
*
*/

	.align 4,0x90
_mathHardAtan2:
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fldl	DARG2(%ebp)
	fpatan

	leave
	ret
    
/*******************************************************************************
*
* mathHardCos - ANSI-compatable hardware floating-point cosine
*
* RETURNS: the cosine of the radian argument dblParam
*
* SEE ALSO: 
* floatLib (1), sin (2), cos (2), tan(2),
* "The C Programming Language - Second Edition"
*
* double mathHardCos (dblParam)
*     double dblParam;	/* angle in radians *
*
*/

	.align 4,0x90
_mathHardCos:
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fcos

	fstl	DARG1(%ebp)
	fwait
	movl	DARG1+4(%ebp),%eax
	movl	DARG1(%ebp),%edx

	leave
	ret

/*******************************************************************************
*
* mathHardExp - hardware floating-point exponential function
*
* RETURNS:
*    Floating-point inverse natural logarithm (e ** (dblExponent)).
*
* SEE ALSO: 
* floatLib (1), "The C Programming Language - Second Edition"
*
* double mathHardExp (dblExponent)
*     double dblExponent;	/* argument *
*
*/

	.align 4,0x90
_powerOftwo:
	/* 2**z. 
	 * z1 is integer of z. z2 is fractal of z.
	 * if z2 is greater than 0.5, z2 -= 0.5, then
	 *     2**z = 2**(z1 + z2 + 0.5) = 2**(z1) * 2**(z2) * 2**(0.5)
	 * if z2 is less than 0.5, then
	 *     2**z = 2**(z1 + z2) = 2**(z1) * 2**(z2)
	 */
	pushl	%ebp
	movl	%esp,%ebp

	/* change the Round Control bits */
	subl	$8,%esp
	fstcw	-4(%ebp)
	fwait
	movw	-4(%ebp),%ax
	andw	$ FPP_ROUND_MASK,%ax
	orw	$ FPP_ROUND_DOWN,%ax
	movw	%ax,-8(%ebp)
	fldcw	-8(%ebp)

	/* get z1(integer part of z) and z2(fractal part of z) */
	fldl	%st
	frndint
	fldcw	-4(%ebp)
	fsub	%st,%st(1)
	fxch	%st(1)
	fchs

	/* get a value 0.5 */
	fld1
	fchs
	fld1
	fscale
	fxch	%st(1)
	fstp	%st

	/* get z2 % 0.5 */
	fxch	%st(1)
	fprem
	fstsw	%ax
	fstp	%st(1)

	/* get A = 2**(z2) */
	f2xm1
	fld1
	faddp	%st,%st(1)
	test	$0x0200,%ax
	jz	powerOftwo0

	/* get A = 2**(z2) * 2**(0.5) */
	fld1
	fadd	%st,%st(0)
	fsqrt
	fmulp	%st,%st(1)

powerOftwo0:
	/* get 2**z = A * 2**(z1) */
	fscale
	fstp	%st(1)
	addl	$8,%esp

	leave
	ret

	.align 4,0x90
_mathHardExp:				/* e**y = 2**(y * log2(e)) */
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fldl2e
	fmulp	%st,%st(1)
	call	_powerOftwo

	leave
	ret

/*******************************************************************************
*
* mathHardFabs - ANSI-compatable hardware floating-point absolute value
*
* RETURNS: The floating-point absolute value of dblParam.
*
* SEE ALSO: 
* floatLib (1), "The C Programming Language - Second Edition"
*
* double mathHardFabs (dblParam)
*     double dblParam;	/* argument *
*
*/

	.align 4,0x90
_mathHardFabs: 
	pushl	%ebp
	movl	%esp,%ebp

	fldl	DARG1(%ebp)
	fabs

	leave
	ret

/*******************************************************************************
*
* mathHardLog - ANSI-compatable hardware floating-point natural logarithm 
*
* RETURNS: The natural logarithm of dblParam.
*
* SEE ALSO: 
* floatLib (1), "The C Programming Language - Second Edition"
*
* double mathHardLog (dblParam)
*     double dblParam;	/* argument *
*
*/

	.align 4,0x90
_mathHardLog:			/* loge(x) = loge(2) * log2(x) */
	pushl	%ebp
	movl	%esp,%ebp

	fldln2			/* st0 = loge(2) */
	fldl	DARG1(%ebp)	/* st0 = x, st1 = loge(2) */
	fyl2x			/* st0 = loge(2) * log2(x) */

	leave
	ret

/*******************************************************************************
*
* mathHardLog10 - ANSI-compatable hardware floating-point base 10 logarithm
*
* RETURNS: The logarithm (base 10) of dblParam.
*
* SEE ALSO: floatLib (1), log2 (2)
*
* double mathHardLog10 (dblParam)
*     double dblParam;	/* argument *
*
*/

	.align 4,0x90
_mathHardLog10:			/* log10(x) = log10(2) * log2(x) */
	pushl	%ebp
	movl	%esp,%ebp

	fldlg2			/* st0 = log10(2) */
	fldl	DARG1(%ebp)	/* st0 = x, st1 = log10(2) */
	fyl2x			/* st0 = log10(2) * log2(x) */

	leave
	ret

/*******************************************************************************
*
* mathHardLog2 - ANSI-compatable hardware floating-point logarithm base 2 
*
* RETURNS: The logarithm (base 2) of dblParam.
*
* SEE ALSO: floatLib (1), log10 (2)
*
* double mathHardLog2 (dblParam)
*     double dblParam;	/* argument *
*
*/

	.align 4,0x90
_mathHardLog2:
	pushl	%ebp
	movl	%esp,%ebp

	fld1
	fldl	DARG1(%ebp)
	fyl2x

	leave
	ret

/*******************************************************************************
*
* mathHardPow - ANSI-compatable hardware floating-point power function
*
* RETURNS: The floating-point value of dblX to the power of dblY.
*
* SEE ALSO: floatLib (1), sqrt (2)
*
* double mathHardPow (dblX, dblY)
*     double dblX;	/* X *
*     double dblY;	/* Y *
*
*/

	.align 4,0x90
_mathHardPow: 				/* x**y = 2**(y * log2(x)) */
	pushl	%ebp
	movl	%esp,%ebp

	fldz
	fldl	DARG1(%ebp)
	fcompp
	fstsw	%ax
	sahf