Double symbolic integral - "Explicit integral could not be found"
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
At i=1 and j=3 the following is displayed "Warning: Explicit integral could not be found." A preliminary evaluation revealed that F(i=1,j=4)is failing to integrate x over the specified range with the result still containing the variable x.
Results:
F(i=1,j=2) "No Problem"
F = -(154946195819313285227594878996707*y*exp(-y^2)*((7186705221432913*2^(1/2)*pi^(1/2)*(erf((2^(1/2)*y)/2) + 1))/36028797018963968 - 1)^2)/81129638414606681695789005144064
F(i=1,j=3) "No Problem"
F = (1113552634535825382050544662406231220558348417491*pi^(1/2)*y*exp(-y^2/2)*(erf(y) + 1)*(7186705221432913*2^(1/2)*pi^(1/2) + 7186705221432913*2^(1/2)*pi^(1/2)*erf((2^(1/2)*y)/2) - 36028797018963968))/52656145834278593348959013841835216159447547700274555627155488768
F(i=1,j=4) "Warning: Explicit integral could not be found."
F = int((154946195819313285227594878996707*x*y*exp(-x^2/2)*exp(-y^2/2)*((7186705221432913*2^(1/2)*pi^(1/2)*(erf((2^(1/2)*x)/2) + 1))/36028797018963968 - (7186705221432913*2^(1/2)*pi^(1/2)*(erf((2^(1/2)*y)/2) + 1))/36028797018963968)^2)/81129638414606681695789005144064, x == -Inf..y)
Code:
clear all
clc
n = 4;
syms i j x y t w;
C2 = (2*pi())^(-0.5);
C3 = factorial(n)/(factorial(i-1)*factorial(j-i-1)*factorial(n-j));
fx = C2*exp(-0.5*x^2);
fy = C2*exp(-0.5*y^2);
ft = C2*exp(-0.5*t^2);
fw = C2*exp(-0.5*w^2);
Fx = int(ft,-inf,x);
Fy = int(fw,-inf,y);
for i=1:n/2
for j=(i+1):3
EC3 = eval(C3);
Fxy = EC3*x*y*fx*fy*Fx^(i-1)*(Fy-Fx)^(j-i-1)*(1-Fy)^(n-j);
F = int(Fxy,x,-inf,y);
Exy(i,j) = double(int(F,y,-inf,inf))
end
end
4 comentarios
Walter Roberson
el 12 de Jun. de 2015
Continuity is a really insufficient criteria for closed form integrals to exist. For example an ellipse is continuous but extensive study over centuries has found no closed form integral for its arc length. http://en.m.wikipedia.org/wiki/Arc_length
Maple for one finds no closed form integral for the case you raise here.
Respuestas (1)
Walter Roberson
el 12 de Jun. de 2015
If that integral has a closed form at all, then it is difficult to find. MATLAB warns you that it has left an integral unevaluated. But it is only a warning. When it comes time to do the double() on the nested unevaluated integrals, a numeric integration will be done to approximate the value. The warning helps you to be aware that the numeric result you get out will be an approximation whereas the other results are exact to within floating point round-off.
If you really do not want to see the warning, then you can use woff(). See lastwarning() to find the identifier of the warning to be able to turn it off.
5 comentarios
Walter Roberson
el 14 de Jun. de 2015
This is potentially a bug in how vpa() does its work. You could try,
G = evalin(symengine, 'numeric::int', F, sym('y=-infinity..infinity'));
Exy = double(G);
but it may use up all of your memory.
Walter Roberson
el 14 de Jun. de 2015
When I ask Maple to do a numeric integration using 13 digits or fewer, it finds a solution in not too long, of about -0.9549296585514. But when I ask it to do 14 or 15 digits, it returns the integral unevaluated. When I ask for 16 digits it runs through all of my memory and demands more.
When I plot I do not see anything about the values that would lead to that behaviour. The peak is at roughly y = -1.2, x = -1.9 and is positive but the values slip negative fairly close-by so it is plausible that the overall integral is negative.
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!