29 views (last 30 days)

Show older comments

I am using Matlab to calculate a double integral. The error is:

Warning: Explicit integral could not be found.

Error using mupadmex

Error in MuPAD command: DOUBLE cannot convert the input expression into a double

array. If the input expression contains a symbolic variable, use the VPA function instead.

Error in sym/double (line 827)

Xstr = mupadmex('symobj::double', S.s, 0);

Error in MeEqu7 (line 20)

out2 =out2 (int(int(yi,u1,0,u2),u2, 0 ,1000))

and my code is as follows:

clear all;

a=4;

la1=1/(pi*500^2); la2= la1*5;

p1=25; p2=p1/25;

sgma2=10^(-11);

index=1;

g=2./a;

syms r u1 u2 u3 aa

powe= 0;

seta= 10^powe;

xNor = (u2./u1).^(-1./g) + 1 ;

aa= -u2.^(-1./g) .*seta.*( pi.*(la1.*p1.^(g) + la2.*p2.^(g)) ).^(1./g) ;

yi = exp(aa./xNor - u2 );

out2 = double (int(int(yi,u1,0,u2),u2, 0 ,1000))

Thanks alot for your help!! I also tried vpa. but the out2 is not a number. Actually, I am expecting the out2 is a number, sothat I can use this number to plot a figure. I did some search here. but their case is not helpful for my case. Could you take a deeper look?

John D'Errico
on 25 Dec 2014

Edited: John D'Errico
on 25 Dec 2014

Look at what you have created.

(int(yi,u1,0,u2)

ans =

int(exp(- u2 - 1934281311383407/(1208925819614629174706176*u2^2*(u1^2/u2^2 + 1))), u1, 0, u2)

This expression has the symbolic variable u2 in it. The error says the integral could not be found, but even if it could have been found, what you are trying to stuff into out2 is a function of u2.

Perhaps you were planning on trying to integrate u2 out given your initial comment, but regardless, you cannot convert a symbolic expression with a variable in it, into a double precision number.

Edit: (in response to the change made to your question.)

You changed one line of your code. The last line of code is now

out2 = double (int(int(yi,u1,0,u2),u2, 0 ,1000))

So, where before, you had only a single integral, you now do have a double nested integral. But, really, nothing material has changed. The double integral has no analytical solution. The inner integral fails, so the outer one fails too.

You cannot convert that result into a double precision number. Wanting something to happen in mathematics is not sufficient. Demanding that the integral exist is not sufficient, by trying to convert it to a double.

int(int(yi,u1,0,u2),u2, 0 ,1000)

ans =

int(int(exp(- u2 - 1934281311383407/(1208925819614629174706176*u2^2*(u1^2/u2^2 + 1))), u1, 0, u2), u2, 0, 1000)

Maybe your goal is to do a numerical integration, and you were thinking that telling MATLAB to convert to double precision that it would magically know what you wanted. MATLAB does not read your mind. If the input to double is not a number or an array of numbers, it cannot somehow know to do a numerical integration instead of the symbolic integration you asked for.

(Note that IF you do decide to solve this using a numerical integration, that integrals of exponentials with extremely wide limits of integration often tend to fail, so you will need to do that integration carefully.) So now lets look at the integrand of this double integral.

exp(- u2 - 1934281311383407./(1208925819614629174706176*u2.^2.*(u1.^2./u2.^2 + 1)));

First of all, the ratio

vpa('1934281311383407/1208925819614629174706176')

ans =

0.0000000016000000000000002650866688430533

is essentially 1.6e-9 in double precision. Next, basic algebra lets us simplify that mess of a kernel a bit more.

exp(- u2 - 1.6e-9./(u1.^2 + u2.^2));

You want to convert this into a double precision result, with the integral on u1 having limits of 0 and u2, and the integral on u2 having limits of 0 and 1000.

So now compare the two parts of what we find inside the exponential. We have a minus u2 term, which varies from 0 to 1000, and a second term which is essentially zero in double precision over the ENTIRE DOMAIN of your integration, compared to a number that varies from 0 to 1000. The result of that double integral will clearly be just slightly less than 1, since we know that

int(exp(-u2),0,1000) == 1 - exp(-1000)

The double integration seems to be doable though.

u1max = @(u2) u2;

fun = @(u2,u1) exp(- u2 - 1.6e-9./(u1.^2 + u2.^2));

integral2(fun,0,1000,0,u1max,'reltol',1e-13,'abstol',1e-13)

Warning: Reached the maximum number of function evaluations (10000). The result fails the global error test.

> In funfun/private/integral2Calc>integral2t at 129

In funfun/private/integral2Calc at 9

In integral2 at 106

ans =

0.999999986764709

So it looks like I cranked down on the tolerances a bit too much, but the result seems reasonable. In the end, if you want to do a numerical integration, use a numerical integration tool, not a symbolic one.

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!