- documented what function does

- cleaned code
- added header about license...


git-svn-id: svn://scm.gforge.inria.fr/svn/mpfr/trunk@823 280ebfd0-de03-0410-8827-d642c229c3f4

1 files changed, 47 insertions, 45 deletions

diff --git a/extract.c b/extract.c
index 4d70095f3..944922f37 100644
--- a/extract.c
+++ b/extract.c
@@ -1,61 +1,63 @@
-#include <stdio.h>
-#include <math.h>
-#include <stdlib.h> 
-#include <strings.h>
+/* mpfr_extract -- bit-extraction function for the binary splitting algorithm
+
+Copyright (C) 2000 Free Software Foundation.
+
+This file is part of the MPFR Library.
+
+The MPFR Library is free software; you can redistribute it and/or modify
+it under the terms of the GNU Library General Public License as published by
+the Free Software Foundation; either version 2 of the License, or (at your
+option) any later version.
+
+The MPFR Library is distributed in the hope that it will be useful, but
+WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Library General Public
+License for more details.
+
+You should have received a copy of the GNU Library General Public License
+along with the MPFR Library; see the file COPYING.LIB.  If not, write to
+the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
+MA 02111-1307, USA. */
+
 #include "gmp.h"
 #include "gmp-impl.h"
 #include "mpfr.h"
-#include "assert.h"
 
+/* given 0 <= |p| < 1, this function extracts limbs of p and puts them in y.
+   It is mainly designed for the "binary splitting" algorithm together with
+   generic.c.
 
-int
+   More precisely, if B = 2^BITS_PER_MP_LIMB:
+   - for i=0, y = floor(p * B)
+   - for i>0, y = (p * B^(2^i)) mod B^(2^(i-1))
+ */
+void
 #if __STDC__
-mpfr_extract(mpz_ptr y, mpfr_srcptr p, int  i)
+mpfr_extract (mpz_ptr y, mpfr_srcptr p, unsigned int i)
 #else
-mpfr_extract(y,p,i)
-mpz_ptr y; 
-mpfr_srcptr p;
-int  i;
+mpfr_extract (y, p, i)
+     mpz_ptr y; 
+     mpfr_srcptr p;
+     unsigned int i;
 #endif
 {
-  int two_i=1 << i;
+  int two_i = 1 << i;
   int two_i_2 = i ? two_i / 2 : 1;
-  int size;
-  int j;
+  mp_size_t size_p = MPFR_ABSSIZE(p);
 
-  /* TODO: gestion de l'infini. */
+  /* as 0 <= |p| < 1, we don't have to care with infinities, NaN, ... */
   
-  if (MPFR_ABSSIZE(p) < two_i){
-    int j;
-    y->_mp_alloc=two_i_2 ;	
-    y->_mp_size=two_i_2 ;	
-    PTR(y) = (*_mp_allocate_func)(two_i_2 * sizeof(mp_limb_t));
+  _mpz_realloc (y, two_i_2);
+  if (size_p < two_i) {
     MPN_ZERO (PTR(y), two_i_2);
-    assert(y->_mp_d!=NULL);
-    for(j= 0; j < MPFR_ABSSIZE(p) -  two_i_2 ; j++){
-      y->_mp_d[j + two_i - MPFR_ABSSIZE(p)] = p->_mp_d[j];
-    }
-  } else
-    {
-      PTR(y) = (*_mp_allocate_func)(two_i_2  * sizeof(mp_limb_t));
-      memcpy(y -> _mp_d, p->_mp_d+MPFR_ABSSIZE(p) - two_i, two_i_2 * sizeof(mp_limb_t));
-      y->_mp_alloc=two_i_2 ;	
-      y->_mp_size=two_i_2 ;
-    }
-  
-  size = MPFR_ABSSIZE(y);
-  for (j = 0; j < size; j++)
-    {
-      if (y->_mp_d[j])
-         { 
-	   if MPFR_ISNEG(p)
-		     mpz_neg(y,y);
-            return 0;
-         }
-    }
-  y->_mp_size=0;	
-  
-  return 0;
+    if (size_p >= two_i_2)
+      MPN_COPY (PTR(y) + two_i - size_p, MPFR_MANT(p), size_p - two_i_2);
+  }
+  else
+    MPN_COPY (PTR(y), MPFR_MANT(p) + size_p - two_i, two_i_2);
+
+  MPN_NORMALIZE (PTR(y), two_i_2);
+  SIZ(y) = (MPFR_ISNEG(p)) ? -two_i_2 : two_i_2;
 }