Negative values while solving equation with gamma functions

조회 수: 8 (최근 30일)
Sreeja
Sreeja 2014년 7월 3일
편집: John D'Errico 2014년 7월 3일
I am trying to solve the following equation which involves a gamma function. vpasolve(gamma(x+1.77)^2==gamma(x)*gamma(x+3.54)-gamma(x+1.77)^2,x), but the answer that I get is negative. Is there a way to solve this?

답변 (1개)

John D'Errico
John D'Errico 2014년 7월 3일
편집: John D'Errico 2014년 7월 3일
Why do people feel the need to use vpasolve?
fzero will be faster, just not quite as accurate. But if your final goal is a double anyway, then why use a mack truck to take a pea to Boston? You can't store 40 digits of precision in a double.
Regardless, PLOT YOUR FUNCTION! Don't just throw up your hands and stop thinking!
ezplot(fun,[3,4])
Had you used a plot to look for a root, you would have found one. Then just give vpasolve a starting value.
vpasolve(gamma(x+1.77)^2==gamma(x)*gamma(x+3.54)-gamma(x+1.77)^2,x,4)
ans =
3.3472921715874174466992415057515
You could also have used a simple scheme to look for zero crossings programmatically, to then give you a bracket for the desired positive root.
  댓글 수: 2
Sreeja
Sreeja 2014년 7월 3일
편집: Sreeja 2014년 7월 3일
Thanks John.
Suppose negative values are obtained while running a for loop, then how can we solve it? For example,
y=[0.4:0.1:1]'; n=numel(y); syms x sym('result',size(y));
for j=1:n, result(j,1)= vpasolve(y(j)*y(j)*gamma(x+(2/1.0720))*gamma(x+(2/1.0720))==(gamma(x)*gamma(x+(4/1.0720))-(gamma(x+(2/1.0720))*gamma(x+(2/1.0720)))),x); end
Also which is the best effective way to solve this equation?
John D'Errico
John D'Errico 2014년 7월 3일
편집: John D'Errico 2014년 7월 3일
As I suggested, generate a list of points of the function, then test for a sign change.
Yes, you DO have a function. Inside your loop, define a function handle as:
fun = @(x) -y(j)*y(j)*gamma(x+(2/1.0720)).*gamma(x+(2/1.0720)) + (gamma(x).*gamma(x+(4/1.0720))-(gamma(x+(2/1.0720)).*gamma(x+(2/1.0720))));
For a given value of y(i), this is a function of x. Evaluate at a list of points. Check to find the first consecutive pair of points such that
fun(x)*fun(x+dx) <= 0
Pass that as a bracket into fzero. Note that the test above can be computed in one line, with no need for a loop. Find will do most of the work.
Whats the problem?

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Loops and Conditional Statements에 대해 자세히 알아보기

Community Treasure Hunt

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

Start Hunting!

Translated by