Back to index

glibc  2.9
Functions | Variables
support.c File Reference
#include "mathimpl.h"

Go to the source code of this file.

Functions

double scalb (double x, int N)
double copysign (double x, double y)
double logb (double x)
 finite (double x)
double drem (double x, double p)
double sqrt (double x)

Variables

static char sccsid [] = "@(#)support.c 8.1 (Berkeley) 6/4/93"
static const unsigned short msign = 0x7fff
static const unsigned short mexp = 0x7ff0
static const short prep1 = 54
static const short gap = 4
static const short bias = 1023
static const double novf = 1.7E308
static const double nunf = 3.0E-308
static const double zero = 0.0

Function Documentation

double copysign ( double  x,
double  y 
)

Definition at line 126 of file support.c.

{
#ifdef national
        unsigned short  *px=(unsigned short *) &x+3,
                        *py=(unsigned short *) &y+3;
#else  /* national */
        unsigned short  *px=(unsigned short *) &x,
                        *py=(unsigned short *) &y;
#endif /* national */

#if defined(vax)||defined(tahoe)
        if ( (*px & mexp) == 0 ) return(x);
#endif /* defined(vax)||defined(tahoe) */

        *px = ( *px & msign ) | ( *py & ~msign );
        return(x);
}

Here is the caller graph for this function:

double drem ( double  x,
double  p 
)

Definition at line 186 of file support.c.

{
        short sign;
        double hp,dp,tmp;
        unsigned short  k;
#ifdef national
        unsigned short
              *px=(unsigned short *) &x  +3,
              *pp=(unsigned short *) &p  +3,
              *pd=(unsigned short *) &dp +3,
              *pt=(unsigned short *) &tmp+3;
#else  /* national */
        unsigned short
              *px=(unsigned short *) &x  ,
              *pp=(unsigned short *) &p  ,
              *pd=(unsigned short *) &dp ,
              *pt=(unsigned short *) &tmp;
#endif /* national */

        *pp &= msign ;

#if defined(vax)||defined(tahoe)
        if( ( *px & mexp ) == ~msign )    /* is x a reserved operand? */
#else  /* defined(vax)||defined(tahoe) */
        if( ( *px & mexp ) == mexp )
#endif /* defined(vax)||defined(tahoe) */
              return  (x-p)-(x-p); /* create nan if x is inf */
       if (p == zero) {
#if defined(vax)||defined(tahoe)
              return(infnan(EDOM));
#else  /* defined(vax)||defined(tahoe) */
              return zero/zero;
#endif /* defined(vax)||defined(tahoe) */
       }

#if defined(vax)||defined(tahoe)
        if( ( *pp & mexp ) == ~msign )    /* is p a reserved operand? */
#else  /* defined(vax)||defined(tahoe) */
        if( ( *pp & mexp ) == mexp )
#endif /* defined(vax)||defined(tahoe) */
              { if (p != p) return p; else return x;}

        else  if ( ((*pp & mexp)>>gap) <= 1 )
                /* subnormal p, or almost subnormal p */
            { double b; b=scalb(1.0,(int)prep1);
              p *= b; x = drem(x,p); x *= b; return(drem(x,p)/b);}
        else  if ( p >= novf/2)
            { p /= 2 ; x /= 2; return(drem(x,p)*2);}
        else
            {
                dp=p+p; hp=p/2;
                sign= *px & ~msign ;
                *px &= msign       ;
                while ( x > dp )
                    {
                        k=(*px & mexp) - (*pd & mexp) ;
                        tmp = dp ;
                        *pt += k ;

#if defined(vax)||defined(tahoe)
                        if( x < tmp ) *pt -= 128 ;
#else  /* defined(vax)||defined(tahoe) */
                        if( x < tmp ) *pt -= 16 ;
#endif /* defined(vax)||defined(tahoe) */

                        x -= tmp ;
                    }
                if ( x > hp )
                    { x -= p ;  if ( x >= hp ) x -= p ; }

#if defined(vax)||defined(tahoe)
              if (x)
#endif /* defined(vax)||defined(tahoe) */
                     *px ^= sign;
                return( x);

            }
}

Here is the call graph for this function:

Here is the caller graph for this function:

finite ( double  x)

Definition at line 172 of file support.c.

{
#if defined(vax)||defined(tahoe)
        return(1);
#else  /* defined(vax)||defined(tahoe) */
#ifdef national
        return( (*((short *) &x+3 ) & mexp ) != mexp );
#else  /* national */
        return( (*((short *) &x ) & mexp ) != mexp );
#endif /* national */
#endif /* defined(vax)||defined(tahoe) */
}

Here is the caller graph for this function:

double logb ( double  x)

Definition at line 145 of file support.c.

{

#ifdef national
        short *px=(short *) &x+3, k;
#else  /* national */
        short *px=(short *) &x, k;
#endif /* national */

#if defined(vax)||defined(tahoe)
        return (int)(((*px&mexp)>>gap)-bias);
#else  /* defined(vax)||defined(tahoe) */
        if( (k= *px & mexp ) != mexp )
            if ( k != 0 )
                return ( (k>>gap) - bias );
            else if( x != zero)
                return ( -1022.0 );
            else
                return(-(1.0/zero));
        else if(x != x)
            return(x);
        else
            {*px &= msign; return(x);}
#endif /* defined(vax)||defined(tahoe) */
}

Here is the caller graph for this function:

double scalb ( double  x,
int  N 
)

Definition at line 85 of file support.c.

{
        int k;

#ifdef national
        unsigned short *px=(unsigned short *) &x + 3;
#else  /* national */
        unsigned short *px=(unsigned short *) &x;
#endif /* national */

        if( x == zero )  return(x);

#if defined(vax)||defined(tahoe)
        if( (k= *px & mexp ) != ~msign ) {
            if (N < -260)
              return(nunf*nunf);
           else if (N > 260) {
              return(copysign(infnan(ERANGE),x));
           }
#else  /* defined(vax)||defined(tahoe) */
        if( (k= *px & mexp ) != mexp ) {
            if( N<-2100) return(nunf*nunf); else if(N>2100) return(novf+novf);
            if( k == 0 ) {
                 x *= scalb(1.0,(int)prep1);  N -= prep1; return(scalb(x,N));}
#endif /* defined(vax)||defined(tahoe) */

            if((k = (k>>gap)+ N) > 0 )
                if( k < (mexp>>gap) ) *px = (*px&~mexp) | (k<<gap);
                else x=novf+novf;               /* overflow */
            else
                if( k > -prep1 )
                                        /* gradual underflow */
                    {*px=(*px&~mexp)|(short)(1<<gap); x *= scalb(1.0,k-1);}
                else
                return(nunf*nunf);
            }
        return(x);
}

Here is the call graph for this function:

Here is the caller graph for this function:

double sqrt ( double  x)

Definition at line 267 of file support.c.

{
        double q,s,b,r;
        double t;
       double const zero=0.0;
        int m,n,i;
#if defined(vax)||defined(tahoe)
        int k=54;
#else  /* defined(vax)||defined(tahoe) */
        int k=51;
#endif /* defined(vax)||defined(tahoe) */

    /* sqrt(NaN) is NaN, sqrt(+-0) = +-0 */
        if(x!=x||x==zero) return(x);

    /* sqrt(negative) is invalid */
        if(x<zero) {
#if defined(vax)||defined(tahoe)
              return (infnan(EDOM));      /* NaN */
#else  /* defined(vax)||defined(tahoe) */
              return(zero/zero);
#endif /* defined(vax)||defined(tahoe) */
       }

    /* sqrt(INF) is INF */
        if(!finite(x)) return(x);

    /* scale x to [1,4) */
        n=logb(x);
        x=scalb(x,-n);
        if((m=logb(x))!=0) x=scalb(x,-m);       /* subnormal number */
        m += n;
        n = m/2;
        if((n+n)!=m) {x *= 2; m -=1; n=m/2;}

    /* generate sqrt(x) bit by bit (accumulating in q) */
            q=1.0; s=4.0; x -= 1.0; r=1;
            for(i=1;i<=k;i++) {
                t=s+1; x *= 4; r /= 2;
                if(t<=x) {
                    s=t+t+2, x -= t; q += r;}
                else
                    s *= 2;
                }

    /* generate the last bit and determine the final rounding */
            r/=2; x *= 4;
            if(x==zero) goto end; 100+r; /* trigger inexact flag */
            if(s<x) {
                q+=r; x -=s; s += 2; s *= 2; x *= 4;
                t = (x-s)-5;
                b=1.0+3*r/4; if(b==1.0) goto end; /* b==1 : Round-to-zero */
                b=1.0+r/4;   if(b>1.0) t=1;      /* b>1 : Round-to-(+INF) */
                if(t>=0) q+=r; }         /* else: Round-to-nearest */
            else {
                s *= 2; x *= 4;
                t = (x-s)-1;
                b=1.0+3*r/4; if(b==1.0) goto end;
                b=1.0+r/4;   if(b>1.0) t=1;
                if(t>=0) q+=r; }

end:        return(scalb(q,n));
}

Here is the call graph for this function:


Variable Documentation

const short bias = 1023 [static]

Definition at line 81 of file support.c.

const short gap = 4 [static]

Definition at line 81 of file support.c.

const unsigned short mexp = 0x7ff0 [static]

Definition at line 80 of file support.c.

const unsigned short msign = 0x7fff [static]

Definition at line 80 of file support.c.

const double novf = 1.7E308 [static]

Definition at line 82 of file support.c.

const double nunf = 3.0E-308 [static]

Definition at line 82 of file support.c.

const short prep1 = 54 [static]

Definition at line 81 of file support.c.

char sccsid[] = "@(#)support.c 8.1 (Berkeley) 6/4/93" [static]

Definition at line 31 of file support.c.

const double zero = 0.0 [static]

Definition at line 82 of file support.c.