GCC Code Coverage Report
Directory: ./ Exec Total Coverage
File: lib/libm/src/s_scalbnl.c Lines: 13 23 56.5 %
Date: 2017-11-13 Branches: 6 16 37.5 %

Line Branch Exec Source
1
/*	$OpenBSD: s_scalbnl.c,v 1.4 2016/09/12 19:47:02 guenther Exp $	*/
2
/* @(#)s_scalbn.c 5.1 93/09/24 */
3
/*
4
 * ====================================================
5
 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6
 *
7
 * Developed at SunPro, a Sun Microsystems, Inc. business.
8
 * Permission to use, copy, modify, and distribute this
9
 * software is freely granted, provided that this notice
10
 * is preserved.
11
 * ====================================================
12
 */
13
14
/*
15
 * scalbnl (long double x, int n)
16
 * scalbnl(x,n) returns x* 2**n  computed by  exponent
17
 * manipulation rather than by actually performing an
18
 * exponentiation or a multiplication.
19
 */
20
21
/*
22
 * We assume that a long double has a 15-bit exponent.  On systems
23
 * where long double is the same as double, scalbnl() is an alias
24
 * for scalbn(), so we don't use this routine.
25
 */
26
27
#include <sys/types.h>
28
#include <machine/ieee.h>
29
#include <float.h>
30
#include <math.h>
31
32
#if LDBL_MAX_EXP != 0x4000
33
#error "Unsupported long double format"
34
#endif
35
36
static const long double
37
huge = 0x1p16000L,
38
tiny = 0x1p-16000L;
39
40
long double
41
scalbnl (long double x, int n)
42
{
43
1636
	union {
44
		long double e;
45
		struct ieee_ext bits;
46
	} u;
47
	int k;
48
818
	u.e = x;
49
818
        k = u.bits.ext_exp;			/* extract exponent */
50
818
        if (k==0) {				/* 0 or subnormal x */
51
118
            if ((u.bits.ext_frach
52
#ifdef EXT_FRACHMBITS
53
		| u.bits.ext_frachm
54
#endif /* EXT_FRACHMBITS */
55
#ifdef EXT_FRACLMBITS
56
		| u.bits.ext_fraclm
57
#endif /* EXT_FRACLMBITS */
58
118
		| u.bits.ext_fracl)==0) return x;	/* +-0 */
59
	    u.e *= 0x1p+128;
60
	    k = u.bits.ext_exp - 128;
61
            if (n< -50000) return tiny*x; 	/*underflow*/
62
	    }
63
759
        if (k==0x7fff) return x+x;		/* NaN or Inf */
64
759
        k = k+n;
65
759
        if (k >= 0x7fff) return huge*copysignl(huge,x); /* overflow  */
66
759
        if (k > 0) 				/* normal result */
67
759
	    {u.bits.ext_exp = k; return u.e;}
68
        if (k <= -128) {
69
            if (n > 50000) 	/* in case integer overflow in n+k */
70
		return huge*copysign(huge,x);	/*overflow*/
71
	    else return tiny*copysign(tiny,x); 	/*underflow*/
72
	}
73
        k += 128;				/* subnormal result */
74
	u.bits.ext_exp = k;
75
        return u.e*0x1p-128;
76
818
}
77
DEF_STD(scalbnl);
78
79
long double
80
ldexpl(long double x, int n)
81
{
82
1076
	return scalbnl(x, n);
83
}
84
DEF_STD(ldexpl);