You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@nuttx.apache.org by xi...@apache.org on 2022/01/01 12:38:43 UTC
[incubator-nuttx] branch master updated: libc/math: fix log and logf calculations on ARMv7 (and maybe others)
This is an automated email from the ASF dual-hosted git repository.
xiaoxiang pushed a commit to branch master
in repository https://gitbox.apache.org/repos/asf/incubator-nuttx.git
The following commit(s) were added to refs/heads/master by this push:
new b859f27 libc/math: fix log and logf calculations on ARMv7 (and maybe others)
b859f27 is described below
commit b859f2750bbc46d30188af73fefaa8c9f6143b0e
Author: Petro Karashchenko <pe...@gmail.com>
AuthorDate: Thu Dec 30 19:08:47 2021 +0200
libc/math: fix log and logf calculations on ARMv7 (and maybe others)
Probably this is a bug of a GCC, but on AMRv7 the code "if (relax_factor > 1)"
generates "bne.n" instruction if "relax_factor" is "int". Few lines above
"relax_factor *= LOG_RELAX_MULTIPLIER;" is done without overflow check hence
at some moment overflow occurs and "relax_factor" becomes a zero and condition
"if (relax_factor > 1)" becomes always evaluated to "true" hence "epsilon"
becomes zero always.
Probably this is not the best way to fix the bug (The best way is to report it
to GCC), but this change allows to get correct behavior of "log" and "logf" for
ARMv7 based MCUs
Signed-off-by: Petro Karashchenko <pe...@gmail.com>
---
libs/libc/math/lib_log.c | 32 +++++++++++++++++---------------
libs/libc/math/lib_logf.c | 8 ++++----
libs/libc/math/lib_logl.c | 24 ++++++++++++++++--------
3 files changed, 37 insertions(+), 27 deletions(-)
diff --git a/libs/libc/math/lib_log.c b/libs/libc/math/lib_log.c
index 638ed62..24edee3 100644
--- a/libs/libc/math/lib_log.c
+++ b/libs/libc/math/lib_log.c
@@ -39,13 +39,15 @@
* Pre-processor Definitions
****************************************************************************/
-/* To avoid looping forever in particular corner cases, every LOGF_MAX_ITER
- * the error criteria is relaxed by a factor LOGF_RELAX_MULTIPLIER.
+#define DBL_MAX_EXP_X 700.0
+
+/* To avoid looping forever in particular corner cases, every LOG_MAX_ITER
+ * the error criteria is relaxed by a factor LOG_RELAX_MULTIPLIER.
* todo: might need to adjust the double floating point version too.
*/
-#define LOGF_MAX_ITER 10
-#define LOGF_RELAX_MULTIPLIER 2
+#define LOG_MAX_ITER 10
+#define LOG_RELAX_MULTIPLIER 2.0
/****************************************************************************
* Public Functions
@@ -62,7 +64,7 @@ double log(double x)
double y_old;
double ey;
double epsilon;
- int relax_factor;
+ double relax_factor;
int iter;
y = 0.0;
@@ -70,7 +72,7 @@ double log(double x)
epsilon = DBL_EPSILON;
iter = 0;
- relax_factor = 1;
+ relax_factor = 1.0;
while (y > y_old + epsilon || y < y_old - epsilon)
{
@@ -78,36 +80,36 @@ double log(double x)
ey = exp(y);
y -= (ey - x) / ey;
- if (y > 700.0)
+ if (y > DBL_MAX_EXP_X)
{
- y = 700.0;
+ y = DBL_MAX_EXP_X;
}
- if (y < -700.0)
+ if (y < -DBL_MAX_EXP_X)
{
- y = -700.0;
+ y = -DBL_MAX_EXP_X;
}
epsilon = (fabs(y) > 1.0) ? fabs(y) * DBL_EPSILON : DBL_EPSILON;
- if (++iter >= LOGF_MAX_ITER)
+ if (++iter >= LOG_MAX_ITER)
{
- relax_factor *= LOGF_RELAX_MULTIPLIER;
+ relax_factor *= LOG_RELAX_MULTIPLIER;
iter = 0;
}
- if (relax_factor > 1)
+ if (relax_factor > 1.0)
{
epsilon *= relax_factor;
}
}
- if (y == 700.0)
+ if (y == DBL_MAX_EXP_X)
{
return INFINITY;
}
- if (y == -700.0)
+ if (y == -DBL_MAX_EXP_X)
{
return INFINITY;
}
diff --git a/libs/libc/math/lib_logf.c b/libs/libc/math/lib_logf.c
index b53805d..efaac61 100644
--- a/libs/libc/math/lib_logf.c
+++ b/libs/libc/math/lib_logf.c
@@ -44,7 +44,7 @@
*/
#define LOGF_MAX_ITER 10
-#define LOGF_RELAX_MULTIPLIER 2
+#define LOGF_RELAX_MULTIPLIER 2.0F
/****************************************************************************
* Public Functions
@@ -60,7 +60,7 @@ float logf(float x)
float y_old;
float ey;
float epsilon;
- int relax_factor;
+ float relax_factor;
int iter;
y = 0.0F;
@@ -68,7 +68,7 @@ float logf(float x)
epsilon = FLT_EPSILON;
iter = 0;
- relax_factor = 1;
+ relax_factor = 1.0F;
while (y > y_old + epsilon || y < y_old - epsilon)
{
@@ -94,7 +94,7 @@ float logf(float x)
iter = 0;
}
- if (relax_factor > 1)
+ if (relax_factor > 1.0F)
{
epsilon *= relax_factor;
}
diff --git a/libs/libc/math/lib_logl.c b/libs/libc/math/lib_logl.c
index 2bab0b2..0250a2a 100644
--- a/libs/libc/math/lib_logl.c
+++ b/libs/libc/math/lib_logl.c
@@ -36,6 +36,12 @@
#include <float.h>
/****************************************************************************
+ * Pre-processor Definitions
+ ****************************************************************************/
+
+#define LDBL_MAX_EXP_X 700.0
+
+/****************************************************************************
* Public Functions
****************************************************************************/
@@ -47,8 +53,8 @@ long double logl(long double x)
long double ey;
long double epsilon;
- y = 0.0;
- y_old = 1.0;
+ y = 0.0;
+ y_old = 1.0;
epsilon = LDBL_EPSILON;
while (y > y_old + epsilon || y < y_old - epsilon)
@@ -57,23 +63,25 @@ long double logl(long double x)
ey = expl(y);
y -= (ey - x) / ey;
- if (y > 700.0)
+ if (y > LDBL_MAX_EXP_X)
{
- y = 700.0;
+ y = LDBL_MAX_EXP_X;
}
- if (y < -700.0)
+ if (y < -LDBL_MAX_EXP_X)
{
- y = -700.0;
+ y = -LDBL_MAX_EXP_X;
}
+
+ epsilon = (fabsl(y) > 1.0) ? fabsl(y) * LDBL_EPSILON : LDBL_EPSILON;
}
- if (y == 700.0)
+ if (y == LDBL_MAX_EXP_X)
{
return INFINITY;
}
- if (y == -700.0)
+ if (y == -LDBL_MAX_EXP_X)
{
return INFINITY;
}