From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: from us-smtp-delivery-124.mimecast.com (us-smtp-delivery-124.mimecast.com [170.10.133.124]) by sourceware.org (Postfix) with ESMTPS id AAC2C38582A1 for ; Sun, 13 Nov 2022 20:39:13 +0000 (GMT) DMARC-Filter: OpenDMARC Filter v1.4.1 sourceware.org AAC2C38582A1 Authentication-Results: sourceware.org; dmarc=pass (p=none dis=none) header.from=redhat.com Authentication-Results: sourceware.org; spf=pass smtp.mailfrom=redhat.com DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=redhat.com; s=mimecast20190719; t=1668371953; h=from:from:reply-to:reply-to:subject:subject:date:date: message-id:message-id:to:to:cc:cc:mime-version:mime-version: content-type:content-type:in-reply-to:in-reply-to: references:references; bh=SLRnCZqLJHlWq7DIWQ1fxtLKtq+9taz5tH6t0q9LHnk=; b=EpoaIyailwb24d8zVfbCJLCuzKpo/GyIyBHYi9xRMp2M9rIQY6G7jn6khivISBI4ChAqwh wXwS2VQ3y7pm0QFMuO2PkSMRQag+b6oKt1UsZ/DDVULnSjO3/khoCNRQlHPDC3QU3xX8CY qycEDQ1m/hXmPfrkIhVNMdMa1d79+0E= Received: from mimecast-mx02.redhat.com (mx3-rdu2.redhat.com [66.187.233.73]) by relay.mimecast.com with ESMTP with STARTTLS (version=TLSv1.2, cipher=TLS_ECDHE_RSA_WITH_AES_256_GCM_SHA384) id us-mta-478-pKCaydkZMY2ZqHrTpSYdxA-1; Sun, 13 Nov 2022 15:39:09 -0500 X-MC-Unique: pKCaydkZMY2ZqHrTpSYdxA-1 Received: from smtp.corp.redhat.com (int-mx10.intmail.prod.int.rdu2.redhat.com [10.11.54.10]) (using TLSv1.2 with cipher AECDH-AES256-SHA (256/256 bits)) (No client certificate requested) by mimecast-mx02.redhat.com (Postfix) with ESMTPS id A478329ABA3F; Sun, 13 Nov 2022 20:39:09 +0000 (UTC) Received: from tucnak.zalov.cz (unknown [10.39.192.38]) by smtp.corp.redhat.com (Postfix) with ESMTPS id 472E7492B08; Sun, 13 Nov 2022 20:39:09 +0000 (UTC) Received: from tucnak.zalov.cz (localhost [127.0.0.1]) by tucnak.zalov.cz (8.17.1/8.17.1) with ESMTPS id 2ADKd4ZF2804762 (version=TLSv1.3 cipher=TLS_AES_256_GCM_SHA384 bits=256 verify=NOT); Sun, 13 Nov 2022 21:39:04 +0100 Received: (from jakub@localhost) by tucnak.zalov.cz (8.17.1/8.17.1/Submit) id 2ADKd3NP2804761; Sun, 13 Nov 2022 21:39:03 +0100 Date: Sun, 13 Nov 2022 21:39:03 +0100 From: Jakub Jelinek To: Aldy Hernandez , "Joseph S. Myers" Cc: GCC patches , Andrew MacLeod Subject: Re: [PATCH] [range-ops] Implement sqrt. Message-ID: Reply-To: Jakub Jelinek References: <20221113200553.440728-1-aldyh@redhat.com> MIME-Version: 1.0 In-Reply-To: <20221113200553.440728-1-aldyh@redhat.com> X-Scanned-By: MIMEDefang 3.1 on 10.11.54.10 X-Mimecast-Spam-Score: 0 X-Mimecast-Originator: redhat.com Content-Type: text/plain; charset=us-ascii Content-Disposition: inline X-Spam-Status: No, score=-3.9 required=5.0 tests=BAYES_00,DKIMWL_WL_HIGH,DKIM_SIGNED,DKIM_VALID,DKIM_VALID_AU,DKIM_VALID_EF,RCVD_IN_DNSWL_NONE,RCVD_IN_MSPIKE_H2,SPF_HELO_NONE,SPF_NONE,TXREP autolearn=ham autolearn_force=no version=3.4.6 X-Spam-Checker-Version: SpamAssassin 3.4.6 (2021-04-09) on server2.sourceware.org List-Id: On Sun, Nov 13, 2022 at 09:05:53PM +0100, Aldy Hernandez wrote: > It seems SQRT is relatively straightforward, and it's something Jakub > wanted for this release. > > Jakub, what do you think? > > p.s. Too tired to think about op1_range. That would be multiplication of the same value twice, i.e. fop_mult with trio that has op1_op2 () == VREL_EQ? But see below, as sqrt won't be always precise, we need to account for some errors. > gcc/ChangeLog: > > * gimple-range-op.cc (class cfn_sqrt): New. > (gimple_range_op_handler::maybe_builtin_call): Add cases for sqrt. Yes, I'd like to see SQRT support in. The only thing I'm worried is that unlike {+,-,*,/}, negation etc. typically implemented in hardware or precise soft-float, sqrt is often implemented in library using multiple floating point arithmetic functions. And different implementations have different accuracy. So, I wonder if we don't need to add a target hook where targets will be able to provide upper bound on error for floating point functions for different floating point modes and some way to signal unknown accuracy/can't be trusted, in which case we would give up or return just the range for VARYING. Then, we could write some tests that say in a loop constructs random floating point values (perhaps sanitized to be non-NAN), calls libm function and the same mpfr one and return maximum error in ulps. And then record those, initially for glibc and most common targets and gradually maintainers could supply more. If we add an infrastructure for that within a few days, then we could start filling the details. One would hope that sqrt has < 10ulps accuracy if not already the 0.5ulp one, but for various other functions I think it can be much more. Oh, nanq ib libquadmath has terrible accuracy, but that one fortunately is not builtin... If we have some small integer for ulps accuracy of calls (we could use 0 for 0.5ulps accuracy aka precise), wonder if we'd handle it just as a loop of doing n times frange_nextafter or something smarter. > --- a/gcc/gimple-range-op.cc > +++ b/gcc/gimple-range-op.cc > @@ -43,6 +43,7 @@ along with GCC; see the file COPYING3. If not see > #include "range.h" > #include "value-query.h" > #include "gimple-range.h" > +#include "fold-const-call.h" > > // Given stmt S, fill VEC, up to VEC_SIZE elements, with relevant ssa-names > // on the statement. For efficiency, it is an error to not pass in enough > @@ -301,6 +302,41 @@ public: > } > } op_cfn_constant_p; > > +// Implement range operator for SQRT. > +class cfn_sqrt : public range_operator_float > +{ > + using range_operator_float::fold_range; > +private: > + REAL_VALUE_TYPE real_sqrt (const REAL_VALUE_TYPE &arg, tree type) const > + { > + tree targ = build_real (type, arg); > + tree res = fold_const_call (as_combined_fn (BUILT_IN_SQRT), type, targ); > + return *TREE_REAL_CST_PTR (res); > + } > + void rv_fold (REAL_VALUE_TYPE &lb, REAL_VALUE_TYPE &ub, bool &maybe_nan, > + tree type, > + const REAL_VALUE_TYPE &lh_lb, > + const REAL_VALUE_TYPE &lh_ub, > + const REAL_VALUE_TYPE &, > + const REAL_VALUE_TYPE &, > + relation_kind) const final override > + { > + if (real_compare (LT_EXPR, &lh_ub, &dconst0)) > + { > + real_nan (&lb, "", 0, TYPE_MODE (type)); > + ub = lb; > + maybe_nan = true; > + return; > + } > + lb = real_sqrt (lh_lb, type); > + ub = real_sqrt (lh_ub, type); > + if (real_compare (GE_EXPR, &lh_lb, &dconst0)) > + maybe_nan = false; > + else > + maybe_nan = true; Doesn't this for say VARYING range result in [NAN, +INF] range? We want [-0.0, +INF]. So perhaps the real_compare should be done before doing the real_sqrt calls and for the maybe_nan case use hardcoded -0.0 as lb? BTW, as for the ulps, another thing to test is whether even when the library has certain number of ulps error worst case whether it still obeys the basic math properties of the function or not. Say for sqrt that it always fits into [-0.0, +INF] (guess because of the flush denormals to zero we wouldn't have a problem here for say 30ulps sqrt that [nextafter (0.0, 1.0) * 16.0, 64.0] wouldn't be considered [-nextafter (0.0, 1.0) * 16.0, 8.0 + 30ulps] but just [-0.0, 8.0 + 30ulps], but later on say sin/cos, which mathematically should have result always in [-1.0, 1.0] +-NAN, it would be interesting to see if there aren't some implementations that would happily return 1.0 + 15ulps or -1.0 - 20ulps. And unrelated thought about reverse y = x * x; if we know y's range - op1_range/op2_range in that case could be handled as sqrt without the library ulps treatment (if we assume that multiplication is always precise), but the question is if op1_range or op2_range is called at all in those cases and whether we could similarly use trio to derive in that case the x's range. Jakub