Symbolic integration in matlab

5 次查看(过去 30 天)
nathan khosla
nathan khosla 2011-10-2
I have: function problemOne
%declare symbolic Variables
N=6.023e23;
k=1.3806503e-23; %units (m^2)kg/(s^2)K
e=k*213.96; %units of k * 1/K
r0 = 426; %units of pm
syms r T;
%set y to be function inside integral in B
y=exp(((-4*e)*(((r0/r)^12)-((r0/r)^6))/k*T)-1)*(4*pi*(r^2));
%define B in terms of y
B= (-1/2)*N*int(y,r,0,Inf);
%solve for B=0;
solve(B,T)
end
This outputs:
Warning: Explicit solution could not be found.
> In solve at 81
In problemOne at 21
ans =
[ empty sym ]
I think the problem is in the integral, since
int(y,r,0,Inf)
yields
Inf - limit(4*pi*r^3*exp(-1)*exp(4005180207889834777988627621265408/783018253290173*T/r^6)*exp(-23937547659622311621519348000222281120506532855808/783018253290173*T/r^12), r = 0, Right)
What do I do? I thought I had a basic knowledge of Matlab (not an expert by any means), but this is confusing me to no end.
Thanks, Nkk

回答(3 个)

William
William 2011-10-2
When you're trying to run the whole thing at once debugging can be next to impossible. My guess is you're getting infinity somewhere in the problem or you have a syntax error within the code.
Try this.
Break it into pieces. have the program print out all the pieces and check the answers to see if they mathematically make sense.
I cannot see what you have specifically written on line 21 which makes it difficult to see what is going on.
Also see if you can plot the integral before you attempt to integrate it. If it is going to infinity than integrating it will lead to infinity.
try the
ezplot() function.

Walter Roberson
Walter Roberson 2011-10-2
y goes to infinity as r goes to infinity. The exp() part has constants divided by r^6 and r^12, so those divisions are going to go to 0 quickly (regardless of T), leading to exp(0). The r^2 term at the end of y therefor dominates, and r^2 strictly increases as r increases so the behavior of y going to infinity is not "removable".
The value at the other end, near r = 0, depends a great deal upon the value of T, especially when r is near 426 and T is near 0.

nathan khosla
nathan khosla 2011-10-2
I know...I think that is my problem. However, I have no idea how to solve it. In the end, I really only want to know when B=0. Is there any easier way to do this? Do I need symbolics? I have a friend who is a bit better at Matlab than I am (read:he is very good) and he said that I may not need them (it was over the phone, so he never saw the code, thus the uncertainty). So...do I need them?
Thanks, Nkk
  1 个评论
Walter Roberson
Walter Roberson 2011-10-2
As B is a constant multiplied by the integral, B will be 0 if and only if the integral is 0. As indicated above, though, the integral goes to infinity as r increases, and as that is not over an infinitely small region, the integral will always be infinite for every value of T, even for T = 0 or T negative.
As at least some points are going to come out positive, in order to get 0 as the integral you would have to have some points come out negative. exp() is, however, never negative unless T is complex. (And even if T is complex, the division by r^6 or r^12 sends the exponential quickly.)
The only remaining possibility is if T was complex and the integral in the region near 0 was so large as to cancel out the remaining infinite tail. Is a complex T plausible in the situation?

请先登录,再进行评论。

Community Treasure Hunt

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

Start Hunting!

Translated by