xref: /onnv-gate/usr/src/lib/libast/common/comp/frexpl.c (revision 12068:08a39a083754)
14887Schin /***********************************************************************
24887Schin *                                                                      *
34887Schin *               This software is part of the ast package               *
4*12068SRoger.Faulkner@Oracle.COM *          Copyright (c) 1985-2010 AT&T Intellectual Property          *
54887Schin *                      and is licensed under the                       *
64887Schin *                  Common Public License, Version 1.0                  *
78462SApril.Chin@Sun.COM *                    by AT&T Intellectual Property                     *
84887Schin *                                                                      *
94887Schin *                A copy of the License is available at                 *
104887Schin *            http://www.opensource.org/licenses/cpl1.0.txt             *
114887Schin *         (with md5 checksum 059e8cd6165cb4c31e351f2b69388fd9)         *
124887Schin *                                                                      *
134887Schin *              Information and Software Systems Research               *
144887Schin *                            AT&T Research                             *
154887Schin *                           Florham Park NJ                            *
164887Schin *                                                                      *
174887Schin *                 Glenn Fowler <gsf@research.att.com>                  *
184887Schin *                  David Korn <dgk@research.att.com>                   *
194887Schin *                   Phong Vo <kpv@research.att.com>                    *
204887Schin *                                                                      *
214887Schin ***********************************************************************/
224887Schin #pragma prototyped
234887Schin 
244887Schin /*
254887Schin  * frexpl/ldexpl implementation
264887Schin  */
274887Schin 
284887Schin #include <ast.h>
294887Schin 
304887Schin #include "FEATURE/float"
314887Schin 
324887Schin #if _lib_frexpl && _lib_ldexpl
334887Schin 
344887Schin NoN(frexpl)
354887Schin 
364887Schin #else
374887Schin 
384887Schin #ifndef LDBL_MAX_EXP
394887Schin #define LDBL_MAX_EXP	DBL_MAX_EXP
404887Schin #endif
414887Schin 
424887Schin #if defined(_ast_fltmax_exp_index) && defined(_ast_fltmax_exp_shift)
434887Schin 
444887Schin #define INIT()		_ast_fltmax_exp_t _pow_
454887Schin #define pow2(i)		(_pow_.f=1,_pow_.e[_ast_fltmax_exp_index]+=((i)<<_ast_fltmax_exp_shift),_pow_.f)
464887Schin 
474887Schin #else
484887Schin 
494887Schin static _ast_fltmax_t	pow2tab[LDBL_MAX_EXP + 1];
504887Schin 
514887Schin static int
524887Schin init(void)
534887Schin {
544887Schin 	register int		x;
554887Schin 	_ast_fltmax_t		g;
564887Schin 
574887Schin 	g = 1;
584887Schin 	for (x = 0; x < elementsof(pow2tab); x++)
594887Schin 	{
604887Schin 		pow2tab[x] = g;
614887Schin 		g *= 2;
624887Schin 	}
634887Schin 	return 0;
644887Schin }
654887Schin 
664887Schin #define INIT()		(pow2tab[0]?0:init())
674887Schin 
684887Schin #define pow2(i)		(pow2tab[i])
694887Schin 
704887Schin #endif
714887Schin 
724887Schin #if !_lib_frexpl
734887Schin 
744887Schin #undef	frexpl
754887Schin 
764887Schin extern _ast_fltmax_t
774887Schin frexpl(_ast_fltmax_t f, int* p)
784887Schin {
794887Schin 	register int		k;
804887Schin 	register int		x;
814887Schin 	_ast_fltmax_t		g;
824887Schin 
834887Schin 	INIT();
844887Schin 
854887Schin 	/*
864887Schin 	 * normalize
874887Schin 	 */
884887Schin 
894887Schin 	x = k = LDBL_MAX_EXP / 2;
904887Schin 	if (f < 1)
914887Schin 	{
924887Schin 		g = 1.0L / f;
934887Schin 		for (;;)
944887Schin 		{
954887Schin 			k = (k + 1) / 2;
964887Schin 			if (g < pow2(x))
974887Schin 				x -= k;
984887Schin 			else if (k == 1 && g < pow2(x+1))
994887Schin 				break;
1004887Schin 			else
1014887Schin 				x += k;
1024887Schin 		}
1034887Schin 		if (g == pow2(x))
1044887Schin 			x--;
1054887Schin 		x = -x;
1064887Schin 	}
1074887Schin 	else if (f > 1)
1084887Schin 	{
1094887Schin 		for (;;)
1104887Schin 		{
1114887Schin 			k = (k + 1) / 2;
1124887Schin 			if (f > pow2(x))
1134887Schin 				x += k;
1144887Schin 			else if (k == 1 && f > pow2(x-1))
1154887Schin 				break;
1164887Schin 			else
1174887Schin 				x -= k;
1184887Schin 		}
1194887Schin 		if (f == pow2(x))
1204887Schin 			x++;
1214887Schin 	}
1224887Schin 	else
1234887Schin 		x = 1;
1244887Schin 	*p = x;
1254887Schin 
1264887Schin 	/*
1274887Schin 	 * shift
1284887Schin 	 */
1294887Schin 
1304887Schin 	x = -x;
1314887Schin 	if (x < 0)
1324887Schin 		f /= pow2(-x);
1334887Schin 	else if (x < LDBL_MAX_EXP)
1344887Schin 		f *= pow2(x);
1354887Schin 	else
1364887Schin 		f = (f * pow2(LDBL_MAX_EXP - 1)) * pow2(x - (LDBL_MAX_EXP - 1));
1374887Schin 	return f;
1384887Schin }
1394887Schin 
1404887Schin #endif
1414887Schin 
1424887Schin #if !_lib_ldexpl
1434887Schin 
1444887Schin #undef	ldexpl
1454887Schin 
1464887Schin extern _ast_fltmax_t
1474887Schin ldexpl(_ast_fltmax_t f, register int x)
1484887Schin {
1494887Schin 	INIT();
1504887Schin 	if (x < 0)
1514887Schin 		f /= pow2(-x);
1524887Schin 	else if (x < LDBL_MAX_EXP)
1534887Schin 		f *= pow2(x);
1544887Schin 	else
1554887Schin 		f = (f * pow2(LDBL_MAX_EXP - 1)) * pow2(x - (LDBL_MAX_EXP - 1));
1564887Schin 	return f;
1574887Schin }
1584887Schin 
1594887Schin #endif
1604887Schin 
1614887Schin #endif
162