No plot in double integral
Mostrar comentarios más antiguos
Hi! I try to solve the double integral, but there is an empty plot. Could you please tell me, what is wrong? In MathCad this integral is solved well.
n=1;
t=1;
r=1;
s=0:0.1:5;
fun = @(x,z,k) (x.*exp(2*n*t.*x.^2)./sqrt(1-x.^2)).*(exp(-2*t.*z.^2).*besselj(0,z.*k.*x).*(1+(sqrt(-pi)./(z*r)).*(1+((z*r).^2)/2).*(erfi(z*r/2)).*(exp(-((z*r).^2)/4))));
f3 = arrayfun(@(k)integral2(@(x,z)fun(x,z,k),0,Inf,0,1),s)
Cor = ((-2*r*sqrt(2*n*t)*exp(-2*n*t))/(pi*erf(sqrt(2*n*t))*(atan(r/(2*sqrt(2*t)))-sqrt(2*t)/(2*r*(1/4+2*t/(r^2))))))*f3;
plot(s,Cor,'g-');
6 comentarios
See above. Integral2 returns that the integral does not exist because Inf or NaN values were returned from fun.
Is it correct that you want to integrate fun with respect to x from 0 to Inf and with respect to z from 0 to 1 ? Or vice versa ?
the cyclist
el 29 de En. de 2023
Editada: the cyclist
el 29 de En. de 2023
Can you upload the formula you are trying to integrate, in math notation? Maybe a mistake was made in transcribing it into MATLAB.
Are you expecting complex values?
Maybe you wanted x to be limited to 0 1 and z to be 0 to infinity? You coded with x being 0 to infinity and z being 0 to 1
n=1;
t=1;
r=1;
syms x z k s real
fun = (x.*exp(2*n*t.*x.^2)./sqrt(1-x.^2)).*(exp(-2*t.*z.^2).*besselj(0,z.*k.*x).*(1+(sqrt(-pi)./(z*r)).*(1+((z*r).^2)/2).*(erfi(z*r/2)).*(exp(-((z*r).^2)/4))));
f3 = subs(int(int(fun,x,0,inf), z, 0, 1), k, s)
fun_zhalf = limit(subs(fun, k, 1), z, 1/2)
xlimited = limit(fun_zhalf, x, 1)
case_to_plot = children(children(xlimited, 3),1)
fplot(case_to_plot, [0 2])
xlimited = limit(fun_zhalf, x, 2)
vpa(xlimited)
the cyclist
el 29 de En. de 2023
I was surprised to see the explicitly imaginary
sqrt(-pi)
in the expression, which is what prompted me to ask to see the math notation version.
I did not notice that. I wonder what it would look like if you used sqrt(pi) instead?
n=1;
t=1;
r=1;
syms x z k s real
Pi = sym(pi);
fun = (x.*exp(2*n*t.*x.^2)./sqrt(1-x.^2)).*(exp(-2*t.*z.^2).*besselj(0,z.*k.*x).*(1+(sqrt(Pi)./(z*r)).*(1+((z*r).^2)/2).*(erfi(z*r/2)).*(exp(-((z*r).^2)/4))));
f3 = subs(int(int(fun,x,0,inf), z, 0, 1), k, s)
fun_zhalf = limit(subs(fun, k, 1), z, 1/2)
xlimited = limit(fun_zhalf, x, 1)
case_to_plot = children(children(xlimited, 3),1)
fplot(case_to_plot, [0 2])
xlimited = limit(fun_zhalf, x, 2)
vpa(xlimited)
You get rid of the real-valued portion of the function (at least at some values), but the complex portion remains.
Hexe
el 30 de En. de 2023
Respuestas (0)
Categorías
Más información sobre Programming en Centro de ayuda y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!













