|
|
1.1 ! root 1: /* ! 2: * Copyright (c) 1985 Regents of the University of California. ! 3: * ! 4: * Use and reproduction of this software are granted in accordance with ! 5: * the terms and conditions specified in the Berkeley Software License ! 6: * Agreement (in particular, this entails acknowledgement of the programs' ! 7: * source, and inclusion of this notice) with the additional understanding ! 8: * that all recipients should regard themselves as participants in an ! 9: * ongoing research project and hence should feel obligated to report ! 10: * their experiences (good or bad) with these elementary function codes, ! 11: * using "sendbug 4bsd-bugs@BERKELEY", to the authors. ! 12: */ ! 13: ! 14: #ifndef lint ! 15: static char sccsid[] = "@(#)log10.c 1.2 (Berkeley) 8/21/85"; ! 16: #endif not lint ! 17: ! 18: /* LOG10(X) ! 19: * RETURN THE BASE 10 LOGARITHM OF x ! 20: * DOUBLE PRECISION (VAX D format 56 bits, IEEE DOUBLE 53 BITS) ! 21: * CODED IN C BY K.C. NG, 1/20/85; ! 22: * REVISED BY K.C. NG on 1/23/85, 3/7/85, 4/16/85. ! 23: * ! 24: * Required kernel function: ! 25: * log(x) ! 26: * ! 27: * Method : ! 28: * log(x) ! 29: * log10(x) = --------- or [1/log(10)]*log(x) ! 30: * log(10) ! 31: * ! 32: * Note: ! 33: * [log(10)] rounded to 56 bits has error .0895 ulps, ! 34: * [1/log(10)] rounded to 53 bits has error .198 ulps; ! 35: * therefore, for better accuracy, in VAX D format, we divide ! 36: * log(x) by log(10), but in IEEE Double format, we multiply ! 37: * log(x) by [1/log(10)]. ! 38: * ! 39: * Special cases: ! 40: * log10(x) is NaN with signal if x < 0; ! 41: * log10(+INF) is +INF with no signal; log10(0) is -INF with signal; ! 42: * log10(NaN) is that NaN with no signal. ! 43: * ! 44: * Accuracy: ! 45: * log10(X) returns the exact log10(x) nearly rounded. In a test run ! 46: * with 1,536,000 random arguments on a VAX, the maximum observed ! 47: * error was 1.74 ulps (units in the last place). ! 48: * ! 49: * Constants: ! 50: * The hexadecimal values are the intended ones for the following constants. ! 51: * The decimal values may be used, provided that the compiler will convert ! 52: * from decimal to binary accurately enough to produce the hexadecimal values ! 53: * shown. ! 54: */ ! 55: ! 56: #ifdef VAX /* VAX D format (56 bits) */ ! 57: /* static double */ ! 58: /* ln10hi = 2.3025850929940456790E0 ; Hex 2^ 2 * .935D8DDDAAA8AC */ ! 59: static long ln10hix[] = { 0x5d8d4113, 0xa8acddaa}; ! 60: #define ln10hi (*(double*)ln10hix) ! 61: #else /* IEEE double */ ! 62: static double ! 63: ivln10 = 4.3429448190325181667E-1 ; /*Hex 2^ -2 * 1.BCB7B1526E50E */ ! 64: #endif ! 65: ! 66: double log10(x) ! 67: double x; ! 68: { ! 69: double log(); ! 70: ! 71: #ifdef VAX ! 72: return(log(x)/ln10hi); ! 73: #else /* IEEE double */ ! 74: return(ivln10*log(x)); ! 75: #endif ! 76: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.