summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorOri Bernstein <ori@eigenstate.org>2021-01-30 09:19:57 -0800
committerOri Bernstein <ori@eigenstate.org>2021-01-30 09:19:57 -0800
commitbb250c4c3b8a0add2a7ec2b972f327b0a3db03d4 (patch)
tree80c77fc5c4f304272664c06124507b407411df4e
parent273c4bff7ae59ff33e0b58c95967397090fd3864 (diff)
downloadplan9front-bb250c4c3b8a0add2a7ec2b972f327b0a3db03d4.tar.xz
ape/libm: implement log2 in libc
-rw-r--r--sys/include/ape/math.h1
-rw-r--r--sys/src/ape/lib/ap/math/log.c40
2 files changed, 26 insertions, 15 deletions
diff --git a/sys/include/ape/math.h b/sys/include/ape/math.h
index 3bd3daade..87dd49a2e 100644
--- a/sys/include/ape/math.h
+++ b/sys/include/ape/math.h
@@ -25,6 +25,7 @@ extern double exp(double);
extern double frexp(double, int *);
extern double ldexp(double, int);
extern double log(double);
+extern double log2(double);
extern double log10(double);
extern double modf(double, double *);
extern double pow(double, double);
diff --git a/sys/src/ape/lib/ap/math/log.c b/sys/src/ape/lib/ap/math/log.c
index bc1c2e1cb..2c7eac276 100644
--- a/sys/src/ape/lib/ap/math/log.c
+++ b/sys/src/ape/lib/ap/math/log.c
@@ -10,16 +10,17 @@
#include <math.h>
#include <errno.h>
-#define log2 0.693147180559945309e0
-#define ln10o1 .4342944819032518276511
-#define sqrto2 0.707106781186547524e0
-#define p0 -.240139179559210510e2
-#define p1 0.309572928215376501e2
-#define p2 -.963769093377840513e1
-#define p3 0.421087371217979714e0
-#define q0 -.120069589779605255e2
-#define q1 0.194809660700889731e2
-#define q2 -.891110902798312337e1
+#define Log2 0.693147180559945309e0
+#define Ln2o1 1.4426950408889634073599
+#define Ln10o1 .4342944819032518276511
+#define Sqrto2 0.707106781186547524e0
+#define P0 -.240139179559210510e2
+#define P1 0.309572928215376501e2
+#define P2 -.963769093377840513e1
+#define P3 0.421087371217979714e0
+#define Q0 -.120069589779605255e2
+#define Q1 0.194809660700889731e2
+#define Q2 -.891110902798312337e1
double
log(double arg)
@@ -36,7 +37,7 @@ log(double arg)
x *= 2;
exp--;
}
- if(x < sqrto2) {
+ if(x < Sqrto2) {
x *= 2;
exp--;
}
@@ -44,19 +45,28 @@ log(double arg)
z = (x-1) / (x+1);
zsq = z*z;
- temp = ((p3*zsq + p2)*zsq + p1)*zsq + p0;
- temp = temp/(((zsq + q2)*zsq + q1)*zsq + q0);
- temp = temp*z + exp*log2;
+ temp = ((P3*zsq + P2)*zsq + P1)*zsq + P0;
+ temp = temp/(((zsq + Q2)*zsq + Q1)*zsq + Q0);
+ temp = temp*z + exp*Log2;
return temp;
}
double
log10(double arg)
{
+ if(arg <= 0) {
+ errno = (arg==0)? ERANGE : EDOM;
+ return -HUGE_VAL;
+ }
+ return log(arg) * Ln10o1;
+}
+double
+log2(double arg)
+{
if(arg <= 0) {
errno = (arg==0)? ERANGE : EDOM;
return -HUGE_VAL;
}
- return log(arg) * ln10o1;
+ return log(arg) * Ln2o1;
}