Bug in diff

[Richard B. Kreckel]

Hi,

On Mon, 10 Mar 2003, Ola Skavhaug wrote:
> I think I've encountered a bug in GiNaC-1.0.13 (on a Debian Linux machine)
> In ginsh:
>
> > diff(sqrt((-sin((x)))^(2.0)),x);
> (1.0+5.0165576136843360246E-20*I)*sin(x)*cos(x)*((-sin(x))^(2.0))^(-1/2)
>
> The answer should be:
> sin(x)*cos(x)*((-sin(x))^(2.0))^(-1/2)
>
> I have not tried to figure out why this happens (perhaps an uninitialized
> variable?)
> I've compiled ginac by myself, using cln-1.1.5 (gcc 3.2 compiler)

I haven't any idea why you think this is due to an uninitialized variable,
but, well, it isn't a bug, strictly speaking.  Such behavior is to be
expected when dealing with floating points.  You should use an integer
exponent.

In your example (or the slightly simpler diff(sqrt((-x)^(2.)),x)), some
stages after power::diff() a mul::combine_overall_coeff(ex, ex) happens,
and therin numeric::power() is being called, which in turn calls
cln::expt(-1,1.0).  That is not guaranteed to be equal -1 by CLTL.

However, I do think that we should be somewhat more careful trying to
check for unit exponents, but for another reason: When
numeric::power_dyn() returns a reference, it can return a pointer to an
*existing* object instead of creating a new one on the heap.  This will
account for some more comparisons by pointer later on and hence be a
performance win.  For reasons of symmetry, the same check should be done
in numeric::power().

I have committed such a patch to the 1.1 and the 1.2 branches.  But
surprises are still gonna happen with floating point examples unless we
start using BCD which I'm sure you don't want.  ;-)

Cheers
   -richy.
-- 
Richard B. Kreckel
<Richard.Kreckel at GiNaC.DE>
<http://www.ginac.de/~kreckel/>