Symbolic integration inside numerical integration

2018 3 17
1 个回答
回答已采纳
更新时间：2018 3 20
3 次查看（30 天）

0 个投票

I would like to solve the following integral numerically:
As far as I can tell, the built-in function integral2 is not applicable here due to the nature of the expression. Instead, from Evaluating Double Integrals, I arrive at the following solution which takes forever to compute due to the third line:
syms x r
firstint=int(1./(1+x.^2.5),x,r,Inf)
answer=int(r.*exp(-r.^2).*firstint,r,0,Inf),
double(answer)
The last line is just to get a numerical answer out of the symbolic integration. To speed up things, I am inclined to replace the third line with the built-in function integral, but it requires me to make firstint into a handle first, for example by using another built-in function called matlabFunction. Unfortunately, I get an error when doing the following:
syms x r
firstint=int(1./(1+x.^2.5),x,r,Inf),
firstintHandle = matlabFunction(firstint),
answer=integral(r.*exp(-r.^2).*firstintHandle,0,Inf),
Any ideas on what to do? Please be very specific. Your help is greatly appreciated.

请先登录，再回答此问题。

 采纳的回答

Star Strider
Star Strider 2018-3-17

0 个投票

See if this does what you want:
firstint = @(r) integral(@(x) 1./(1+x.^2.5), r, Inf, 'ArrayValued',1);
secondint = @(r) r.*exp(-r.^2).*exp(-firstint(r));
answer = integral(secondint, 0, Inf, 'ArrayValued',1)
answer =
0.274214226644371
The code appears to me to be correct. I defer to you to be certain it is.

6 个评论
显示 4更早的评论 隐藏 4更早的评论

Gunnar
Gunnar 2018-3-17
Wow, that integration took less than an instance. Many, many thanks!
Gunnar
Gunnar 2018-3-20
A quick follow-up if I may? It is a variant of the first integral:
integral(@(x)1./(1+x.^1.1),1,Inf,'ArrayValued',1);
but where I this time get the following error message:
Warning: Reached the limit on the maximum number of intervals in use. Approximate bound on error is 1.1e-01. The integral may not exist, or it may be difficult to approximate numerically to the requested accuracy.
Any ideas on what I'm doing wrong here?
Star Strider
Star Strider 2018-3-20
You are not doing anything wrong that I can see. If you set it without the 'ArrayValued' flag, it works, although it still complains with the Warning:
Q = integral(@(x)1./(1+x.^1.1),1,Inf)
Q =
9.424915202068878
It also has a closed-form solution:
syms x
f(x) = 1/(1 + x^1.1);
I = int(f, x, 1, Inf)
vpaI = vpa(I)
I =
(pi*100i)/121 + (10*log(20000000000))/11 + log((1771561*10^(7/11)*11^(4/11))/10000000) - (pi*exp((pi*1i)/11)*90i)/121 + (pi*exp((pi*2i)/11)*80i)/121 - (pi*exp((pi*3i)/11)*70i)/121 + (pi*exp((pi*4i)/11)*60i)/121 - (pi*exp((pi*5i)/11)*50i)/121 - (40*pi*exp((pi*1i)/22))/121 + (30*pi*exp((pi*3i)/22))/121 - (20*pi*exp((pi*5i)/22))/121 + (10*pi*exp((pi*7i)/22))/121 - (10*log(19487171000*exp((pi*1i)/11) - 19487171000)*exp((pi*1i)/11))/11 + (10*log(19487171000*exp((pi*2i)/11) + 19487171000)*exp((pi*2i)/11))/11 - (10*log(19487171000*exp((pi*3i)/11) - 19487171000)*exp((pi*3i)/11))/11 + (10*log(19487171000*exp((pi*4i)/11) + 19487171000)*exp((pi*4i)/11))/11 - (10*log(19487171000*exp((pi*5i)/11) - 19487171000)*exp((pi*5i)/11))/11 + (10*log(19487171000*exp((pi*6i)/11) + 19487171000)*exp((pi*6i)/11))/11 - (10*log(19487171000*exp((pi*7i)/11) - 19487171000)*exp((pi*7i)/11))/11 + (10*log(19487171000*exp((pi*8i)/11) + 19487171000)*exp((pi*8i)/11))/11 - (10*log(19487171000*exp((pi*9i)/11) - 19487171000)*exp((pi*9i)/11))/11 + (10*log(19487171000*exp((pi*10i)/11) + 19487171000)*exp((pi*10i)/11))/11
vpaI =
9.4316568296129897202042695044568 - 1.1754943508222875079687365372222e-38i
The result is complex, although the imaginary part is vanishingly small.
Gunnar
Gunnar 2018-3-20
编辑：Gunnar 2018-3-20
Thanks for your replies, Star Strider. I also managed to get a numeric answer based on the first approach, but the error message indicates that it might not be accurate. In comparison, your second solution takes much longer to compute due to the need to do symbolic integration, but on the other hand there is no error message. It is a tradeoff, and it is difficult to know what is the best course of action. Indeed, MATLAB seems to offer multiple ways to numerically evaluate integrals (see for example function 'trapz'), and some integrals may work with one method but not the other. With that said, many thanks for the help so far.

请先登录，再进行评论。

更多回答（0 个）

类别

在 帮助中心 和 File Exchange 中查找有关 Numerical Integration and Differentiation 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by