- bugfix: 2^(3/2) does not throw any more. This was now secured in...

- ...a new function paranoia_check10.

diff --git a/check/paranoia_check.cpp b/check/paranoia_check.cpp
index b7552dcdf8508fd7707f5710670f5367cc77889b..c5f3557337c003ae636ee4a132a72c673cf1ec76 100644 (file)
--- a/check/paranoia_check.cpp
+++ b/check/paranoia_check.cpp
@@ -218,7 +218,7 @@ static unsigned paranoia_check8(void)
     ex f = e.normal();
 
     // The bug caused a division by zero in normal(), so the following
     ex f = e.normal();
 
     // The bug caused a division by zero in normal(), so the following

-    // check is actually bogus...
+    // check is actually quite bogus...

     if (!f.is_equal(e)) {
         clog << "normal(-x/(x+1)) returns " << f << " instead of -x/(x+1)\n";
         ++result;
     if (!f.is_equal(e)) {
         clog << "normal(-x/(x+1)) returns " << f << " instead of -x/(x+1)\n";
         ++result;


@@ -229,7 +229,7 @@ static unsigned paranoia_check8(void)
 // This one was a result of a modification to frac_cancel() & Co. to avoid
 // expanding the numerator and denominator when bringing them from Q[X] to
 // Z[X]. multiply_lcm() forgot to multiply the x-linear term with the LCM of
 // This one was a result of a modification to frac_cancel() & Co. to avoid
 // expanding the numerator and denominator when bringing them from Q[X] to
 // Z[X]. multiply_lcm() forgot to multiply the x-linear term with the LCM of

-// the coefficient's denominators (2 in this case). Introduced on Jan 25th
+// the coefficient's denominators (2 in this case).  Introduced on Jan 25th

 // 2000 and fixed on Jan 31th.
 static unsigned paranoia_check9(void)
 {
 // 2000 and fixed on Jan 31th.
 static unsigned paranoia_check9(void)
 {


@@ -246,6 +246,31 @@ static unsigned paranoia_check9(void)
     return result;
 }
 
     return result;
 }
 

+// I have no idea when this broke.  It has been working long ago, before 0.4.0
+// and on Feb 13th 2000 I found out that things like 2^(3/2) throw an exception
+// "power::eval(): pow(0,0) is undefined" instead of simplifying to 2*2^(1/2).
+// It was fixed that same day.
+static unsigned paranoia_check10(void)
+{
+    unsigned result = 0;
+    
+    ex b = numeric(2);
+    ex e = numeric(3,2);
+    ex r;
+    
+    try {
+        r = pow(b,e).eval();
+        if (!(r-2*sqrt(ex(2))).is_zero()) {
+            clog << "2^(3/2) erroneously returned " << r << " instead of 2*sqrt(2)" << endl;
+            ++result;
+        }
+    } catch (const exception &e) {
+        clog << "2^(3/2) throws " << e.what() << endl;
+        ++result;
+    }
+    return result;
+}
+

 unsigned paranoia_check(void)
 {
     unsigned result = 0;
 unsigned paranoia_check(void)
 {
     unsigned result = 0;


@@ -262,6 +287,7 @@ unsigned paranoia_check(void)
     result += paranoia_check7();
     result += paranoia_check8();
     result += paranoia_check9();
     result += paranoia_check7();
     result += paranoia_check8();
     result += paranoia_check9();

+    result += paranoia_check10();

 
     if (!result) {
         cout << " passed ";
 
     if (!result) {
         cout << " passed ";

diff --git a/ginac/expairseq.cpp b/ginac/expairseq.cpp
index fe910d0316c6c39502d460929deb090f0c0c74b7..7da18758c2ee3c2e14b11e30462a114e3d65b80f 100644 (file)
--- a/ginac/expairseq.cpp
+++ b/ginac/expairseq.cpp
@@ -953,7 +953,6 @@ void expairseq::make_flat(const exvector & v)
         (*cit).printraw(cout);
     }
     cout << endl;
         (*cit).printraw(cout);
     }
     cout << endl;

-    cout.flush();

     */
 }
 
     */
 }
 


@@ -965,11 +964,12 @@ void expairseq::make_flat(const epvector & v)
     // and their cumulative number of operands
     int nexpairseqs=0;
     int noperands=0;
     // and their cumulative number of operands
     int nexpairseqs=0;
     int noperands=0;

-    cit=v.begin();
+
+    cit = v.begin();

     while (cit!=citend) {
         if (cit->rest.bp->tinfo()==tinfo()) {
             nexpairseqs++;
     while (cit!=citend) {
         if (cit->rest.bp->tinfo()==tinfo()) {
             nexpairseqs++;

-            noperands+=ex_to_expairseq((*cit).rest).seq.size();
+            noperands += ex_to_expairseq((*cit).rest).seq.size();

         }
         ++cit;
     }
         }
         ++cit;
     }


@@ -978,7 +978,7 @@ void expairseq::make_flat(const epvector & v)
     seq.reserve(v.size()+noperands-nexpairseqs);
 
     // copy elements and split off numerical part
     seq.reserve(v.size()+noperands-nexpairseqs);
 
     // copy elements and split off numerical part

-    cit=v.begin();
+    cit = v.begin();

     while (cit!=citend) {
         if ((cit->rest.bp->tinfo()==tinfo())&&can_make_flat(*cit)) {
             const expairseq & subseqref=ex_to_expairseq((*cit).rest);
     while (cit!=citend) {
         if ((cit->rest.bp->tinfo()==tinfo())&&can_make_flat(*cit)) {
             const expairseq & subseqref=ex_to_expairseq((*cit).rest);

diff --git a/ginac/expairseq.h b/ginac/expairseq.h
index 9c32dbdf5fb8ec7dd446c10fc146fe617500c9a7..bcd2f66fa8896c7f322eece6c3172a8379cf98d9 100644 (file)
--- a/ginac/expairseq.h
+++ b/ginac/expairseq.h
@@ -165,7 +165,8 @@ protected:
                             vector<bool> & touched,
                             unsigned & number_of_zeroes);
     bool has_coeff_0(void) const;
                             vector<bool> & touched,
                             unsigned & number_of_zeroes);
     bool has_coeff_0(void) const;

-    void add_numerics_to_hashtab(epvector::iterator first_numeric,                                              epvector::const_iterator last_non_zero);
+    void add_numerics_to_hashtab(epvector::iterator first_numeric,
+                                 epvector::const_iterator last_non_zero);

 #endif // def EXPAIRSEQ_USE_HASHTAB
     bool is_canonical() const;
     epvector * expandchildren(unsigned options) const;
 #endif // def EXPAIRSEQ_USE_HASHTAB
     bool is_canonical() const;
     epvector * expandchildren(unsigned options) const;

diff --git a/ginac/mul.cpp b/ginac/mul.cpp
index 6a32fb58033ebb604efaa5d49a9efcaa1e0eea03..054b2131637f01051a3b70c4d0b4abe133d63f5e 100644 (file)
--- a/ginac/mul.cpp
+++ b/ginac/mul.cpp
@@ -92,7 +92,7 @@ mul::mul(const ex & lh, const ex & rh)
 {
     debugmsg("mul constructor from ex,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
 {
     debugmsg("mul constructor from ex,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;

-    overall_coeff=_ex1();
+    overall_coeff = _ex1();

     construct_from_2_ex(lh,rh);
     GINAC_ASSERT(is_canonical());
 }
     construct_from_2_ex(lh,rh);
     GINAC_ASSERT(is_canonical());
 }


@@ -101,33 +101,16 @@ mul::mul(const exvector & v)
 {
     debugmsg("mul constructor from exvector",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
 {
     debugmsg("mul constructor from exvector",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;

-    overall_coeff=_ex1();
+    overall_coeff = _ex1();

     construct_from_exvector(v);
     GINAC_ASSERT(is_canonical());
 }
 
     construct_from_exvector(v);
     GINAC_ASSERT(is_canonical());
 }
 

-/*
-mul::mul(const epvector & v, bool do_not_canonicalize)
-{
-    debugmsg("mul constructor from epvector,bool",LOGLEVEL_CONSTRUCT);
-    tinfo_key = TINFO_mul;
-    if (do_not_canonicalize) {
-        seq=v;
-#ifdef EXPAIRSEQ_USE_HASHTAB
-        combine_same_terms(); // to build hashtab
-#endif // def EXPAIRSEQ_USE_HASHTAB
-    } else {
-        construct_from_epvector(v);
-    }
-    GINAC_ASSERT(is_canonical());
-}
-*/
-

 mul::mul(const epvector & v)
 {
     debugmsg("mul constructor from epvector",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
 mul::mul(const epvector & v)
 {
     debugmsg("mul constructor from epvector",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;

-    overall_coeff=_ex1();
+    overall_coeff = _ex1();

     construct_from_epvector(v);
     GINAC_ASSERT(is_canonical());
 }
     construct_from_epvector(v);
     GINAC_ASSERT(is_canonical());
 }


@@ -136,7 +119,7 @@ mul::mul(const epvector & v, const ex & oc)
 {
     debugmsg("mul constructor from epvector,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
 {
     debugmsg("mul constructor from epvector,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;

-    overall_coeff=oc;
+    overall_coeff = oc;

     construct_from_epvector(v);
     GINAC_ASSERT(is_canonical());
 }
     construct_from_epvector(v);
     GINAC_ASSERT(is_canonical());
 }


@@ -146,7 +129,7 @@ mul::mul(epvector * vp, const ex & oc)
     debugmsg("mul constructor from epvector *,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
     GINAC_ASSERT(vp!=0);
     debugmsg("mul constructor from epvector *,ex",LOGLEVEL_CONSTRUCT);
     tinfo_key = TINFO_mul;
     GINAC_ASSERT(vp!=0);

-    overall_coeff=oc;
+    overall_coeff = oc;

     construct_from_epvector(*vp);
     delete vp;
     GINAC_ASSERT(is_canonical());
     construct_from_epvector(*vp);
     delete vp;
     GINAC_ASSERT(is_canonical());


@@ -161,7 +144,7 @@ mul::mul(const ex & lh, const ex & mh, const ex & rh)
     factors.push_back(lh);
     factors.push_back(mh);
     factors.push_back(rh);
     factors.push_back(lh);
     factors.push_back(mh);
     factors.push_back(rh);

-    overall_coeff=_ex1();
+    overall_coeff = _ex1();

     construct_from_exvector(factors);
     GINAC_ASSERT(is_canonical());
 }
     construct_from_exvector(factors);
     GINAC_ASSERT(is_canonical());
 }

diff --git a/ginac/power.cpp b/ginac/power.cpp
index a85bc79f4fb2e366ec30773f83a593d2521cdec2..060d5a73c09109bf3d9a462df54ca0d56d724578 100644 (file)
--- a/ginac/power.cpp
+++ b/ginac/power.cpp
@@ -338,18 +338,18 @@ ex power::eval(int level) const
     const ex & ebasis    = level==1 ? basis    : basis.eval(level-1);
     const ex & eexponent = level==1 ? exponent : exponent.eval(level-1);
 
     const ex & ebasis    = level==1 ? basis    : basis.eval(level-1);
     const ex & eexponent = level==1 ? exponent : exponent.eval(level-1);
 

-    bool basis_is_numerical=0;
-    bool exponent_is_numerical=0;
+    bool basis_is_numerical = 0;
+    bool exponent_is_numerical = 0;

     numeric * num_basis;
     numeric * num_exponent;
 
     if (is_exactly_of_type(*ebasis.bp,numeric)) {
     numeric * num_basis;
     numeric * num_exponent;
 
     if (is_exactly_of_type(*ebasis.bp,numeric)) {

-        basis_is_numerical=1;
-        num_basis=static_cast<numeric *>(ebasis.bp);
+        basis_is_numerical = 1;
+        num_basis = static_cast<numeric *>(ebasis.bp);

     }
     if (is_exactly_of_type(*eexponent.bp,numeric)) {
     }
     if (is_exactly_of_type(*eexponent.bp,numeric)) {

-        exponent_is_numerical=1;
-        num_exponent=static_cast<numeric *>(eexponent.bp);
+        exponent_is_numerical = 1;
+        num_exponent = static_cast<numeric *>(eexponent.bp);

     }
 
     // ^(x,0) -> 1 (0^0 also handled here)
     }
 
     // ^(x,0) -> 1 (0^0 also handled here)


@@ -391,10 +391,10 @@ ex power::eval(int level) const
         if (basis_is_crational && exponent_is_crational
             && num_exponent->is_real()
             && !num_exponent->is_integer()) {
         if (basis_is_crational && exponent_is_crational
             && num_exponent->is_real()
             && !num_exponent->is_integer()) {

-            numeric r, q, n, m;
-            n = num_exponent->numer();
-            m = num_exponent->denom();
-            q = iquo(n, m, r);
+            numeric n = num_exponent->numer();
+            numeric m = num_exponent->denom();
+            numeric r;
+            numeric q = iquo(n, m, r);

             if (r.is_negative()) {
                 r = r.add(m);
                 q = q.sub(_num1());
             if (r.is_negative()) {
                 r = r.add(m);
                 q = q.sub(_num1());


@@ -402,29 +402,23 @@ ex power::eval(int level) const
             if (q.is_zero())  // the exponent was in the allowed range 0<(n/m)<1
                 return this->hold();
             else {
             if (q.is_zero())  // the exponent was in the allowed range 0<(n/m)<1
                 return this->hold();
             else {

-                epvector res(2);
+                epvector res;
+                res.reserve(1);

                 res.push_back(expair(ebasis,r.div(m)));
                 res.push_back(expair(ebasis,r.div(m)));

-                res.push_back(expair(ex(num_basis->power(q)),_ex1()));
-                return (new mul(res))->setflag(status_flags::dynallocated | status_flags::evaluated);
-                /*return mul(num_basis->power(q),
-                           power(ex(*num_basis),ex(r.div(m)))).hold();
-                */
-                /* return (new mul(num_basis->power(q),
-                   power(*num_basis,r.div(m)).hold()))->setflag(status_flags::dynallocated | status_flags::evaluated);
-                */
+                return (new mul(res,ex(num_basis->power(q))))->setflag(status_flags::dynallocated | status_flags::evaluated);

             }
         }
     }
 
     // ^(^(x,c1),c2) -> ^(x,c1*c2)
     // (c1, c2 numeric(), c2 integer or -1 < c1 <= 1,
             }
         }
     }
 
     // ^(^(x,c1),c2) -> ^(x,c1*c2)
     // (c1, c2 numeric(), c2 integer or -1 < c1 <= 1,

-    // case c1=1 should not happen, see below!)
+    // case c1==1 should not happen, see below!)

     if (exponent_is_numerical && is_ex_exactly_of_type(ebasis,power)) {
     if (exponent_is_numerical && is_ex_exactly_of_type(ebasis,power)) {

-        const power & sub_power=ex_to_power(ebasis);
-        const ex & sub_basis=sub_power.basis;
-        const ex & sub_exponent=sub_power.exponent;
+        const power & sub_power = ex_to_power(ebasis);
+        const ex & sub_basis = sub_power.basis;
+        const ex & sub_exponent = sub_power.exponent;

         if (is_ex_exactly_of_type(sub_exponent,numeric)) {
         if (is_ex_exactly_of_type(sub_exponent,numeric)) {

-            const numeric & num_sub_exponent=ex_to_numeric(sub_exponent);
+            const numeric & num_sub_exponent = ex_to_numeric(sub_exponent);

             GINAC_ASSERT(num_sub_exponent!=numeric(1));
             if (num_exponent->is_integer() || abs(num_sub_exponent)<1) {
                 return power(sub_basis,num_sub_exponent.mul(*num_exponent));
             GINAC_ASSERT(num_sub_exponent!=numeric(1));
             if (num_exponent->is_integer() || abs(num_sub_exponent)<1) {
                 return power(sub_basis,num_sub_exponent.mul(*num_exponent));