From d647127d14bd624a461a22609838e79279e69e7d Mon Sep 17 00:00:00 2001 From: Rajan Gill Date: Sun, 16 Sep 2018 15:14:15 -0600 Subject: [PATCH] libs/libc/math: Small change to log() and logf() that improves accuracy and convergence time --- libs/libc/math/lib_log.c | 19 +++++++++++-------- libs/libc/math/lib_logf.c | 19 +++++++++++-------- 2 files changed, 22 insertions(+), 16 deletions(-) diff --git a/libs/libc/math/lib_log.c b/libs/libc/math/lib_log.c index c2d414eb8e..59d7bff62b 100644 --- a/libs/libc/math/lib_log.c +++ b/libs/libc/math/lib_log.c @@ -60,24 +60,24 @@ double log(double x) { double y; double y_old; - double ey; + double ney; double epsilon; int relax_factor; - int i; + int iter; y = 0.0; y_old = 1.0; epsilon = DBL_EPSILON; - i = 0; + iter = 0; relax_factor = 1; while (y > y_old + epsilon || y < y_old - epsilon) { y_old = y; - ey = exp(y); - y -= (ey - x) / ey; + ney = exp(-y); + y -= 1 - x * ney; if (y > 700.0) { @@ -91,13 +91,16 @@ double log(double x) epsilon = (fabs(y) > 1.0) ? fabs(y) * DBL_EPSILON : DBL_EPSILON; - if (++i >= LOGF_MAX_ITER) + if (++iter >= LOGF_MAX_ITER) { relax_factor *= LOGF_RELAX_MULTIPLIER; - i = 0; + iter = 0; } - epsilon *= relax_factor; + if (relax_factor > 1) + { + epsilon *= relax_factor; + } } if (y == 700.0) diff --git a/libs/libc/math/lib_logf.c b/libs/libc/math/lib_logf.c index 0faaf3b55e..1ae9579166 100644 --- a/libs/libc/math/lib_logf.c +++ b/libs/libc/math/lib_logf.c @@ -58,23 +58,23 @@ float logf(float x) { float y; float y_old; - float ey; + float ney; float epsilon; int relax_factor; - int i; + int iter; y = 0.0F; y_old = 1.0F; epsilon = FLT_EPSILON; - i = 0; + iter = 0; relax_factor = 1; while (y > y_old + epsilon || y < y_old - epsilon) { y_old = y; - ey = expf(y); - y -= (ey - x) / ey; + ney = expf(-y); + y -= 1 - x * ney; if (y > FLT_MAX_EXP_X) { @@ -88,13 +88,16 @@ float logf(float x) epsilon = (fabsf(y) > 1.0F) ? fabsf(y) * FLT_EPSILON : FLT_EPSILON; - if (++i >= LOGF_MAX_ITER) + if (++iter >= LOGF_MAX_ITER) { relax_factor *= LOGF_RELAX_MULTIPLIER; - i = 0; + iter = 0; } - epsilon *= relax_factor; + if (relax_factor > 1) + { + epsilon *= relax_factor; + } } if (y == FLT_MAX_EXP_X)