[scilab-Users] Strange problem with "integrate" : "Convergence problem..." error 24 (Scilab 4.1.2)

[Michaël Baudin]

Hi Remi,

I inquired your problem and it seems to me that all is fine,
that is, the problem you face is a limitation of the computationnal
method and not a bug.

The "integrate" command is a Scilab macro which is
based on the "intg" Scilab primitive.

Now try the following Scilab script :

function y=f(x)
  y=sin(x)
endfunction
a = 0;
b=2*%pi;
computed = intg(a,b,f);

The expected result in that case is 0., since -cos(x) is periodic
on the [0,2 PI] range.
The same error message appears.
The following is the call stack, which shows how the quapro library
is used to compute the integral based on a "21-point gauss-kronrod 
approximation".

 differential_equations/sci_gateway/c/sci_intg
 > differential_equations/sci_gateway/fortran/intg.f
 > differential_equations/src/fortran/dqag0.f
 > dqags > quarul

By doing a step-by-step check of the process, it
appears that the numerical integration is correctly computed (near the
machine's epsilon), but a secondary check sets a boolean flag to a
value indicating that the error for that computation is too large.

A simple workaround for your problem may be to decompose the
range in two parts :

i1 = intg(0,%pi,f);
i2 = intg(%pi,2*%pi,f);

then the expected result is the sum of i1 and i2.
Other solutions include traditionnal ode solver (like dassl ) or 
Monte-Carlo methods.

Best regards,

Michaël Baudin

Remi G. a écrit :
> Hello!
>
> I'm facing a strange problem with the function "integrate" (and also 
> with "intg") with Scilab 4.1.2 :
>
> If I want to integrate x^5 between -%pi and %pi, it works (hopefully) :
>  integrate('t^5','t',-%pi,%pi)
>
> But if replace 5 by a variable k (that I've defined before of course) 
> it will not work for k>=5!
>  i=5
>  integrate('t^i','t',-%pi,%pi)
>
> Here is the error message I get :
>    !--error 24
>  convergence problem...
>  At line      26 of function integrate called by :
>  integrate('t^i','t',-%pi,%pi)
>
> That's really annoying because I can't code my least square 
> approximation procedure which is... for tomorrow!
> If someone could give me a way to go through this issue!
> Thank you very much!
>
> Rémi.