From nobody Sun Jun 7 21:00:58 2026 X-Original-To: dev-commits-src-main@mlmmj.nyi.freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2610:1c1:1:606c::19:1]) by mlmmj.nyi.freebsd.org (Postfix) with ESMTP id 4gYSJC0zyGz6g5CL for ; Sun, 07 Jun 2026 21:00:59 +0000 (UTC) (envelope-from git@FreeBSD.org) Received: from mxrelay.nyi.freebsd.org (mxrelay.nyi.freebsd.org [IPv6:2610:1c1:1:606c::19:3]) (using TLSv1.3 with cipher TLS_AES_256_GCM_SHA384 (256/256 bits) key-exchange X25519 server-signature RSA-PSS (4096 bits) server-digest SHA256 client-signature RSA-PSS (4096 bits) client-digest SHA256) (Client CN "mxrelay.nyi.freebsd.org", Issuer "R13" (not verified)) by mx1.freebsd.org (Postfix) with ESMTPS id 4gYSJB4j5yz3SmB for ; Sun, 07 Jun 2026 21:00:58 +0000 (UTC) (envelope-from git@FreeBSD.org) DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=freebsd.org; s=dkim; t=1780866058; h=from:from:reply-to:subject:subject:date:date:message-id:message-id: to:to:cc:cc:mime-version:mime-version:content-type:content-type: content-transfer-encoding:content-transfer-encoding; bh=jMBWLqfrsmdp/y+RXcKA5pdRGrsmPP4qm9JzgS8y5w0=; b=lFUr/pVPYTppr/M+9RyX2iMzhAeylHYCnWD7WQTdYy1PF3rsYVIcdyf4RrelPWJYaMSW8f C7h6E8YcaR6lptq7gozAcJUJf+yGkV8liWw6rlts6AFUONNFvwmfCMqM1b4mbULHeQoDtK 7JrIRtUX5z+7CwsPZ2ln7t2Ji5t3SR1U7PdYthRW4bRy8oyPTfViZrmb5z279ExNcOm+4s YDNM8WbRsooiBgrihwHKrmY/Eeon9dxkrpWLoBV4TvNY8I0vjQD5vB5WEjUzVHXiX8qjq4 SoqzI/8lAJpuVqYE9zShRbokoUL32ARnyW+q/o3MuLFc6OtCcGTwNZnmP8MhXQ== ARC-Seal: i=1; s=dkim; d=freebsd.org; t=1780866058; a=rsa-sha256; cv=none; b=Bl+kFmsYapdGDlPA0tz4/KelqVs5Rt0EvjuVNEtACj+c0bhq49Q5oGwrpE76wftsR9j/6K 3UWmpFZQJXmOM9nJjq9PtpgmWOsIQdX8uxJr3/ElCijnMMRtRRv6AqxuUhPxIeB1YocYKO M5FIhLLgMusNkz/96RhJpivKrwJDcy1vNLEZwiljonegtyd9j2ciLIWkFBpCqATC2K6OR1 HMv1+5RFfy0KQdR65x+w9xUM5RY93itjtisF/h+fjw4o7e3Ez1st4NRfKZE9jXwVv5GO1+ 8uqZAQblYhxlpghto16ZnSgNsm3cgIy7Huc1S8/0FG7soBZYomjVV/fSI9dEMQ== ARC-Authentication-Results: i=1; mx1.freebsd.org; none ARC-Message-Signature: i=1; a=rsa-sha256; c=relaxed/relaxed; d=freebsd.org; s=dkim; t=1780866058; h=from:from:reply-to:subject:subject:date:date:message-id:message-id: to:to:cc:cc:mime-version:mime-version:content-type:content-type: content-transfer-encoding:content-transfer-encoding; bh=jMBWLqfrsmdp/y+RXcKA5pdRGrsmPP4qm9JzgS8y5w0=; b=wYb1TGIx2IAi4FrxaWv6lBYsymuTwCjPPB4g7IaPsDtpYWw9g4LdY2EExTvCbLEpagsor8 nVD4iEwyU9/G9pmpplIvFplW8jYOnAptMROYz0fBKXG5dw3GPIAIkXiTFnjRDS2bY/zdQ/ t249TP1GGpm+hG2eWSeWGxMrtvUPqb0IeDRjmJXn63+p3bu3pfNO4A8jO7WMxzlGRiggy7 E+R08AJw1sAZEQIIcmjgrbtmZvOZau+WXlG83Pqsl6E9ldJAkzl3vvgGUHjwGrVb4f6/a7 RP9fKuCrDz9SOjkSSDDc6im7P/YsrMmfahUmiFqDq4PROgawJaU/2WcLcWV38Q== Received: from gitrepo.freebsd.org (gitrepo.freebsd.org [IPv6:2610:1c1:1:6068::e6a:5]) by mxrelay.nyi.freebsd.org (Postfix) with ESMTP id 4gYSJB40tbz3Ck for ; Sun, 07 Jun 2026 21:00:58 +0000 (UTC) (envelope-from git@FreeBSD.org) Received: from git (uid 1279) (envelope-from git@FreeBSD.org) id 47a4d by gitrepo.freebsd.org (DragonFly Mail Agent v0.13+ on gitrepo.freebsd.org); Sun, 07 Jun 2026 21:00:58 +0000 To: src-committers@FreeBSD.org, dev-commits-src-all@FreeBSD.org, dev-commits-src-main@FreeBSD.org Cc: Steve Kargl From: Robert Clausecker Subject: git: c3f6dcea1992 - main - msun: Fix up for recent rsqrt[fl] functions List-Id: Commit messages for the main branch of the src repository List-Archive: https://lists.freebsd.org/archives/dev-commits-src-main List-Help: List-Post: List-Subscribe: List-Unsubscribe: X-BeenThere: dev-commits-src-main@freebsd.org Sender: owner-dev-commits-src-main@FreeBSD.org List-Id: List-Post: List-Help: List-Subscribe: List-Unsubscribe: List-Owner: Precedence: list MIME-Version: 1.0 Content-Type: text/plain; charset=utf-8 Content-Transfer-Encoding: 8bit X-Git-Committer: fuz X-Git-Repository: src X-Git-Refname: refs/heads/main X-Git-Reftype: branch X-Git-Commit: c3f6dcea199289329c1d3b91b69e5a4821fc3dff Auto-Submitted: auto-generated Date: Sun, 07 Jun 2026 21:00:58 +0000 Message-Id: <6a25dc0a.47a4d.198123cf@gitrepo.freebsd.org> The branch main has been updated by fuz: URL: https://cgit.FreeBSD.org/src/commit/?id=c3f6dcea199289329c1d3b91b69e5a4821fc3dff commit c3f6dcea199289329c1d3b91b69e5a4821fc3dff Author: Steve Kargl AuthorDate: 2026-06-07 19:12:16 +0000 Commit: Robert Clausecker CommitDate: 2026-06-07 20:59:19 +0000 msun: Fix up for recent rsqrt[fl] functions Paul Zimmermann (of Core-Math and MPFR fame) graciously tested the recently committed rsqrt[fl]() functions. He identified 127 incorrectly rounded values for rsqrtf() in round-to-nearest mode. This patch fixes the rounding in RN. Exhaustive testing now shows that rsqrtf() is corrected rounded for RN. He also tested rsqrt() and rsqrtl() in the interval [1,4). Both appear to be correctly rounded. Finally, the patch includes small changes to comments. A concise list of changes is * lib/msun/src/s_rsqrt.c: . Fix comments. * lib/msun/src/s_rsqrtf.c . Fix comments. . Exhaustive testing by Paul Zimmermann found 127 incorrectly rounded values in round-to-nearests. These gave have the form 0x1.13e07pN with N an odd integer. With this patch, all values are now correctly rounded in round-to-nearest. * lib/msun/src/s_rsqrtl.c . Fix comments. . Move all variable declarations to top of function and sort. PR: 295706 MFC after: 1 week Fixes: 3085fc9d97bd83785ba3ba43e0378d7d67987d1f --- lib/msun/src/s_rsqrt.c | 4 ++-- lib/msun/src/s_rsqrtf.c | 37 ++++++++++++++++++++++--------------- lib/msun/src/s_rsqrtl.c | 8 +++----- 3 files changed, 27 insertions(+), 22 deletions(-) diff --git a/lib/msun/src/s_rsqrt.c b/lib/msun/src/s_rsqrt.c index 0a513afe8ed2..cac16f837a05 100644 --- a/lib/msun/src/s_rsqrt.c +++ b/lib/msun/src/s_rsqrt.c @@ -108,14 +108,14 @@ rsqrt(double x) hx = (hx & 0x000fffff) | 0x3fe00000; SET_HIGH_WORD(x, hx); - /* m is odd. Put x into [0.25,5) and increase m. */ + /* m is odd. Put x into [0.25,0.5) and increase m. */ if (m & 1) { x /= 2; m += 1; } m = -(m >> 1); /* Prepare for 2^(-m/2). */ - y = 1 / sqrt(x); /* ~52-bit estimate. */ + y = 1 / sqrt(x); /* ~52-bit estimate. */ h = y / 2; diff --git a/lib/msun/src/s_rsqrtf.c b/lib/msun/src/s_rsqrtf.c index b71f7baf5657..faeb1c300e18 100644 --- a/lib/msun/src/s_rsqrtf.c +++ b/lib/msun/src/s_rsqrtf.c @@ -108,7 +108,7 @@ rsqrtf(float x) ix = (ix & 0x007fffff) | 0x3f000000; SET_FLOAT_WORD(x, ix); /* x is in [0.5,1). */ - /* m is odd. Put x into [0.25,5) and increase m. */ + /* m is odd. Put x into [0.25,0.5) and increase m. */ if (m & 1) { x /= 2; m += 1; @@ -122,33 +122,40 @@ rsqrtf(float x) * a max ulp of ~1.49. */ y = 1 / sqrtf(x); - h = y / 2; /* - * For values of x with a representation of 0x1.fffffcpN with - * N an odd integer, the computed rsqrtf() is not correctly - * rounded in round-to-nearest without toggling the rounding - * mode to FE_TOWARDZERO. Note, FE_DOWNWARD also works. - * However, messing with the rounding mode is expensive, so - * only do it when necessary. Example, x = 0x1.fffffcp3 gives - * y --> 0x3f800001. + * For values of x with representations of the forms 0x1.fffffcpN + * or 0x1.13e07pN with N an odd integer, the computed rsqrtf() is + * not correctly rounded in round-to-nearest without fiddling with + * the rounding mode and/or bits of x. */ GET_FLOAT_WORD(ix, y); - if ((ix & 0x000fffff) == 1) { + if ((ix & 0x000fffff) == 1) { /* 0x1.fffffcpN */ + rnd = fegetround(); + fesetround(FE_DOWNWARD); + _MUL(x, y, zh, zl); + _XMUL(zh, zl, h, 0, ph, pl); + fesetround(rnd); + _XADD(-ph, -pl, half, 0, rh, rl); + y = h + h * rh; + } else if((ix & 0x00ffffff) == 0x00ae6055) { /* 0x1.13e07pN */ + GET_FLOAT_WORD(ix, x); + SET_FLOAT_WORD(x, ix - 1); rnd = fegetround(); - fesetround(FE_TOWARDZERO); + fesetround(FE_DOWNWARD); _MUL(x, y, zh, zl); _XMUL(zh, zl, h, 0, ph, pl); - fesetround(rnd); + _XADD(-ph, -pl, half, 0, rh, rl); + y = h + h * rh; + fesetround(rnd); } else { _MUL(x, y, zh, zl); _XMUL(zh, zl, h, 0, ph, pl); + _XADD(-ph, -pl, half, 0, rh, rl); + y = h * rh + h; } - _XADD(-ph, -pl, half, 0, rh, rl); - y = h * rh + h; - ix = (uint32_t)(m + 128) << 23; SET_FLOAT_WORD(x, ix); return (y *= x); diff --git a/lib/msun/src/s_rsqrtl.c b/lib/msun/src/s_rsqrtl.c index 34d5cdc48173..bf3d32481e3f 100644 --- a/lib/msun/src/s_rsqrtl.c +++ b/lib/msun/src/s_rsqrtl.c @@ -108,7 +108,7 @@ rsqrtl(long double x) } u.bits.exp = 0x3ffe; - /* m is odd. Put x into [0.25,5) and increase m. */ + /* m is odd. Put x into [0.25,0.5) and increase m. */ if (m & 1) { u.e /= 2; m += 1; @@ -141,8 +141,9 @@ long double rsqrtl(long double x) { volatile static const double vzero = 0; + static const double half = 0.5; int hx, m, rnd; - long double y; + long double h, ph, pl, rh, rl, y, zh, zl; /* x = +-0. Raise exception. */ if (x == 0) @@ -183,9 +184,6 @@ rsqrtl(long double x) y = 1 / sqrt((double)x); /* ~52-bit estimate. */ y -= y * (x * y * y - 1) / 2; /* ~104-bit estimate. */ - static const double half = 0.5; - long double h, ph, pl, rh, rl, zh, zl; - h = y / 2; rnd = fegetround();