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