From mboxrd@z Thu Jan 1 00:00:00 1970 X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on inbox.vuxu.org X-Spam-Level: X-Spam-Status: No, score=-3.3 required=5.0 tests=MAILING_LIST_MULTI, RCVD_IN_DNSWL_MED,RCVD_IN_MSPIKE_H3,RCVD_IN_MSPIKE_WL autolearn=ham autolearn_force=no version=3.4.4 Received: (qmail 28202 invoked from network); 10 May 2021 20:53:17 -0000 Received: from mother.openwall.net (195.42.179.200) by inbox.vuxu.org with ESMTPUTF8; 10 May 2021 20:53:17 -0000 Received: (qmail 22358 invoked by uid 550); 10 May 2021 20:53:15 -0000 Mailing-List: contact musl-help@lists.openwall.com; run by ezmlm Precedence: bulk List-Post: List-Help: List-Unsubscribe: List-Subscribe: List-ID: Reply-To: musl@lists.openwall.com Received: (qmail 22325 invoked from network); 10 May 2021 20:53:15 -0000 Date: Mon, 10 May 2021 22:53:03 +0200 From: Szabolcs Nagy To: musl@lists.openwall.com Message-ID: <20210510205303.GO2799122@port70.net> Mail-Followup-To: musl@lists.openwall.com MIME-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline Subject: [musl] [PATCH] math: avoid internal overflow in expm1l.c the ld80 expm1l implementation in c tries to compute 2^k q + 2^k - 1 but 2^k can overflow (when k==16384) while q is small and negative. then the result should be finite, but the code gives nan. to avoid the overflow use 2 (2^(k-1) q + 2^(k-1) - 0.5). note: this implementation is not used on x86 currently. --- src/math/expm1l.c | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/math/expm1l.c b/src/math/expm1l.c index d1715078..79795dc0 100644 --- a/src/math/expm1l.c +++ b/src/math/expm1l.c @@ -109,9 +109,9 @@ long double expm1l(long double x) /* exp(x) = exp(k ln 2) exp(remainder ln 2) = 2^k exp(remainder ln 2). We have qx = exp(remainder ln 2) - 1, so - exp(x) - 1 = 2^k (qx + 1) - 1 = 2^k qx + 2^k - 1. */ - px = scalbnl(1.0, k); - x = px * qx + (px - 1.0); + exp(x) - 1 = 2^k (qx + 1) - 1 = 2 (2^(k-1) qx + 2^(k-1) - 0.5). */ + px = scalbnl(1.0, k-1); + x = 2*(px * qx + (px - 0.5)); return x; } #elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384 -- 2.29.2