Changeset 236a4c in git

Timestamp:

Feb 19, 2013, 5:54:48 PM (10 years ago)

Author:

Yue Ren <ren@…>

Branches:

(u'jengelh-datetime', 'ceac47cbc86fe4a15902392bdbb9bd2ae0ea02c6')(u'spielwiese', 'f875bbaccd0831e36aaed09ff6adeb3eb45aeb94')

Children:

3693cf538753b44f4862e8f7dc384771128b32d0

Parents:

a81ddb0212ef423502e854619e48d78e1477e806

git-author:

Yue Ren <ren@mathematik.uni-kl.de>2013-02-19 17:54:48+01:00

git-committer:

Yue Ren <ren@mathematik.uni-kl.de>2013-02-19 20:07:24+01:00

Message:

fix: number to float conversion on 64bit machines

File:

• libpolys/coeffs/shortfl.cc (modified) (6 diffs)

libpolys/coeffs/shortfl.cc

@@ -6 +6 @@
 * ABSTRACT:
 */
-
+#include <iostream>
 #include <coeffs/shortfl.h>
 
 #include <string.h>
+#include <math.h>
 #include <coeffs/coeffs.h>
 #include <coeffs/numbers.h>
@@ -373 +374 @@
   assume( getCoeffType(r) == ID );
   assume( getCoeffType(aRing) ==  n_Zp );
-  
+
   int i = (int)((long)from);
   float f = (float)i;
@@ -389 +390 @@
 
 static number nrMapC(number from, const coeffs aRing, const coeffs r)
-{  
+{
   assume( getCoeffType(r) == ID );
   assume( getCoeffType(aRing) == n_long_C );
-  
+
   gmp_float h = ((gmp_complex*)from)->real();
   float t =(float)mpf_get_d((mpf_srcptr)&h);
@@ -406 +407 @@
 */
 #define SR_HDL(A) ((long)(A))
-#define mpz_isNeg(A) ((A)->_mp_size<0)
-#define mpz_limb_size(A) ((A)->_mp_size)
-#define mpz_limb_d(A) ((A)->_mp_d)
-#define MPZ_DIV(A,B,C) mpz_tdiv_q((A),(B),(C))
 #define IS_INT(A) ((A)->s==3)
 #define IS_IMM(A) (SR_HDL(A)&SR_INT)
 #define GET_NOM(A) ((A)->z)
 #define GET_DENOM(A) ((A)->n)
-#define MPZ_INIT mpz_init
-#define MPZ_CLEAR mpz_clear
 
   assume( getCoeffType(r) == ID );
   assume( getCoeffType(aRing) == n_Q );
-
-  mpz_t h;
-  mpz_ptr g,z,n;
-  int i,j,t,s;
-  float ba,rr,rn,y;
 
   if (IS_IMM(from))
     return nf((float)nlInt(from,NULL /* dummy for nlInt*/)).N();
-  z=GET_NOM(from);
-  s=0X10000;
-  ba=(float)s;
-  ba*=ba;
-  rr=0.0;
-  i=mpz_size1(z);
+
+  /* read out the enumerator */
+  mpz_ptr z=GET_NOM(from);
+  if(mpz_size1(z)>4)
+  {
+    WerrorS("float overflow");
+    return nf(0.0).N();
+  }
+  mpf_t e;
+  mpf_init(e);
+  mpf_set_z(e,z);
+
+  /* if number was an integer, we are done*/
   if(IS_INT(from))
   {
-    if(i>4)
-    {
-      WerrorS("float overflow");
-      return nf(rr).N();
-    }
-    i--;
-    rr=(float)mpz_limb_d(z)[i];
-    while(i>0)
-    {
-      i--;
-      y=(float)mpz_limb_d(z)[i];
-      rr=rr*ba+y;
-    }
-    if(mpz_isNeg(z))
-      rr=-rr;
-    return nf(rr).N();
-  }
-  n=GET_DENOM(from);
-  j=s=mpz_limb_size(n);
-  if(j>i)
-  {
-    g=n; n=z; z=g;
-    t=j; j=i; i=t;
-  }
-  t=i-j;
-  if(t>4)
-  {
-    if(j==s)
-      WerrorS("float overflow");
-    return nf(rr).N();
-  }
-  if(t>1)
-  {
-    g=h;
-    MPZ_INIT(g);
-    MPZ_DIV(g,z,n);
-    t=mpz_size1(g);
-    if(t>4)
-    {
-      MPZ_CLEAR(g);
-      if(j==s)
-        WerrorS("float overflow");
-      return nf(rr).N();
-    }
-    t--;
-    rr=(float)mpz_limb_d(g)[t];
-    while(t)
-    {
-      t--;
-      y=(float)mpz_limb_d(g)[t];
-      rr=rr*ba+y;
-    }
-    MPZ_CLEAR(g);
-    if(j!=s)
-      rr=1.0/rr;
-    if(mpz_isNeg(z))
-      rr=-rr;
-    return nf(rr).N();
-  }
-  rn=(float)mpz_limb_d(n)[j-1];
-  rr=(float)mpz_limb_d(z)[i-1];
-  if(j>1)
-  {
-    rn=rn*ba+(float)mpz_limb_d(n)[j-2];
-    rr=rr*ba+(float)mpz_limb_d(z)[i-2];
-    i--;
-  }
-  if(t!=0)
-    rr=rr*ba+(float)mpz_limb_d(z)[i-2];
-  if(j==s)
-    rr=rr/rn;
-  else
-    rr=rn/rr;
-  if(mpz_isNeg(z))
-    rr=-rr;
-  return nf(rr).N();
-}
-
+    double basis;
+    signed long int exp;
+    basis = mpf_get_d_2exp(&exp, e);
+    float f = ldexp(basis,exp);
+    mpf_clear(e);
+    return nf(f).N();
+  }
+
+  /* else read out the denominator */
+  mpz_ptr n = GET_DENOM(from);
+  if(mpz_size1(n)>4)
+  {
+    WerrorS("float overflow");
+    mpf_clear(e);
+    return nf(0.0).N();
+  }
+  mpf_t d;
+  mpf_init(d);
+  mpf_set_z(d,n);
+
+  /* and compute the quotient */
+  mpf_t q;
+  mpf_init(q);
+  mpf_div(q,e,d);
+
+  double basis;
+  signed long int exp;
+  basis = mpf_get_d_2exp(&exp, q);
+  float f = ldexp(basis,exp);
+  mpf_clear(e);
+  mpf_clear(d);
+  mpf_clear(q);
+  return nf(f).N();
+}
 
 nMapFunc nrSetMap(const coeffs src, const coeffs dst)
@@ -544 +497 @@
 {
   assume( getCoeffType(n) == ID );
-   
+
   assume( p == NULL );
-   
+
   n->cfKillChar = ndKillChar; /* dummy */
   n->ch = 0;
@@ -578 +531 @@
   n->cfDBTest=ndDBTest; // not yet implemented: nrDBTest;
 #endif
-  
+
   n->nCoeffIsEqual = ndCoeffIsEqual;
 
   n->float_len = SHORT_REAL_LENGTH;
   n->float_len2 = SHORT_REAL_LENGTH;
-   
+
   // TODO: Any variables?
   return FALSE;