[Maxima-bugs] [ maxima-Bugs-1965640 ] Problems with $specint

[SourceForge.net <no...@so...>]

Bugs item #1965640, was opened at 2008-05-16 22:04
Message generated for change (Comment added) made by crategus
You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1965640&group_id=4933

Please note that this message will contain a full copy of the comment thread,
including the initial issue submission, for this request,
not just the latest update.
Category: Lisp Core - Integration
Group: Includes proposed fix
Status: Open
Resolution: None
Priority: 5
Private: No
Submitted By: Crategus (crategus)
Assigned to: Nobody/Anonymous (nobody)
Summary: Problems with $specint

Initial Comment:
I would like to open a bug report to collect the results of investigations of the function $specint. To get an overview of the problems with $specint I have taken the 99 tabulated Laplace transforms from the website of EqWorld.

55 of the 99 examples fail with the original code. I have divided the problems in the following cases:

1. Maxima has no algorithm:

In most cases Maxima gives an internal symbol like other-defint-to-follow-negtest or arbpow-failed. This is a known problem. In some cases we get a correct noun form of the unevaluated integral. At last there are many problems which gives a wrong answer. 

A lot of examples include terms like t^(-1) ... or t^(-1/2) ... Maxima can't calculat these integrals but we know solutions.  

A simple type of integral Maxima can't evaluate is the division by the sum of constants:

(%i7) specint(%e^(-s*t)/(x+y),t);
(%o7) other-defint-to-follow-negtest

In this cases the exponential function is hidden in a summation. I have found a correction which works generally and gives the correct result:

(%i7) specint(%e^(-s*t)/(x+y),t);
(%o7) (1/(x+y)*s)

2. Maxima has an algorithm for a special function but don't give the correct result:

We get no results for functions like bessel_k, bessel_y, log, erf, erfc etc. For all these functions Laplace transforms are tabulated. Beside the test of EqWorld I tried to get results for the internal functions %l[n,a](x) - the Laguerre function - or %he[n](x) - the Hermite function. But I dont' get the expected result.

Here the example for the Laguerre function:

(%i6) kill(all);
(%o0) done

(%i1) assume(s>0,n>0),declare(n,integer);
(%o1) [s > 0,n > 0]

(%i2) specint(%e^(-2*t)*%l[n,0](t),t);
Maxima encountered a Lisp error:
Error in MACSYMA-TOP-LEVEL [or a callee]: $N is not of type NUMBER.
Automatically continuing.
To reenable the Lisp debugger set *debugger-hook* to nil.

After correction of the code (the correction is not included in the diff appended to the next post):

(%i6) kill(all);
(%o0) done
(%i1) assume(s>0,n>0),declare(n,integer);
(%o1) [s > 0,n > 0]
(%i2) specint(%e^(-s*t)*%l[n,0](t),t);
(%o2) (1-1/s)^n/s

In the case of the Laguerre function I have found a bug in the transformation. After correction, Maxima gives the expected result shown above for the Laguerre function. There may be further bugs. Or limitations of the algorithmen prevend the calculation of results. At least, these limitations should be documented for the user.

3. Maxima gets extra factors or terms in the result

A simple example is the bessel_i function. Here we get an additional phase factor: 

  %e^(%i*%pi*v/(-1)^(v2/2) 

in all calculations. This factor vanish when we introduce a small correction to the code. In other cases the problem seems to be more difficult.

Dieter Kaiser



----------------------------------------------------------------------

>Comment By: Crategus (crategus)
Date: 2008-05-24 23:16

Message:
Logged In: YES 
user_id=2039760
Originator: YES

I have attached the test file for the routine LT-EXP.
File Added: test_hypgeo_LT-EXP.mac

----------------------------------------------------------------------

Comment By: Crategus (crategus)
Date: 2008-05-24 23:10

Message:
Logged In: YES 
user_id=2039760
Originator: YES

I have collected examples to test only the routine LT-EXP and the code for
doing the Laplace transforms. I have found 30 examples tabulated by EqWorld
and A&S testing this routine.

With the original code 17 examples fail. The first changes showed in this
thread reduce this to 6 examples. Correcting a bug in LT-EXP all examples
work.

8. Bug in LT-EXP and in the following routines F29P149, F36P147 and
F37P147

The reason for the 6 remaining errors is a bug in the calculation of the
Laplace transform in the routines LT-EXP, F29P146, F36P147 and F37P147. The
constant term is missing. Additionally it is useful to check v=0 before
calling F36P147 and F37P147. For v <>0 we have no Laplace transform.
Without this check we can get wrong results.


BUT: 

For 3 examples I have an additional factor in the result of Maxima which
is not present in tabulated result of EqWorld. In one case the factor is
%e^(s^2/(s*a)) and in two cases the factor is 4. Furthermore, there is one
case Maxima get a different sign of a term. Because the algorithm works
very well in all other cases, I think there are errors in the table of
EqWorld. But, it is better to verify thesw cases independly.

Second:

I have changed code in the routine LT-EXP and F35P147 to get the results
for expression with sqrt(t) in the exponent (Point 6). The changed code
works very well and gives a lot of additional results involving
trigonometric and hyperbolic functions with sqrt(t) as argument. But, I
haven't verified why this does work mathematically.

I have attached a file with an update of the changes.

Dieter Kaiser

File Added: diff_hypgeo.txt

----------------------------------------------------------------------

Comment By: Crategus (crategus)
Date: 2008-05-16 22:11

Message:
Logged In: YES 
user_id=2039760
Originator: YES

As a first step to improve the code of $specint I would like to present 7
changes:

1. Function DEFEXEC

If we cant't find a parameter, we apply $factor to the expression. Now
Maxima finds the result for expressions like  t/(x+y) --> 1/(s^2(x+y))
where x and y are free of the integration variable.

2. Function ARBPOW1

I have specialized the pattern match to be sure that in the expression
c*t^v the parameter c is free of the integration variable. This condition
now will fail if we enter $specint with expression like u(t) or
t^(1/2)*(a+t)^(-1) and with the changes below we get nice and correct noun
forms.

3. Function LT-SF-LOG

Because we have specialized the pattern match we add at the end of the
function as return value a noun form.

4. Function LT-ARBPOW

A lot of integrals fail at this point. We add as the return value a noun
form. 

5. Function LT-SF-LOG, Condition ONEI

This is an example how we can avoid additional phase factors. If we use %i
directly in the calculation all additional phase factors in the
calculations vanish and the results are correct. There are more places we
can apply this change to obtain easier results.

6. Bug in LT-EXP and F35P147

I have found a bug in the routine lt-exp and f35p147. This bug prevents
the calculation of integrals with e.g. sin(2*sqrt(a*t)). $SPECINT gives the
result 0.
Here the output of Maxima after correction:

(%i6) radcan(specint(%e^(-s*t)*sin(2*sqrt(a*t)),t));
(%o6) sqrt(%pi)*sqrt(a)*%e^-(a/s)/s^(3/2)

That is perfectly the tabulated expression and $SPECINT now works for a
lot of other integrals too.

7. Extension of the algorithm

To show how we can extend the algorithm of $SPECINT to calculate further
integrals, I have added code to calculate integrals of the form
t^-1*(%e^(-a*t)-%e^(-b*t)). The code works also for integrals like
t^-1*sin(a*t). Here an example

(%i4) specint(%e^(-s*t)*t^-1*sin(a*t),t);
(%o4) %i*(log(s-%i*a)-log(s+%i*a))/2

That's equivalent to the tabulated answer atan(a/s).


With this changes problems 55, 150 and 157 of rtest14.mac will produce
different results. In all cases the noun form is improved and now more
correct.

The numbers of correct results of the test file test_eqworld.mac is
increased to 87. When Maxima can't evaluate the integral but returns a
correct noun form I declared the test as "(OK noun form)". There are 12
remaining examples which fails. This examples mostly include the log or erf
function. I think there is something wrong with the mathematic.

I have added a diff to show the above described changes to the code. 

Hint: 

The test file will stop 5 times and ask for the sign of the internal
variable psey. That's a known bug. I have this bug not remarked as an
error, because the results are correct. I try to find the reason of the
bug.

Dieter Kaiser

File Added: diff_hypgeo.txt

----------------------------------------------------------------------

You can respond by visiting: 
https://sourceforge.net/tracker/?func=detail&atid=104933&aid=1965640&group_id=4933