Thank you for reminding, I almost forget to add this to toolbox - will add right now.

Yes, absolutely! Just forgot to add it to the list.

Actually I saw the syntax in GCARE too :).

Not sure if I get the meaning of the output - how re-ordering could fail? Probably will just return 1.

The issue was probably related to the fact that we didn't digitally sign all components of toolbox.

Now we do sign all the components and probably this is resolved.

Let me know if it still occurs.

This bug has been resolved in Windows version. You can download it from our website.

GNU Linux & Mac OSX will follow today.

In power function we handle inputs like complex numbers by default. Semantic for various singular cases in complex arithmetic is not defined uniformly and implementation dependent (e.g. see here: http://www.advanpix.com/2015/10/19/devnotes-1-multiplicative-operations-with-complex-infinity/)

That is the reason for such misbehavior. We will fix this today.

Thank you for reporting the issue

I am not sure, probably you need to use better precision (from numerical point of view).

But why don't you compute the integral analytically?

Try using Maple or Mathematica. I have just tried Maple - it gave long but sensible expression for the any r,p,q,s.

Hi again,

I couldn't find your last comment (did you delete it?).

Just few comments on the code.

The function should be able to accept the vector-parameters (quadqk is a vectorized method). So instead of

f = @(x)(cos(p*x)+cosh(p*x)+q*sin(p*x)+q*sinh(p*x))*(cos(r*x)+cosh(r*x)+s*sin(r*x)+s*sinh(r*x));

use element-wise multiplication:

f = @(x)(cos(p*x)+cosh(p*x)+q*sin(p*x)+q*sinh(p*x)).*(cos(r*x)+cosh(r*x)+s*sin(r*x)+s*sinh(r*x));

Also be careful with precision and interval counts. My example shows extreme values, please use more reasonable at first, e.g.:

mp.Digits(34);
...
[q,errbnd] = quadgk(f,mp(0),L,'RelTol',100*mp.eps(),'AbsTol',100*mp.eps(), 'MaxIntervalCount', 650)

Quadgk will warn you if it needs more levels of recursion (MaxIntervalCount)

Yes, of course, toolbox follows standard Matlab syntax and rules.

If you can do something with standard Matlab - most likely it can be done with toolbox too.

Here is simple example for demonstration:

mp.Digits(50);
f = @(x)sin(x);
[q,errbnd] = quadgk(f,mp(0),mp(1),'RelTol',100*mp.eps,'AbsTol',100*mp.eps, 'MaxIntervalCount', 1000000)
q = 
    0.45969769413186028259906339255702339626768957938206
errbnd = 
    7.872028207137840807482477381844110449433981844857e-51

Define your integrand as a function, call quadgk with corresponding parameters.

Make sure the values 0, L and tolerance are also of mp-type. Toolbox will detect the mp-objects and will call its own multiprecision version of quadgk.

Let me know how it works.