*
*   $VER: tanh.s 33.1 (22.1.97)
*
*   hyperbolic tangent
*
*   Version history:
*
*   33.1    22.1.97 (c) Motorola
*
*           - snipped from M68060SP sources
*

    machine 68040
    fpu     1

    XDEF    _tanh
    XDEF    @tanh

    XREF    @exp
    XREF    @expm1

*************************************************************************
* tanh():  computes the hyperbolic tangent of a normalized input        *
*                                                                       *
* INPUT *************************************************************** *
*       fp0 = extended precision input                                  *
*                                                                       *
* OUTPUT ************************************************************** *
*       fp0 = tanh(X)                                                   *
*                                                                       *
* ACCURACY and MONOTONICITY ******************************************* *
*       The returned result is within 3 ulps in 64 significant bit,     *
*       i.e. within 0.5001 ulp to 53 bits if the result is subsequently *
*       rounded to double precision. The result is provably monotonic   *
*       in double precision.                                            *
*                                                                       *
* ALGORITHM *********************************************************** *
*                                                                       *
*       TANH                                                            *
*       1. If |X| >= (5/2) log2 or |X| <= 2**(-40), go to 3.            *
*                                                                       *
*       2. (2**(-40) < |X| < (5/2) log2) Calculate tanh(X) by           *
*               sgn := sign(X), y := 2|X|, z := expm1(Y), and           *
*               tanh(X) = sgn*( z/(2+z) ).                              *
*               Exit.                                                   *
*                                                                       *
*       3. (|X| <= 2**(-40) or |X| >= (5/2) log2). If |X| < 1,          *
*               go to 7.                                                *
*                                                                       *
*       4. (|X| >= (5/2) log2) If |X| >= 50 log2, go to 6.              *
*                                                                       *
*       5. ((5/2) log2 <= |X| < 50 log2) Calculate tanh(X) by           *
*               sgn := sign(X), y := 2|X|, z := exp(Y),                 *
*               tanh(X) = sgn - [ sgn*2/(1+z) ].                        *
*               Exit.                                                   *
*                                                                       *
*       6. (|X| >= 50 log2) Tanh(X) = +-1 (round to nearest). Thus, we  *
*               calculate Tanh(X) by                                    *
*               sgn := sign(X), Tiny := 2**(-126),                      *
*               tanh(X) := sgn - sgn*Tiny.                              *
*               Exit.                                                   *
*                                                                       *
*       7. (|X| < 2**(-40)). Tanh(X) = X.       Exit.                   *
*                                                                       *
*************************************************************************


X       EQU             -12
XFRAC   EQU             X+4
V       EQU             -12
SGN     EQU             -16

TEMP_SIZE EQU           16

_tanh
        fmove.d         (4,sp),fp0
@tanh
        link            a0,#-TEMP_SIZE
        fmove.x         fp0,(X,a0)
        move.l          (X,a0),d1
        move.w          (XFRAC,a0),d1
        move.l          d1,(X,a0)
        and.l           #$7FFFFFFF,d1
        cmp.l           #$3fd78000,d1           ; is |X| < 2^(-40)?
        blt.w           .TANHBORS               ; yes
        cmp.l           #$3fffddce,d1           ; is |X| > (5/2)LOG2?
        bgt.w           .TANHBORS               ; yes

;--THIS IS THE USUAL CASE
;--Y = 2|X|, Z = EXPM1(Y), TANH(X) = SIGN(X) * Z / (Z+2).

        move.l          (X,a0),d1
        move.l          d1,(SGN,a0)
        and.l           #$7FFF0000,d1
        add.l           #$00010000,d1           ; EXPONENT OF 2|X|
        move.l          d1,(X,a0)
        and.l           #$80000000,(SGN,a0)
        fmove.x         (X,a0),fp0              ; FP0 IS Y = 2|X|

        jsr             @expm1                  ; FP0 IS Z = EXPM1(Y)

        fmove.x         fp0,fp1
        fadd.s          #$40000000,fp1          ; Z+2
        move.l          (SGN,a0),d1
        fmove.x         fp1,(V,a0)
        eor.l           d1,(V,a0)

        fdiv.x          (V,a0),fp0
        unlk            a0
        rts

.TANHBORS
        cmp.l           #$3FFF8000,d1
        blt.w           .TANHSM

        cmp.l           #$40048AA1,d1
        bgt.w           .TANHHUGE

;-- (5/2) LOG2 < |X| < 50 LOG2,
;--TANH(X) = 1 - (2/[EXP(2X)+1]). LET Y = 2|X|, SGN = SIGN(X),
;--TANH(X) = SGN -      SGN*2/[EXP(Y)+1].

        move.l          (X,a0),d1
        move.l          d1,(SGN,a0)
        and.l           #$7FFF0000,d1
        add.l           #$00010000,d1           ; EXPO OF 2|X|
        move.l          d1,(X,a0)               ; Y = 2|X|
        and.l           #$80000000,(SGN,a0)
        move.l          (SGN,a0),d1
        fmove.x         (X,a0),fp0              ; Y = 2|X|

        jsr             @exp                    ; FP0 IS EXP(Y)

        move.l          (SGN,a0),d1
        fadd.s          #$3F800000,fp0          ; EXP(Y)+1

        eor.l           #$C0000000,d1           ; -SIGN(X)*2
        fmove.s         d1,fp1                  ; -SIGN(X)*2 IN SGL FMT
        fdiv.x          fp0,fp1                 ; -SIGN(X)2 / [EXP(Y)+1 ]

        move.l          (SGN,a0),d1
        or.l            #$3F800000,d1           ; SGN
        fmove.s         d1,fp0                  ; SGN IN SGL FMT

        fadd.x          fp1,fp0
.TANHSM
        unlk            a0
        rts

;---RETURN SGN(X) - SGN(X)EPS
.TANHHUGE
        move.l          (X,a0),d1
        and.l           #$80000000,d1
        or.l            #$3F800000,d1
        fmove.s         d1,fp0
        and.l           #$80000000,d1
        eor.l           #$80800000,d1           ; -SIGN(X)*EPS

        fadd.s          d1,fp0
        unlk            a0
        rts

