Home | History | Annotate | Line # | Download | only in spmath
dbl_float.h revision 1.2.14.1
      1  1.2.14.1       snj /*	$NetBSD: dbl_float.h,v 1.2.14.1 2009/02/02 00:13:47 snj Exp $	*/
      2       1.1  fredette 
      3  1.2.14.1       snj /*	$OpenBSD: dbl_float.h,v 1.10 2004/01/02 14:39:01 mickey Exp $	*/
      4       1.1  fredette 
      5       1.1  fredette /*
      6       1.1  fredette  * Copyright 1996 1995 by Open Software Foundation, Inc.
      7       1.1  fredette  *              All Rights Reserved
      8       1.1  fredette  *
      9       1.1  fredette  * Permission to use, copy, modify, and distribute this software and
     10       1.1  fredette  * its documentation for any purpose and without fee is hereby granted,
     11       1.1  fredette  * provided that the above copyright notice appears in all copies and
     12       1.1  fredette  * that both the copyright notice and this permission notice appear in
     13       1.1  fredette  * supporting documentation.
     14       1.1  fredette  *
     15       1.1  fredette  * OSF DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE
     16       1.1  fredette  * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
     17       1.1  fredette  * FOR A PARTICULAR PURPOSE.
     18       1.1  fredette  *
     19       1.1  fredette  * IN NO EVENT SHALL OSF BE LIABLE FOR ANY SPECIAL, INDIRECT, OR
     20       1.1  fredette  * CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM
     21       1.1  fredette  * LOSS OF USE, DATA OR PROFITS, WHETHER IN ACTION OF CONTRACT,
     22       1.1  fredette  * NEGLIGENCE, OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION
     23       1.1  fredette  * WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
     24       1.1  fredette  */
     25       1.1  fredette /*
     26       1.1  fredette  * pmk1.1
     27       1.1  fredette  */
     28       1.1  fredette /*
     29       1.1  fredette  * (c) Copyright 1986 HEWLETT-PACKARD COMPANY
     30       1.1  fredette  *
     31       1.1  fredette  * To anyone who acknowledges that this file is provided "AS IS"
     32       1.1  fredette  * without any express or implied warranty:
     33       1.1  fredette  *     permission to use, copy, modify, and distribute this file
     34       1.1  fredette  * for any purpose is hereby granted without fee, provided that
     35       1.1  fredette  * the above copyright notice and this notice appears in all
     36       1.1  fredette  * copies, and that the name of Hewlett-Packard Company not be
     37       1.1  fredette  * used in advertising or publicity pertaining to distribution
     38       1.1  fredette  * of the software without specific, written prior permission.
     39       1.1  fredette  * Hewlett-Packard Company makes no representations about the
     40       1.1  fredette  * suitability of this software for any purpose.
     41       1.1  fredette  */
     42       1.1  fredette 
     43       1.1  fredette #include <sys/cdefs.h>
     44       1.1  fredette 
     45       1.1  fredette /**************************************
     46       1.1  fredette  * Declare double precision functions *
     47       1.1  fredette  **************************************/
     48       1.1  fredette 
     49       1.1  fredette /* 32-bit word grabing functions */
     50       1.1  fredette #define Dbl_firstword(value) Dallp1(value)
     51       1.1  fredette #define Dbl_secondword(value) Dallp2(value)
     52       1.1  fredette #define Dbl_thirdword(value) dummy_location
     53       1.1  fredette #define Dbl_fourthword(value) dummy_location
     54       1.1  fredette 
     55       1.1  fredette #define Dbl_sign(object) Dsign(object)
     56       1.1  fredette #define Dbl_exponent(object) Dexponent(object)
     57       1.1  fredette #define Dbl_signexponent(object) Dsignexponent(object)
     58       1.1  fredette #define Dbl_mantissap1(object) Dmantissap1(object)
     59       1.1  fredette #define Dbl_mantissap2(object) Dmantissap2(object)
     60       1.1  fredette #define Dbl_exponentmantissap1(object) Dexponentmantissap1(object)
     61       1.1  fredette #define Dbl_allp1(object) Dallp1(object)
     62       1.1  fredette #define Dbl_allp2(object) Dallp2(object)
     63       1.1  fredette 
     64       1.1  fredette /* dbl_and_signs ands the sign bits of each argument and puts the result
     65       1.1  fredette  * into the first argument. dbl_or_signs ors those same sign bits */
     66       1.1  fredette #define Dbl_and_signs( src1dst, src2)		\
     67       1.1  fredette     Dallp1(src1dst) = (Dallp1(src2)|~(1<<31)) & Dallp1(src1dst)
     68       1.1  fredette #define Dbl_or_signs( src1dst, src2)		\
     69       1.1  fredette     Dallp1(src1dst) = (Dallp1(src2)&(1<<31)) | Dallp1(src1dst)
     70       1.1  fredette 
     71       1.1  fredette /* The hidden bit is always the low bit of the exponent */
     72       1.1  fredette #define Dbl_clear_exponent_set_hidden(srcdst) Deposit_dexponent(srcdst,1)
     73       1.1  fredette #define Dbl_clear_signexponent_set_hidden(srcdst) \
     74       1.1  fredette     Deposit_dsignexponent(srcdst,1)
     75       1.1  fredette #define Dbl_clear_sign(srcdst) Dallp1(srcdst) &= ~(1<<31)
     76       1.1  fredette #define Dbl_clear_signexponent(srcdst) \
     77       1.1  fredette     Dallp1(srcdst) &= Dmantissap1((unsigned)-1)
     78       1.1  fredette 
     79       1.1  fredette /* Exponent field for doubles has already been cleared and may be
     80       1.1  fredette  * included in the shift.  Here we need to generate two double width
     81       1.1  fredette  * variable shifts.  The insignificant bits can be ignored.
     82       1.1  fredette  *      MTSAR f(varamount)
     83       1.1  fredette  *      VSHD	srcdst.high,srcdst.low => srcdst.low
     84       1.1  fredette  *	VSHD	0,srcdst.high => srcdst.high
     85       1.1  fredette  * This is very difficult to model with C expressions since the shift amount
     86       1.1  fredette  * could exceed 32.  */
     87       1.1  fredette /* varamount must be less than 64 */
     88       1.1  fredette #define Dbl_rightshift(srcdstA, srcdstB, varamount)			\
     89       1.1  fredette     {if((varamount) >= 32) {						\
     90       1.1  fredette 	Dallp2(srcdstB) = Dallp1(srcdstA) >> (varamount-32);		\
     91       1.1  fredette 	Dallp1(srcdstA)=0;						\
     92       1.1  fredette     }									\
     93       1.1  fredette     else if(varamount > 0) {						\
     94       1.1  fredette 	Variable_shift_double(Dallp1(srcdstA), Dallp2(srcdstB),		\
     95       1.1  fredette 	  (varamount), Dallp2(srcdstB));				\
     96       1.1  fredette 	Dallp1(srcdstA) >>= varamount;					\
     97       1.1  fredette     } }
     98       1.1  fredette /* varamount must be less than 64 */
     99       1.1  fredette #define Dbl_rightshift_exponentmantissa(srcdstA, srcdstB, varamount)	\
    100       1.1  fredette     {if((varamount) >= 32) {						\
    101       1.1  fredette 	Dallp2(srcdstB) = Dexponentmantissap1(srcdstA) >> ((varamount)-32); \
    102       1.1  fredette 	Dallp1(srcdstA) &= (1<<31);  /* clear exponentmantissa field */ \
    103       1.1  fredette     }									\
    104       1.1  fredette     else if(varamount > 0) {						\
    105       1.1  fredette 	Variable_shift_double(Dexponentmantissap1(srcdstA), Dallp2(srcdstB), \
    106       1.1  fredette 	(varamount), Dallp2(srcdstB));					\
    107       1.1  fredette 	Deposit_dexponentmantissap1(srcdstA,				\
    108       1.1  fredette 	    (Dexponentmantissap1(srcdstA)>>(varamount)));			\
    109       1.1  fredette     } }
    110       1.1  fredette /* varamount must be less than 64 */
    111       1.1  fredette #define Dbl_leftshift(srcdstA, srcdstB, varamount)			\
    112       1.1  fredette     {if((varamount) >= 32) {						\
    113       1.1  fredette 	Dallp1(srcdstA) = Dallp2(srcdstB) << (varamount-32);		\
    114       1.1  fredette 	Dallp2(srcdstB)=0;						\
    115       1.1  fredette     }									\
    116       1.1  fredette     else {								\
    117       1.1  fredette 	if ((varamount) > 0) {						\
    118       1.1  fredette 	    Dallp1(srcdstA) = (Dallp1(srcdstA) << (varamount)) |	\
    119       1.1  fredette 		(Dallp2(srcdstB) >> (32-(varamount)));			\
    120       1.1  fredette 	    Dallp2(srcdstB) <<= varamount;				\
    121       1.1  fredette 	}								\
    122       1.1  fredette     } }
    123       1.1  fredette #define Dbl_leftshiftby1_withextent(lefta,leftb,right,resulta,resultb)	\
    124       1.1  fredette     Shiftdouble(Dallp1(lefta), Dallp2(leftb), 31, Dallp1(resulta));	\
    125       1.1  fredette     Shiftdouble(Dallp2(leftb), Extall(right), 31, Dallp2(resultb))
    126       1.1  fredette 
    127       1.1  fredette #define Dbl_rightshiftby1_withextent(leftb,right,dst)		\
    128       1.1  fredette     Extall(dst) = (Dallp2(leftb) << 31) | ((unsigned)Extall(right) >> 1) | \
    129       1.1  fredette 		  Extlow(right)
    130       1.1  fredette 
    131       1.1  fredette #define Dbl_arithrightshiftby1(srcdstA,srcdstB)			\
    132       1.1  fredette     Shiftdouble(Dallp1(srcdstA),Dallp2(srcdstB),1,Dallp2(srcdstB));\
    133       1.1  fredette     Dallp1(srcdstA) = (int)Dallp1(srcdstA) >> 1
    134       1.1  fredette 
    135       1.1  fredette /* Sign extend the sign bit with an integer destination */
    136       1.1  fredette #define Dbl_signextendedsign(value)  Dsignedsign(value)
    137       1.1  fredette 
    138       1.1  fredette #define Dbl_isone_hidden(dbl_value) (Is_dhidden(dbl_value)!=0)
    139       1.1  fredette /* Singles and doubles may include the sign and exponent fields.  The
    140       1.1  fredette  * hidden bit and the hidden overflow must be included. */
    141       1.1  fredette #define Dbl_increment(dbl_valueA,dbl_valueB) \
    142       1.1  fredette     if( (Dallp2(dbl_valueB) += 1) == 0 )  Dallp1(dbl_valueA) += 1
    143       1.1  fredette #define Dbl_increment_mantissa(dbl_valueA,dbl_valueB) \
    144       1.1  fredette     if( (Dmantissap2(dbl_valueB) += 1) == 0 )  \
    145       1.1  fredette     Deposit_dmantissap1(dbl_valueA,dbl_valueA+1)
    146       1.1  fredette #define Dbl_decrement(dbl_valueA,dbl_valueB) \
    147       1.1  fredette     if( Dallp2(dbl_valueB) == 0 )  Dallp1(dbl_valueA) -= 1; \
    148       1.1  fredette     Dallp2(dbl_valueB) -= 1
    149       1.1  fredette 
    150       1.1  fredette #define Dbl_isone_sign(dbl_value) (Is_dsign(dbl_value)!=0)
    151       1.1  fredette #define Dbl_isone_hiddenoverflow(dbl_value) (Is_dhiddenoverflow(dbl_value)!=0)
    152       1.1  fredette #define Dbl_isone_lowmantissap1(dbl_valueA) (Is_dlowp1(dbl_valueA)!=0)
    153       1.1  fredette #define Dbl_isone_lowmantissap2(dbl_valueB) (Is_dlowp2(dbl_valueB)!=0)
    154       1.1  fredette #define Dbl_isone_signaling(dbl_value) (Is_dsignaling(dbl_value)!=0)
    155       1.1  fredette #define Dbl_is_signalingnan(dbl_value) (Dsignalingnan(dbl_value)==0xfff)
    156       1.1  fredette #define Dbl_isnotzero(dbl_valueA,dbl_valueB) \
    157       1.1  fredette     (Dallp1(dbl_valueA) || Dallp2(dbl_valueB))
    158       1.1  fredette #define Dbl_isnotzero_hiddenhigh7mantissa(dbl_value) \
    159       1.1  fredette     (Dhiddenhigh7mantissa(dbl_value)!=0)
    160       1.1  fredette #define Dbl_isnotzero_exponent(dbl_value) (Dexponent(dbl_value)!=0)
    161       1.1  fredette #define Dbl_isnotzero_mantissa(dbl_valueA,dbl_valueB) \
    162       1.1  fredette     (Dmantissap1(dbl_valueA) || Dmantissap2(dbl_valueB))
    163       1.1  fredette #define Dbl_isnotzero_mantissap1(dbl_valueA) (Dmantissap1(dbl_valueA)!=0)
    164       1.1  fredette #define Dbl_isnotzero_mantissap2(dbl_valueB) (Dmantissap2(dbl_valueB)!=0)
    165       1.1  fredette #define Dbl_isnotzero_exponentmantissa(dbl_valueA,dbl_valueB) \
    166       1.1  fredette     (Dexponentmantissap1(dbl_valueA) || Dmantissap2(dbl_valueB))
    167       1.1  fredette #define Dbl_isnotzero_low4p2(dbl_value) (Dlow4p2(dbl_value)!=0)
    168       1.1  fredette #define Dbl_iszero(dbl_valueA,dbl_valueB) (Dallp1(dbl_valueA)==0 && \
    169       1.1  fredette     Dallp2(dbl_valueB)==0)
    170       1.1  fredette #define Dbl_iszero_allp1(dbl_value) (Dallp1(dbl_value)==0)
    171       1.1  fredette #define Dbl_iszero_allp2(dbl_value) (Dallp2(dbl_value)==0)
    172       1.1  fredette #define Dbl_iszero_hidden(dbl_value) (Is_dhidden(dbl_value)==0)
    173       1.1  fredette #define Dbl_iszero_hiddenoverflow(dbl_value) (Is_dhiddenoverflow(dbl_value)==0)
    174       1.1  fredette #define Dbl_iszero_hiddenhigh3mantissa(dbl_value) \
    175       1.1  fredette     (Dhiddenhigh3mantissa(dbl_value)==0)
    176       1.1  fredette #define Dbl_iszero_hiddenhigh7mantissa(dbl_value) \
    177       1.1  fredette     (Dhiddenhigh7mantissa(dbl_value)==0)
    178       1.1  fredette #define Dbl_iszero_sign(dbl_value) (Is_dsign(dbl_value)==0)
    179       1.1  fredette #define Dbl_iszero_exponent(dbl_value) (Dexponent(dbl_value)==0)
    180       1.1  fredette #define Dbl_iszero_mantissa(dbl_valueA,dbl_valueB) \
    181       1.1  fredette     (Dmantissap1(dbl_valueA)==0 && Dmantissap2(dbl_valueB)==0)
    182       1.1  fredette #define Dbl_iszero_exponentmantissa(dbl_valueA,dbl_valueB) \
    183       1.1  fredette     (Dexponentmantissap1(dbl_valueA)==0 && Dmantissap2(dbl_valueB)==0)
    184       1.1  fredette #define Dbl_isinfinity_exponent(dbl_value)		\
    185       1.1  fredette     (Dexponent(dbl_value)==DBL_INFINITY_EXPONENT)
    186       1.1  fredette #define Dbl_isnotinfinity_exponent(dbl_value)		\
    187       1.1  fredette     (Dexponent(dbl_value)!=DBL_INFINITY_EXPONENT)
    188       1.1  fredette #define Dbl_isinfinity(dbl_valueA,dbl_valueB)			\
    189       1.1  fredette     (Dexponent(dbl_valueA)==DBL_INFINITY_EXPONENT &&	\
    190       1.1  fredette     Dmantissap1(dbl_valueA)==0 && Dmantissap2(dbl_valueB)==0)
    191       1.1  fredette #define Dbl_isnan(dbl_valueA,dbl_valueB)		\
    192       1.1  fredette     (Dexponent(dbl_valueA)==DBL_INFINITY_EXPONENT &&	\
    193       1.1  fredette     (Dmantissap1(dbl_valueA)!=0 || Dmantissap2(dbl_valueB)!=0))
    194       1.1  fredette #define Dbl_isnotnan(dbl_valueA,dbl_valueB)		\
    195       1.1  fredette     (Dexponent(dbl_valueA)!=DBL_INFINITY_EXPONENT ||	\
    196       1.1  fredette     (Dmantissap1(dbl_valueA)==0 && Dmantissap2(dbl_valueB)==0))
    197       1.1  fredette 
    198       1.1  fredette #define Dbl_islessthan(dbl_op1a,dbl_op1b,dbl_op2a,dbl_op2b)	\
    199       1.1  fredette     (Dallp1(dbl_op1a) < Dallp1(dbl_op2a) ||			\
    200       1.1  fredette      (Dallp1(dbl_op1a) == Dallp1(dbl_op2a) &&			\
    201       1.1  fredette       Dallp2(dbl_op1b) < Dallp2(dbl_op2b)))
    202       1.1  fredette #define Dbl_isgreaterthan(dbl_op1a,dbl_op1b,dbl_op2a,dbl_op2b)	\
    203       1.1  fredette     (Dallp1(dbl_op1a) > Dallp1(dbl_op2a) ||			\
    204       1.1  fredette      (Dallp1(dbl_op1a) == Dallp1(dbl_op2a) &&			\
    205       1.1  fredette       Dallp2(dbl_op1b) > Dallp2(dbl_op2b)))
    206       1.1  fredette #define Dbl_isnotlessthan(dbl_op1a,dbl_op1b,dbl_op2a,dbl_op2b)	\
    207       1.1  fredette     (Dallp1(dbl_op1a) > Dallp1(dbl_op2a) ||			\
    208       1.1  fredette      (Dallp1(dbl_op1a) == Dallp1(dbl_op2a) &&			\
    209       1.1  fredette       Dallp2(dbl_op1b) >= Dallp2(dbl_op2b)))
    210       1.1  fredette #define Dbl_isnotgreaterthan(dbl_op1a,dbl_op1b,dbl_op2a,dbl_op2b) \
    211       1.1  fredette     (Dallp1(dbl_op1a) < Dallp1(dbl_op2a) ||			\
    212       1.1  fredette      (Dallp1(dbl_op1a) == Dallp1(dbl_op2a) &&			\
    213       1.1  fredette       Dallp2(dbl_op1b) <= Dallp2(dbl_op2b)))
    214       1.1  fredette #define Dbl_isequal(dbl_op1a,dbl_op1b,dbl_op2a,dbl_op2b)	\
    215       1.1  fredette      ((Dallp1(dbl_op1a) == Dallp1(dbl_op2a)) &&			\
    216       1.1  fredette       (Dallp2(dbl_op1b) == Dallp2(dbl_op2b)))
    217       1.1  fredette 
    218       1.1  fredette #define Dbl_leftshiftby8(dbl_valueA,dbl_valueB) \
    219       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),24,Dallp1(dbl_valueA)); \
    220       1.1  fredette     Dallp2(dbl_valueB) <<= 8
    221       1.1  fredette #define Dbl_leftshiftby7(dbl_valueA,dbl_valueB) \
    222       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),25,Dallp1(dbl_valueA)); \
    223       1.1  fredette     Dallp2(dbl_valueB) <<= 7
    224       1.1  fredette #define Dbl_leftshiftby4(dbl_valueA,dbl_valueB) \
    225       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),28,Dallp1(dbl_valueA)); \
    226       1.1  fredette     Dallp2(dbl_valueB) <<= 4
    227       1.1  fredette #define Dbl_leftshiftby3(dbl_valueA,dbl_valueB) \
    228       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),29,Dallp1(dbl_valueA)); \
    229       1.1  fredette     Dallp2(dbl_valueB) <<= 3
    230       1.1  fredette #define Dbl_leftshiftby2(dbl_valueA,dbl_valueB) \
    231       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),30,Dallp1(dbl_valueA)); \
    232       1.1  fredette     Dallp2(dbl_valueB) <<= 2
    233       1.1  fredette #define Dbl_leftshiftby1(dbl_valueA,dbl_valueB) \
    234       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),31,Dallp1(dbl_valueA)); \
    235       1.1  fredette     Dallp2(dbl_valueB) <<= 1
    236       1.1  fredette 
    237       1.1  fredette #define Dbl_rightshiftby8(dbl_valueA,dbl_valueB) \
    238       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),8,Dallp2(dbl_valueB)); \
    239       1.1  fredette     Dallp1(dbl_valueA) >>= 8
    240       1.1  fredette #define Dbl_rightshiftby4(dbl_valueA,dbl_valueB) \
    241       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),4,Dallp2(dbl_valueB)); \
    242       1.1  fredette     Dallp1(dbl_valueA) >>= 4
    243       1.1  fredette #define Dbl_rightshiftby2(dbl_valueA,dbl_valueB) \
    244       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),2,Dallp2(dbl_valueB)); \
    245       1.1  fredette     Dallp1(dbl_valueA) >>= 2
    246       1.1  fredette #define Dbl_rightshiftby1(dbl_valueA,dbl_valueB) \
    247       1.1  fredette     Shiftdouble(Dallp1(dbl_valueA),Dallp2(dbl_valueB),1,Dallp2(dbl_valueB)); \
    248       1.1  fredette     Dallp1(dbl_valueA) >>= 1
    249       1.1  fredette 
    250       1.1  fredette /* This magnitude comparison uses the signless first words and
    251       1.1  fredette  * the regular part2 words.  The comparison is graphically:
    252       1.1  fredette  *
    253       1.1  fredette  *       1st greater?  -------------
    254       1.1  fredette  *				   |
    255       1.1  fredette  *       1st less?-----------------+---------
    256       1.1  fredette  *				   |	    |
    257       1.1  fredette  *       2nd greater or equal----->|	    |
    258       1.1  fredette  *				 False     True
    259       1.1  fredette  */
    260       1.1  fredette #define Dbl_ismagnitudeless(leftB,rightB,signlessleft,signlessright)	\
    261       1.1  fredette       ((signlessleft <= signlessright) &&				\
    262       1.1  fredette        ( (signlessleft < signlessright) || (Dallp2(leftB)<Dallp2(rightB)) ))
    263       1.1  fredette 
    264       1.1  fredette #define Dbl_copytoint_exponentmantissap1(src,dest) \
    265       1.1  fredette     dest = Dexponentmantissap1(src)
    266       1.1  fredette 
    267       1.1  fredette /* A quiet NaN has the high mantissa bit clear and at least on other (in this
    268       1.1  fredette  * case the adjacent bit) bit set. */
    269       1.1  fredette #define Dbl_set_quiet(dbl_value) Deposit_dhigh2mantissa(dbl_value,1)
    270       1.1  fredette #define Dbl_set_exponent(dbl_value, exp) Deposit_dexponent(dbl_value,exp)
    271       1.1  fredette 
    272       1.1  fredette #define Dbl_set_mantissa(desta,destb,valuea,valueb)	\
    273       1.1  fredette     Deposit_dmantissap1(desta,valuea);			\
    274       1.1  fredette     Dmantissap2(destb) = Dmantissap2(valueb)
    275       1.1  fredette #define Dbl_set_mantissap1(desta,valuea)		\
    276       1.1  fredette     Deposit_dmantissap1(desta,valuea)
    277       1.1  fredette #define Dbl_set_mantissap2(destb,valueb)		\
    278       1.1  fredette     Dmantissap2(destb) = Dmantissap2(valueb)
    279       1.1  fredette 
    280       1.1  fredette #define Dbl_set_exponentmantissa(desta,destb,valuea,valueb)	\
    281       1.1  fredette     Deposit_dexponentmantissap1(desta,valuea);			\
    282       1.1  fredette     Dmantissap2(destb) = Dmantissap2(valueb)
    283       1.1  fredette #define Dbl_set_exponentmantissap1(dest,value)			\
    284       1.1  fredette     Deposit_dexponentmantissap1(dest,value)
    285       1.1  fredette 
    286       1.1  fredette #define Dbl_copyfromptr(src,desta,destb) \
    287       1.1  fredette     Dallp1(desta) = src->wd0;		\
    288       1.1  fredette     Dallp2(destb) = src->wd1
    289       1.1  fredette #define Dbl_copytoptr(srca,srcb,dest)	\
    290       1.1  fredette     dest->wd0 = Dallp1(srca);		\
    291       1.1  fredette     dest->wd1 = Dallp2(srcb)
    292       1.1  fredette 
    293       1.1  fredette /*  An infinity is represented with the max exponent and a zero mantissa */
    294       1.1  fredette #define Dbl_setinfinity_exponent(dbl_value) \
    295       1.1  fredette     Deposit_dexponent(dbl_value,DBL_INFINITY_EXPONENT)
    296       1.1  fredette #define Dbl_setinfinity_exponentmantissa(dbl_valueA,dbl_valueB)	\
    297       1.1  fredette     Deposit_dexponentmantissap1(dbl_valueA,			\
    298       1.1  fredette     (DBL_INFINITY_EXPONENT << (32-(1+DBL_EXP_LENGTH))));	\
    299       1.1  fredette     Dmantissap2(dbl_valueB) = 0
    300       1.1  fredette #define Dbl_setinfinitypositive(dbl_valueA,dbl_valueB)		\
    301       1.1  fredette     Dallp1(dbl_valueA)						\
    302       1.1  fredette 	= (DBL_INFINITY_EXPONENT << (32-(1+DBL_EXP_LENGTH)));	\
    303       1.1  fredette     Dmantissap2(dbl_valueB) = 0
    304       1.1  fredette #define Dbl_setinfinitynegative(dbl_valueA,dbl_valueB)		\
    305       1.1  fredette     Dallp1(dbl_valueA) = (1<<31) |				\
    306       1.1  fredette 	(DBL_INFINITY_EXPONENT << (32-(1+DBL_EXP_LENGTH)));	\
    307       1.1  fredette     Dmantissap2(dbl_valueB) = 0
    308       1.1  fredette #define Dbl_setinfinity(dbl_valueA,dbl_valueB,sign)		\
    309       1.1  fredette     Dallp1(dbl_valueA) = (sign << 31) |				\
    310       1.1  fredette 	(DBL_INFINITY_EXPONENT << (32-(1+DBL_EXP_LENGTH)));	\
    311       1.1  fredette     Dmantissap2(dbl_valueB) = 0
    312       1.1  fredette 
    313       1.1  fredette #define Dbl_sethigh4bits(dbl_value, extsign) Deposit_dhigh4p1(dbl_value,extsign)
    314       1.1  fredette #define Dbl_set_sign(dbl_value,sign) Deposit_dsign(dbl_value,sign)
    315       1.1  fredette #define Dbl_invert_sign(dbl_value) Deposit_dsign(dbl_value,~Dsign(dbl_value))
    316       1.1  fredette #define Dbl_setone_sign(dbl_value) Deposit_dsign(dbl_value,1)
    317       1.1  fredette #define Dbl_setone_lowmantissap2(dbl_value) Deposit_dlowp2(dbl_value,1)
    318       1.1  fredette #define Dbl_setzero_sign(dbl_value) Dallp1(dbl_value) &= 0x7fffffff
    319       1.1  fredette #define Dbl_setzero_exponent(dbl_value)			\
    320       1.1  fredette     Dallp1(dbl_value) &= 0x800fffff
    321       1.1  fredette #define Dbl_setzero_mantissa(dbl_valueA,dbl_valueB)	\
    322       1.1  fredette     Dallp1(dbl_valueA) &= 0xfff00000;			\
    323       1.1  fredette     Dallp2(dbl_valueB) = 0
    324       1.1  fredette #define Dbl_setzero_mantissap1(dbl_value) Dallp1(dbl_value) &= 0xfff00000
    325       1.1  fredette #define Dbl_setzero_mantissap2(dbl_value) Dallp2(dbl_value) = 0
    326       1.1  fredette #define Dbl_setzero_exponentmantissa(dbl_valueA,dbl_valueB)	\
    327       1.1  fredette     Dallp1(dbl_valueA) &= 0x80000000;		\
    328       1.1  fredette     Dallp2(dbl_valueB) = 0
    329       1.1  fredette #define Dbl_setzero_exponentmantissap1(dbl_valueA)	\
    330       1.1  fredette     Dallp1(dbl_valueA) &= 0x80000000
    331       1.1  fredette #define Dbl_setzero(dbl_valueA,dbl_valueB) \
    332       1.1  fredette     Dallp1(dbl_valueA) = 0; Dallp2(dbl_valueB) = 0
    333       1.1  fredette #define Dbl_setzerop1(dbl_value) Dallp1(dbl_value) = 0
    334       1.1  fredette #define Dbl_setzerop2(dbl_value) Dallp2(dbl_value) = 0
    335       1.1  fredette #define Dbl_setnegativezero(dbl_value) \
    336       1.1  fredette     Dallp1(dbl_value) = 1 << 31; Dallp2(dbl_value) = 0
    337       1.1  fredette #define Dbl_setnegativezerop1(dbl_value) Dallp1(dbl_value) = 1 << 31
    338       1.1  fredette 
    339       1.1  fredette /* Use the following macro for both overflow & underflow conditions */
    340       1.1  fredette #define ovfl -
    341       1.1  fredette #define unfl +
    342       1.1  fredette #define Dbl_setwrapped_exponent(dbl_value,exponent,op) \
    343       1.1  fredette     Deposit_dexponent(dbl_value,(exponent op DBL_WRAP))
    344       1.1  fredette 
    345       1.1  fredette #define Dbl_setlargestpositive(dbl_valueA,dbl_valueB)			\
    346  1.2.14.1       snj     Dallp1(dbl_valueA) = ((DBL_EMAX+DBL_BIAS) << (32-(1+DBL_EXP_LENGTH))) \
    347       1.1  fredette 			| ((1<<(32-(1+DBL_EXP_LENGTH))) - 1 );		\
    348       1.1  fredette     Dallp2(dbl_valueB) = 0xFFFFFFFF
    349       1.1  fredette #define Dbl_setlargestnegative(dbl_valueA,dbl_valueB)			\
    350  1.2.14.1       snj     Dallp1(dbl_valueA) = ((DBL_EMAX+DBL_BIAS) << (32-(1+DBL_EXP_LENGTH))) \
    351       1.1  fredette 			| ((1<<(32-(1+DBL_EXP_LENGTH))) - 1 ) | (1<<31); \
    352       1.1  fredette     Dallp2(dbl_valueB) = 0xFFFFFFFF
    353       1.1  fredette #define Dbl_setlargest_exponentmantissa(dbl_valueA,dbl_valueB)		\
    354       1.1  fredette     Deposit_dexponentmantissap1(dbl_valueA,				\
    355  1.2.14.1       snj 	(((DBL_EMAX+DBL_BIAS) << (32-(1+DBL_EXP_LENGTH)))		\
    356       1.1  fredette 			| ((1<<(32-(1+DBL_EXP_LENGTH))) - 1 )));	\
    357       1.1  fredette     Dallp2(dbl_valueB) = 0xFFFFFFFF
    358       1.1  fredette 
    359       1.1  fredette #define Dbl_setnegativeinfinity(dbl_valueA,dbl_valueB)			\
    360       1.1  fredette     Dallp1(dbl_valueA) = ((1<<DBL_EXP_LENGTH) | DBL_INFINITY_EXPONENT)	\
    361       1.1  fredette 			 << (32-(1+DBL_EXP_LENGTH)) ;			\
    362       1.1  fredette     Dallp2(dbl_valueB) = 0
    363       1.1  fredette #define Dbl_setlargest(dbl_valueA,dbl_valueB,sign)			\
    364       1.1  fredette     Dallp1(dbl_valueA) = (sign << 31) |					\
    365  1.2.14.1       snj 	((DBL_EMAX+DBL_BIAS) << (32-(1+DBL_EXP_LENGTH))) |		\
    366       1.1  fredette 	 ((1 << (32-(1+DBL_EXP_LENGTH))) - 1 );				\
    367       1.1  fredette     Dallp2(dbl_valueB) = 0xFFFFFFFF
    368       1.1  fredette 
    369       1.1  fredette 
    370       1.1  fredette /* The high bit is always zero so arithmetic or logical shifts will work. */
    371       1.1  fredette #define Dbl_right_align(srcdstA,srcdstB,shift,extent)			\
    372       1.1  fredette     if( shift >= 32 )							\
    373       1.1  fredette 	{								\
    374       1.1  fredette 	/* Big shift requires examining the portion shift off		\
    375       1.1  fredette 	the end to properly set inexact.  */				\
    376       1.1  fredette 	if(shift < 64)							\
    377       1.1  fredette 	    {								\
    378       1.1  fredette 	    if(shift > 32)						\
    379       1.1  fredette 		{							\
    380       1.1  fredette 		Variable_shift_double(Dallp1(srcdstA),Dallp2(srcdstB),	\
    381       1.1  fredette 		 shift-32, Extall(extent));				\
    382       1.1  fredette 		if(Dallp2(srcdstB) << (64 - (shift))) Ext_setone_low(extent); \
    383       1.1  fredette 		}							\
    384       1.1  fredette 	    else Extall(extent) = Dallp2(srcdstB);			\
    385       1.1  fredette 	    Dallp2(srcdstB) = Dallp1(srcdstA) >> (shift - 32);		\
    386       1.1  fredette 	    }								\
    387       1.1  fredette 	else								\
    388       1.1  fredette 	    {								\
    389       1.1  fredette 	    Extall(extent) = Dallp1(srcdstA);				\
    390       1.1  fredette 	    if(Dallp2(srcdstB)) Ext_setone_low(extent);			\
    391       1.1  fredette 	    Dallp2(srcdstB) = 0;					\
    392       1.1  fredette 	    }								\
    393       1.1  fredette 	Dallp1(srcdstA) = 0;						\
    394       1.1  fredette 	}								\
    395       1.1  fredette     else								\
    396       1.1  fredette 	{								\
    397       1.1  fredette 	/* Small alignment is simpler.  Extension is easily set. */	\
    398       1.1  fredette 	if (shift > 0)							\
    399       1.1  fredette 	    {								\
    400       1.1  fredette 	    Extall(extent) = Dallp2(srcdstB) << (32 - (shift));		\
    401       1.1  fredette 	    Variable_shift_double(Dallp1(srcdstA),Dallp2(srcdstB),shift, \
    402       1.1  fredette 	     Dallp2(srcdstB));						\
    403       1.1  fredette 	    Dallp1(srcdstA) >>= shift;					\
    404       1.1  fredette 	    }								\
    405       1.1  fredette 	else Extall(extent) = 0;					\
    406       1.1  fredette 	}
    407       1.1  fredette 
    408       1.1  fredette /*
    409       1.1  fredette  * Here we need to shift the result right to correct for an overshift
    410       1.1  fredette  * (due to the exponent becoming negative) during normalization.
    411       1.1  fredette  */
    412       1.1  fredette #define Dbl_fix_overshift(srcdstA,srcdstB,shift,extent)			\
    413       1.1  fredette 	    Extall(extent) = Dallp2(srcdstB) << (32 - (shift));		\
    414       1.1  fredette 	    Dallp2(srcdstB) = (Dallp1(srcdstA) << (32 - (shift))) |	\
    415       1.1  fredette 		(Dallp2(srcdstB) >> (shift));				\
    416       1.1  fredette 	    Dallp1(srcdstA) = Dallp1(srcdstA) >> shift
    417       1.1  fredette 
    418       1.1  fredette #define Dbl_hiddenhigh3mantissa(dbl_value) Dhiddenhigh3mantissa(dbl_value)
    419       1.1  fredette #define Dbl_hidden(dbl_value) Dhidden(dbl_value)
    420       1.1  fredette #define Dbl_lowmantissap2(dbl_value) Dlowp2(dbl_value)
    421       1.1  fredette 
    422       1.1  fredette /* The left argument is never smaller than the right argument */
    423       1.1  fredette #define Dbl_subtract(lefta,leftb,righta,rightb,resulta,resultb)			\
    424       1.1  fredette     if( Dallp2(rightb) > Dallp2(leftb) ) Dallp1(lefta)--;	\
    425       1.1  fredette     Dallp2(resultb) = Dallp2(leftb) - Dallp2(rightb);		\
    426       1.1  fredette     Dallp1(resulta) = Dallp1(lefta) - Dallp1(righta)
    427       1.1  fredette 
    428       1.1  fredette /* Subtract right augmented with extension from left augmented with zeros and
    429       1.1  fredette  * store into result and extension. */
    430       1.1  fredette #define Dbl_subtract_withextension(lefta,leftb,righta,rightb,extent,resulta,resultb)	\
    431       1.1  fredette     Dbl_subtract(lefta,leftb,righta,rightb,resulta,resultb);		\
    432       1.1  fredette     if( (Extall(extent) = 0-Extall(extent)) )				\
    433       1.1  fredette 	{								\
    434       1.1  fredette 	if((Dallp2(resultb)--) == 0) Dallp1(resulta)--;			\
    435       1.1  fredette 	}
    436       1.1  fredette 
    437       1.1  fredette #define Dbl_addition(lefta,leftb,righta,rightb,resulta,resultb)		\
    438       1.1  fredette     /* If the sum of the low words is less than either source, then	\
    439       1.1  fredette      * an overflow into the next word occurred. */			\
    440       1.1  fredette     Dallp1(resulta) = Dallp1(lefta) + Dallp1(righta);			\
    441       1.1  fredette     if((Dallp2(resultb) = Dallp2(leftb) + Dallp2(rightb)) < Dallp2(rightb)) \
    442       1.1  fredette 	Dallp1(resulta)++
    443       1.1  fredette 
    444       1.1  fredette #define Dbl_xortointp1(left,right,result)			\
    445       1.1  fredette     result = Dallp1(left) XOR Dallp1(right)
    446       1.1  fredette 
    447       1.1  fredette #define Dbl_xorfromintp1(left,right,result)			\
    448       1.1  fredette     Dallp1(result) = left XOR Dallp1(right)
    449       1.1  fredette 
    450       1.1  fredette #define Dbl_swap_lower(left,right)				\
    451       1.1  fredette     Dallp2(left)  = Dallp2(left) XOR Dallp2(right);		\
    452       1.1  fredette     Dallp2(right) = Dallp2(left) XOR Dallp2(right);		\
    453       1.1  fredette     Dallp2(left)  = Dallp2(left) XOR Dallp2(right)
    454       1.1  fredette 
    455       1.1  fredette /* Need to Initialize */
    456       1.1  fredette #define Dbl_makequietnan(desta,destb)					\
    457  1.2.14.1       snj     Dallp1(desta) = ((DBL_EMAX+DBL_BIAS)+1)<< (32-(1+DBL_EXP_LENGTH))	\
    458       1.1  fredette 		| (1<<(32-(1+DBL_EXP_LENGTH+2)));			\
    459       1.1  fredette     Dallp2(destb) = 0
    460       1.1  fredette #define Dbl_makesignalingnan(desta,destb)				\
    461  1.2.14.1       snj     Dallp1(desta) = ((DBL_EMAX+DBL_BIAS)+1)<< (32-(1+DBL_EXP_LENGTH))	\
    462       1.1  fredette 		| (1<<(32-(1+DBL_EXP_LENGTH+1)));			\
    463       1.1  fredette     Dallp2(destb) = 0
    464       1.1  fredette 
    465       1.1  fredette #define Dbl_normalize(dbl_opndA,dbl_opndB,exponent)			\
    466       1.1  fredette 	while(Dbl_iszero_hiddenhigh7mantissa(dbl_opndA)) {		\
    467       1.1  fredette 		Dbl_leftshiftby8(dbl_opndA,dbl_opndB);			\
    468       1.1  fredette 		exponent -= 8;						\
    469       1.1  fredette 	}								\
    470       1.1  fredette 	if(Dbl_iszero_hiddenhigh3mantissa(dbl_opndA)) {			\
    471       1.1  fredette 		Dbl_leftshiftby4(dbl_opndA,dbl_opndB);			\
    472       1.1  fredette 		exponent -= 4;						\
    473       1.1  fredette 	}								\
    474       1.1  fredette 	while(Dbl_iszero_hidden(dbl_opndA)) {				\
    475       1.1  fredette 		Dbl_leftshiftby1(dbl_opndA,dbl_opndB);			\
    476       1.1  fredette 		exponent -= 1;						\
    477       1.1  fredette 	}
    478       1.1  fredette 
    479       1.1  fredette #define Twoword_add(src1dstA,src1dstB,src2A,src2B)		\
    480       1.1  fredette 	/*							\
    481       1.1  fredette 	 * want this macro to generate:				\
    482       1.1  fredette 	 *	ADD	src1dstB,src2B,src1dstB;		\
    483       1.1  fredette 	 *	ADDC	src1dstA,src2A,src1dstA;		\
    484       1.1  fredette 	 */							\
    485       1.1  fredette 	if ((src1dstB) + (src2B) < (src1dstB)) Dallp1(src1dstA)++; \
    486       1.1  fredette 	Dallp1(src1dstA) += (src2A);				\
    487       1.1  fredette 	Dallp2(src1dstB) += (src2B)
    488       1.1  fredette 
    489       1.1  fredette #define Twoword_subtract(src1dstA,src1dstB,src2A,src2B)		\
    490       1.1  fredette 	/*							\
    491       1.1  fredette 	 * want this macro to generate:				\
    492       1.1  fredette 	 *	SUB	src1dstB,src2B,src1dstB;		\
    493       1.1  fredette 	 *	SUBB	src1dstA,src2A,src1dstA;		\
    494       1.1  fredette 	 */							\
    495       1.1  fredette 	if ((src1dstB) < (src2B)) Dallp1(src1dstA)--;		\
    496       1.1  fredette 	Dallp1(src1dstA) -= (src2A);				\
    497       1.1  fredette 	Dallp2(src1dstB) -= (src2B)
    498       1.1  fredette 
    499       1.1  fredette #define Dbl_setoverflow(resultA,resultB)				\
    500       1.1  fredette 	/* set result to infinity or largest number */			\
    501       1.1  fredette 	switch (Rounding_mode()) {					\
    502       1.1  fredette 		case ROUNDPLUS:						\
    503       1.1  fredette 			if (Dbl_isone_sign(resultA)) {			\
    504       1.1  fredette 				Dbl_setlargestnegative(resultA,resultB); \
    505       1.1  fredette 			}						\
    506       1.1  fredette 			else {						\
    507       1.1  fredette 				Dbl_setinfinitypositive(resultA,resultB); \
    508       1.1  fredette 			}						\
    509       1.1  fredette 			break;						\
    510       1.1  fredette 		case ROUNDMINUS:					\
    511       1.1  fredette 			if (Dbl_iszero_sign(resultA)) {			\
    512       1.1  fredette 				Dbl_setlargestpositive(resultA,resultB); \
    513       1.1  fredette 			}						\
    514       1.1  fredette 			else {						\
    515       1.1  fredette 				Dbl_setinfinitynegative(resultA,resultB); \
    516       1.1  fredette 			}						\
    517       1.1  fredette 			break;						\
    518       1.1  fredette 		case ROUNDNEAREST:					\
    519       1.1  fredette 			Dbl_setinfinity_exponentmantissa(resultA,resultB); \
    520       1.1  fredette 			break;						\
    521       1.1  fredette 		case ROUNDZERO:						\
    522       1.1  fredette 			Dbl_setlargest_exponentmantissa(resultA,resultB); \
    523       1.1  fredette 	}
    524       1.1  fredette 
    525       1.1  fredette #define Dbl_denormalize(opndp1,opndp2,exponent,guard,sticky,inexact)	\
    526       1.1  fredette     Dbl_clear_signexponent_set_hidden(opndp1);				\
    527       1.1  fredette     if (exponent >= (1-DBL_P)) {					\
    528       1.1  fredette 	if (exponent >= -31) {						\
    529       1.1  fredette 	    guard = (Dallp2(opndp2) >> (-(exponent))) & 1;		\
    530       1.1  fredette 	    if (exponent < 0) sticky |= Dallp2(opndp2) << (32+exponent); \
    531       1.1  fredette 	    if (exponent > -31) {					\
    532       1.1  fredette 		Variable_shift_double(opndp1,opndp2,1-exponent,opndp2);	\
    533       1.1  fredette 		Dallp1(opndp1) >>= 1-exponent;				\
    534       1.1  fredette 	    }								\
    535       1.1  fredette 	    else {							\
    536       1.1  fredette 		Dallp2(opndp2) = Dallp1(opndp1);			\
    537       1.1  fredette 		Dbl_setzerop1(opndp1);					\
    538       1.1  fredette 	    }								\
    539       1.1  fredette 	}								\
    540       1.1  fredette 	else {								\
    541       1.1  fredette 	    guard = (Dallp1(opndp1) >> (-32-(exponent))) & 1;		\
    542       1.1  fredette 	    if (exponent == -32) sticky |= Dallp2(opndp2);		\
    543       1.1  fredette 	    else sticky |= (Dallp2(opndp2) | Dallp1(opndp1) << (64+(exponent))); \
    544       1.1  fredette 	    Dallp2(opndp2) = Dallp1(opndp1) >> (-31-(exponent));	\
    545       1.1  fredette 	    Dbl_setzerop1(opndp1);					\
    546       1.1  fredette 	}								\
    547       1.1  fredette 	inexact = guard | sticky;					\
    548       1.1  fredette     }									\
    549       1.1  fredette     else {								\
    550       1.1  fredette 	guard = 0;							\
    551       1.1  fredette 	sticky |= (Dallp1(opndp1) | Dallp2(opndp2));			\
    552       1.1  fredette 	Dbl_setzero(opndp1,opndp2);					\
    553       1.1  fredette 	inexact = sticky;						\
    554       1.1  fredette     }
    555       1.1  fredette 
    556       1.1  fredette 
    557       1.2     skrll int dbl_fadd(dbl_floating_point *, dbl_floating_point*, dbl_floating_point*, unsigned int *);
    558       1.2     skrll int dbl_fcmp(dbl_floating_point *, dbl_floating_point*, unsigned int, unsigned int *);
    559       1.2     skrll int dbl_fdiv(dbl_floating_point *, dbl_floating_point *, dbl_floating_point *, unsigned int *);
    560       1.2     skrll int dbl_fmpy(dbl_floating_point *, dbl_floating_point *, dbl_floating_point*, unsigned int *);
    561       1.2     skrll int dbl_frem(dbl_floating_point *, dbl_floating_point *, dbl_floating_point*, unsigned int *);
    562       1.2     skrll int dbl_fsqrt(dbl_floating_point *, dbl_floating_point *, unsigned int *);
    563       1.2     skrll int dbl_fsub(dbl_floating_point *, dbl_floating_point *, dbl_floating_point*, unsigned int *);
    564       1.1  fredette 
    565       1.2     skrll dbl_floating_point dbl_setoverflow(unsigned int);
    566       1.1  fredette 
    567       1.2     skrll int sgl_to_dbl_fcnvff(sgl_floating_point *, dbl_floating_point *, unsigned int *);
    568       1.2     skrll int dbl_to_sgl_fcnvff(dbl_floating_point *, sgl_floating_point *, unsigned int *);
    569       1.1  fredette 
    570       1.2     skrll int dbl_frnd(dbl_floating_point *, dbl_floating_point *, unsigned int *);
    571