Negative values while solving equation with gamma functions

8 vues (au cours des 30 derniers jours)
Sreeja
Sreeja le 3 Juil 2014
Modifié(e) : John D'Errico le 3 Juil 2014
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?

Réponses (1)

John D'Errico
John D'Errico le 3 Juil 2014
Modifié(e) : John D'Errico le 3 Juil 2014
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 commentaires
Sreeja
Sreeja le 3 Juil 2014
Modifié(e) : Sreeja le 3 Juil 2014
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 le 3 Juil 2014
Modifié(e) : John D'Errico le 3 Juil 2014
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?

Connectez-vous pour commenter.

Catégories

En savoir plus sur Loops and Conditional Statements dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by