From mboxrd@z Thu Jan 1 00:00:00 1970 Path: main.gmane.org!not-for-mail From: Kevin Ryde Newsgroups: gmane.lisp.guile.devel Subject: inexact->exact bignum rounding Date: Thu, 02 Oct 2003 10:40:53 +1000 Sender: guile-devel-bounces+guile-devel=m.gmane.org@gnu.org Message-ID: <87brt04gkq.fsf@zip.com.au> NNTP-Posting-Host: deer.gmane.org Mime-Version: 1.0 Content-Type: text/plain; charset=us-ascii X-Trace: sea.gmane.org 1065061215 15790 80.91.224.253 (2 Oct 2003 02:20:15 GMT) X-Complaints-To: usenet@sea.gmane.org NNTP-Posting-Date: Thu, 2 Oct 2003 02:20:15 +0000 (UTC) Original-X-From: guile-devel-bounces+guile-devel=m.gmane.org@gnu.org Thu Oct 02 04:20:13 2003 Return-path: Original-Received: from monty-python.gnu.org ([199.232.76.173]) by deer.gmane.org with esmtp (Exim 3.35 #1 (Debian)) id 1A4t4f-0001Wp-00 for ; Thu, 02 Oct 2003 04:20:13 +0200 Original-Received: from localhost ([127.0.0.1] helo=monty-python.gnu.org) by monty-python.gnu.org with esmtp (Exim 4.24) id 1A4t10-0008Tf-MC for guile-devel@m.gmane.org; Wed, 01 Oct 2003 22:16:26 -0400 Original-Received: from list by monty-python.gnu.org with tmda-scanned (Exim 4.24) id 1A4sX6-0001EK-K4 for guile-devel@gnu.org; Wed, 01 Oct 2003 21:45:32 -0400 Original-Received: from mail by monty-python.gnu.org with spam-scanned (Exim 4.24) id 1A4sWO-00014Y-JW for guile-devel@gnu.org; Wed, 01 Oct 2003 21:45:19 -0400 Original-Received: from [61.8.0.36] (helo=snoopy.pacific.net.au) by monty-python.gnu.org with esmtp (Exim 4.24) id 1A4rWl-0006Wm-HE for guile-devel@gnu.org; Wed, 01 Oct 2003 20:41:08 -0400 Original-Received: from mongrel.pacific.net.au (mongrel.pacific.net.au [61.8.0.107]) by snoopy.pacific.net.au (8.12.3/8.12.3/Debian-6.6) with ESMTP id h920exYH008633 for ; Thu, 2 Oct 2003 10:40:59 +1000 Original-Received: from localhost (ppp97.dyn228.pacific.net.au [203.143.228.97]) by mongrel.pacific.net.au (8.12.3/8.12.3/Debian-6.6) with ESMTP id h920bust016664 for ; Thu, 2 Oct 2003 10:37:56 +1000 Original-Received: from gg by localhost with local (Exim 3.35 #1 (Debian)) id 1A4rWY-0002Au-00; Thu, 02 Oct 2003 10:40:54 +1000 Original-To: guile-devel@gnu.org Mail-Copies-To: never User-Agent: Gnus/5.1003 (Gnus v5.10.3) Emacs/21.3 (gnu/linux) X-BeenThere: guile-devel@gnu.org X-Mailman-Version: 2.1.2 Precedence: list List-Id: Developers list for Guile, the GNU extensibility library List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , Errors-To: guile-devel-bounces+guile-devel=m.gmane.org@gnu.org Xref: main.gmane.org gmane.lisp.guile.devel:2834 X-Report-Spam: http://spam.gmane.org/gmane.lisp.guile.devel:2834 I'd like to propose a new style scm_i_big2dbl, below (not really tested yet). This is of course used by exact->inexact, and also a few other places where it makes sense to do the same as exact->inexact (in particular make-rectangular for instance). /* scm_i_big2dbl() rounds to the closest representable double, in accordance with R5RS exact->inexact. In GMP prior to 4.2, mpz_get_d rounding was unspecified. It followed the hardware rounding mode, but on the absolute value actually. We get a 4.2 style truncate to zero by applying mpz_get_d to just the high DBL_MANT_DIG bits (no rounding). This is not fast, but it's only for old GMP. In GMP 4.2 and up, mpz_get_d rounds towards zero. It basically takes the high DBL_MANT_DIG bits, and we adjust to our desired round-to-nearest by looking at the next bit below those. Note that bignums exactly half way between representable doubles are rounded to the next higher absolute value. This seems like an adequate interpretation of "numerically closest", and it's easier to do than a full "nearest-even" style. */ double scm_i_big2dbl (SCM b) { double result; size_t bits; #if __GNU_MP_VERSION < 4 \ || (__GNU_MP_VERSION == 4 && __GNU_MP_VERSION_MINOR < 2) /* GMP prior to 4.2 */ mpz_t tmp; mpz_init (tmp); mpz_tdiv_q_2exp (tmp, SCM_I_BIG_MPZ (b), bits-DBL_MANT_DIG); result = ldexp (mpz_get_d (tmp), bits-DBL_MANT_DIG); mpz_clear (tmp); #else /* GMP 4.2 and up */ result = mpz_get_d (SCM_I_BIG_MPZ (b)); #endif bits = mpz_sizeinbase (SCM_I_BIG_MPZ (b), 2); if (bits > DBL_MANT_DIG) if (mpz_tstbit (SCM_I_BIG_MPZ (b), bits-DBL_MANT_DIG-1)) result += ldexp ((double) mpz_sgn (SCM_I_BIG_MPZ (b)), bits-DBL_MANT_DIG-1); scm_remember_upto_here_1 (b); return result; } _______________________________________________ Guile-devel mailing list Guile-devel@gnu.org http://mail.gnu.org/mailman/listinfo/guile-devel