xref: /illumos-gate/usr/src/contrib/ast/src/lib/libast/comp/frexpl.c (revision b30d193948be5a7794d7ae3ba0ed9c2f72c88e0f)
1 /***********************************************************************
2 *                                                                      *
3 *               This software is part of the ast package               *
4 *          Copyright (c) 1985-2011 AT&T Intellectual Property          *
5 *                      and is licensed under the                       *
6 *                 Eclipse Public License, Version 1.0                  *
7 *                    by AT&T Intellectual Property                     *
8 *                                                                      *
9 *                A copy of the License is available at                 *
10 *          http://www.eclipse.org/org/documents/epl-v10.html           *
11 *         (with md5 checksum b35adb5213ca9657e911e9befb180842)         *
12 *                                                                      *
13 *              Information and Software Systems Research               *
14 *                            AT&T Research                             *
15 *                           Florham Park NJ                            *
16 *                                                                      *
17 *                 Glenn Fowler <gsf@research.att.com>                  *
18 *                  David Korn <dgk@research.att.com>                   *
19 *                   Phong Vo <kpv@research.att.com>                    *
20 *                                                                      *
21 ***********************************************************************/
22 #pragma prototyped
23 
24 /*
25  * frexpl/ldexpl implementation
26  */
27 
28 #include <ast.h>
29 
30 #include "FEATURE/float"
31 
32 #if _lib_frexpl && _lib_ldexpl
33 
34 NoN(frexpl)
35 
36 #else
37 
38 #ifndef LDBL_MAX_EXP
39 #define LDBL_MAX_EXP	DBL_MAX_EXP
40 #endif
41 
42 #if defined(_ast_fltmax_exp_index) && defined(_ast_fltmax_exp_shift)
43 
44 #define INIT()		_ast_fltmax_exp_t _pow_
45 #define pow2(i)		(_pow_.f=1,_pow_.e[_ast_fltmax_exp_index]+=((i)<<_ast_fltmax_exp_shift),_pow_.f)
46 
47 #else
48 
49 static _ast_fltmax_t	pow2tab[LDBL_MAX_EXP + 1];
50 
51 static int
52 init(void)
53 {
54 	register int		x;
55 	_ast_fltmax_t		g;
56 
57 	g = 1;
58 	for (x = 0; x < elementsof(pow2tab); x++)
59 	{
60 		pow2tab[x] = g;
61 		g *= 2;
62 	}
63 	return 0;
64 }
65 
66 #define INIT()		(pow2tab[0]?0:init())
67 
68 #define pow2(i)		(pow2tab[i])
69 
70 #endif
71 
72 #if !_lib_frexpl
73 
74 #undef	frexpl
75 
76 extern _ast_fltmax_t
77 frexpl(_ast_fltmax_t f, int* p)
78 {
79 	register int		k;
80 	register int		x;
81 	_ast_fltmax_t		g;
82 
83 	INIT();
84 
85 	/*
86 	 * normalize
87 	 */
88 
89 	x = k = LDBL_MAX_EXP / 2;
90 	if (f < 1)
91 	{
92 		g = 1.0L / f;
93 		for (;;)
94 		{
95 			k = (k + 1) / 2;
96 			if (g < pow2(x))
97 				x -= k;
98 			else if (k == 1 && g < pow2(x+1))
99 				break;
100 			else
101 				x += k;
102 		}
103 		if (g == pow2(x))
104 			x--;
105 		x = -x;
106 	}
107 	else if (f > 1)
108 	{
109 		for (;;)
110 		{
111 			k = (k + 1) / 2;
112 			if (f > pow2(x))
113 				x += k;
114 			else if (k == 1 && f > pow2(x-1))
115 				break;
116 			else
117 				x -= k;
118 		}
119 		if (f == pow2(x))
120 			x++;
121 	}
122 	else
123 		x = 1;
124 	*p = x;
125 
126 	/*
127 	 * shift
128 	 */
129 
130 	x = -x;
131 	if (x < 0)
132 		f /= pow2(-x);
133 	else if (x < LDBL_MAX_EXP)
134 		f *= pow2(x);
135 	else
136 		f = (f * pow2(LDBL_MAX_EXP - 1)) * pow2(x - (LDBL_MAX_EXP - 1));
137 	return f;
138 }
139 
140 #endif
141 
142 #if !_lib_ldexpl
143 
144 #undef	ldexpl
145 
146 extern _ast_fltmax_t
147 ldexpl(_ast_fltmax_t f, register int x)
148 {
149 	INIT();
150 	if (x < 0)
151 		f /= pow2(-x);
152 	else if (x < LDBL_MAX_EXP)
153 		f *= pow2(x);
154 	else
155 		f = (f * pow2(LDBL_MAX_EXP - 1)) * pow2(x - (LDBL_MAX_EXP - 1));
156 	return f;
157 }
158 
159 #endif
160 
161 #endif
162