On 7/28/09, Rafa? Topolnicki <rtopolnicki at o2.pl> wrote:
>>> %t6 : a[n] = if n=1 then 'limit( subst(n=_n,rhs(%t6)),_n,n ) else
>>> %rhs(%t6);
That should be: if equal(n, 1) then ....
because "=" is a test for literal identity (i.e. same expression)
while equal is a test for equivalence (equal value).
n=1 always yields false when n is unbound, while equal(n, 1) yields
an unevaluated equal(n, 1) (i.e. postpone the evaluation to a Boolean
until n has a value).
> of course there should be
> (%i9) fourexpand(flist,x,%pi,1);
>
> Division by 0
> #0: fourexpand(l=[%t5,%t6,%t7],x=x,p=%pi,nn=1)
> -- an error. To debug this try debugmode(true);
Even with a[n] = if equal(n, 1) then 0 else rhs(...), I still get
"Divison by 0".
It turns out that is a bug in the depths of the sign-testing code,
which is triggered by the summation code which is called in fourexpand.
For the record, the bug looks like this:
if x > 1 then a else b;
=> if x > 1 then a else b; /* OK */
assume (x >= 1);
if x > 1 then a else b;
=> a /* OOPS! should be same as before */
and it originates from here:
?dcomp (1, 1.0);
=> pos
when that should be 'zero.
I have a patch to fix it which I'll commit for the next release (in August).
Agreed that the Fourier code should try to detect problematic
cases and apply a limit or whatever, so the user doesn't have to.
best
Robert Dodier