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) 08/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
log10(x)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 }
77