From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: from mail3-relais-sop.national.inria.fr (mail3-relais-sop.national.inria.fr [192.134.164.104]) by sourceware.org (Postfix) with ESMTPS id A835A3858427 for ; Thu, 10 Feb 2022 16:21:44 +0000 (GMT) DMARC-Filter: OpenDMARC Filter v1.4.1 sourceware.org A835A3858427 Authentication-Results: sourceware.org; dmarc=none (p=none dis=none) header.from=inria.fr Authentication-Results: sourceware.org; spf=pass smtp.mailfrom=inria.fr DKIM-Signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=inria.fr; s=dc; h=date:message-id:from:to:cc:in-reply-to:subject: references; bh=ZhpwvPLAGFL9Rd1A9CNaK4R4SqLUkbC21ewDT7uFZts=; b=oSE20obxHz1IJAMazYT/2SexUwYKDCoqFMm7tqH3pKD5JWclpgb0u8kU +sUunQlRnGLVomyXxaBxmGantzXWPC6PQT6iyEmQJ7ch9cMCG9hKlLNV8 etwQ58I6urTiBObsuJatLKse9yfpaqtj5yfroVBKRVuOgi8r5NnG60QYt s=; X-IronPort-AV: E=Sophos;i="5.88,359,1635199200"; d="scan'208";a="5583824" Received: from tomate.loria.fr (HELO tomate) ([152.81.10.51]) by mail3-relais-sop.national.inria.fr with ESMTP/TLS/ECDHE-RSA-AES256-GCM-SHA384; 10 Feb 2022 17:21:43 +0100 Date: Thu, 10 Feb 2022 17:21:43 +0100 Message-Id: From: Paul Zimmermann To: andoni.arregui@gtd-gmbh.de Cc: newlib@sourceware.org In-Reply-To: <53d624876330f071a9a1588651b5c73db8ca2070.camel@gtd-gmbh.de> (message from Andoni Arregi on Thu, 10 Feb 2022 17:10:35 +0100) Subject: Re: [PATCH 1/2] Improve lgammaf range for very small cases References: <99cdb95d3263da03844fa0343c22d9e4ab87375e.camel@gtd-gmbh.de> <53d624876330f071a9a1588651b5c73db8ca2070.camel@gtd-gmbh.de> X-Spam-Status: No, score=-9.2 required=5.0 tests=BAYES_00, DKIM_SIGNED, DKIM_VALID, DKIM_VALID_AU, DKIM_VALID_EF, GIT_PATCH_0, RCVD_IN_MSPIKE_H3, RCVD_IN_MSPIKE_WL, SPF_HELO_NONE, SPF_PASS, TXREP, T_SCC_BODY_TEXT_LINE autolearn=ham autolearn_force=no version=3.4.4 X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on server2.sourceware.org X-BeenThere: newlib@sourceware.org X-Mailman-Version: 2.1.29 Precedence: list List-Id: Newlib mailing list List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Thu, 10 Feb 2022 16:21:46 -0000 Dear Andoni, I have tested this patch: the maximal error for all binary32 inputs to lgammaf decreases to 7.50e+06 ulps, which is still huge, but less huge than 8.93e+43... Using RedHat newlib MPFR library: 4.1.0 MPFR header: 4.1.0 (based on 4.1.0) Checking function mylgammaf with MPFR_RNDN libm wrong by up to 7.50e+06 ulp(s) [7497618] for x=-0x1.3a7fcap+1 mylgamma gives -0x1p-24 mpfr_mylgamma gives -0x1.e4cf24p-24 Total: errors=509423944 (11.91%) errors2=11684280 maxerr=7.50e+06 ulp(s) Best regards, Paul > From: Andoni Arregi > Date: Thu, 10 Feb 2022 17:10:35 +0100 > Organization: GTD GmbH > User-Agent: Evolution 3.42.3 > > The original cut for small arguments at |x|<2**-70 (copied from the > double version) produces that when computing nadj we get a subnormal > number for t*x and thus, the division of pi/subnormal will be INF and > the logarithm of it too, which is wrong as a result for lgammaf in this > range. > The proposed new limit seems to be safe and has been tested to > produce accurate results. > (Courstesy of Andreas Jung, ESA) > --- > newlib/libm/math/erf_lgamma.c | 2 +- > 1 file changed, 1 insertion(+), 1 deletion(-) > > diff --git a/newlib/libm/math/erf_lgamma.c b/newlib/libm/math/erf_lgamma.c > index f88f63092..84d02159b 100644 > --- a/newlib/libm/math/erf_lgamma.c > +++ b/newlib/libm/math/erf_lgamma.c > @@ -168,7 +168,7 @@ static float zero= 0.0000000000e+00; > *signgamp = -1; > return one/(x-x); > } > - if(ix<0x1c800000) { /* |x|<2**-70, return -log(|x|) */ > + if(ix<0x30800000) { /* |x|<2**-30, return -log(|x|) */ > if(hx<0) { > *signgamp = -1; > return -__ieee754_logf(-x); > -- > 2.35.1 > > >