From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on inbox.vuxu.org X-Spam-Level: X-Spam-Status: No, score=-3.1 required=5.0 tests=HEADER_FROM_DIFFERENT_DOMAINS, MAILING_LIST_MULTI,RCVD_IN_DNSWL_MED,RCVD_IN_MSPIKE_H4, RCVD_IN_MSPIKE_WL,T_SCC_BODY_TEXT_LINE autolearn=ham autolearn_force=no version=3.4.4 Received: from second.openwall.net (second.openwall.net [193.110.157.125]) by inbox.vuxu.org (Postfix) with SMTP id 4279823DAE for ; Tue, 18 Jun 2024 14:24:10 +0200 (CEST) Received: (qmail 1155 invoked by uid 550); 18 Jun 2024 12:24:06 -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 1123 invoked from network); 18 Jun 2024 12:24:06 -0000 Date: Tue, 18 Jun 2024 14:23:57 +0200 From: Szabolcs Nagy To: Damian McGuckin Cc: MUSL Message-ID: <20240618122357.GL3766212@port70.net> Mail-Followup-To: Damian McGuckin , MUSL References: <31679941-f6c2-64fa-7a8d-6b6f7112c31@esi.com.au> MIME-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline In-Reply-To: <31679941-f6c2-64fa-7a8d-6b6f7112c31@esi.com.au> Subject: Re: [musl] roundf() (and round(), and ...) * Damian McGuckin [2024-06-17 11:48:38 +1000]: > Before I submit any suggestion for a change to roundf/round/.., could I get > a critique of the style to make sure I compy with the style guide. > > Also, I have not used anything except for the now default FLT_EVAL_METHOD > so I need guidance as to where to use 'float_t'. However, the only place > where such arithmetic (here) might be affected is the expression 'a + a' > which is exact anyway so I think it is irrelevant. note that it is unspecified if your algorithm raises inexact or not. (iirc musl asm implementations don't raise inexact, the c code does, c23 now requires no inexact which i guess is what you tried to follow, but that is hard to do in c) otherwise i think your style is fine, but i added some comments. FLT_EVAL_METHOD is not relevant for this code. benchmark data may be useful (or code size e.g. on a soft-float target) because that may be a valid justification to use this implementation. > > Thanks - Damian > > Code follows: > > float roundf(float x) > { > static const int b = 0x7f; > const float a = fabsf(x); > union { float f; uint32_t _f; } r = { a }, _x = { x }; > > if (((int) (r._f >> 23)) < b + 23) /* i.e. effectively |x| < 2^(p-1) */ > { > /* > * capture the sign of the argument > */ > const uint32_t s = _x._f - r._f; i used to use explicit unions then switched to helper function asuint(x) because i found that a bit clearer and compiler optimizes it just fine. and some ppl complained that memcpy is more correct than union, so it is better hidden away behind a helper function if somebody wants to switch. nowadays i'd probably write if (asuint(a) >> 23 < asuint(0x1p23f) >> 23) but e.g. i find it just as clear writing const uint32_t abits = asuint(a); if (abits >> 23 < 0x7f + 23) > /* > * this should achieve rounding to nearest with any > * ties (half-way cases) being rounded away-from-zero. > * (is it wise to use uint32_t instead of int32_t here?) > */ > const uint32_t rf = ((uint32_t) (a + a)) - (uint32_t) a; > > x = (r.f = (float) rf, r._f |= s, r.f); it is isa dependent if int32_t or uint32_t is better, but i'd expect signed convert is more widely supported than unsigned. so i'd write this as const float r = (float)((int32_t)(a + a) - (int32_t)a); x = asfloat(asuint(r) | s); > } > return x; > } > > Note that as per the latest IEEE-754 standard, the above does not raise an > exception in the event of the rounding being inexact. This is not backwards it is unspecified if (uint32_t)a raises inexact exception, so this will be target and compiler dependent in practice. ieee-754-2008: 5.8 Details of conversions from floating-point to integer formats ... A language standard that permits implicit conversions or expressions involving mixed types should require that these be implemented with the inexact-signaling conversion operations below. c23: F.4 Floating to integer conversion ... ... whether conversion of a non-integral floating value raises the "inexact" floating-point exception is unspecified.[note] ... note: In those cases where it matters, library functions can be used to effect such conversions with or without raising the "inexact" floating-point exception. > compatible with the existing MUSL equivalents. While it is another issue > altogether, this latest standard also drops the nearbyint() in favour of > routines called (as per the C standard) 'roundevenf()/roundeven()/...'.