Which is the best method to find zero ? and choose the max value and plot a graph comproving the zeros.
2 次查看(过去 30 天)
显示 更早的评论
For a given parameters
ss=0.6;
T=0.1;
g=1;
b=0.5;
f1 = @(r) 1 - (b./r).^2 - (g^-2)*((2/15)*(ss)^9 *(1./(r - 1).^9 - 1./(r + 1).^9 - 9./(8*r).*(1./(r - 1).^8 - 1./(r + 1).^8)) -(ss)^3 *(1./(r-1).^3 - 1./(r+1).^3 - 3./(2*r).*(1./(r-1).^2 - 1./(r+1).^2)));
I want to find the r value when
f1 == 0
So, which is the best method ? findzero, fminbnd, newtonraphson ? fsolve ?
i try
a1= fzero(f1,1.00000001)
assume(r, 'real');
t=vpasolve(1 - (b/r)^2 - (g^-2)*(2/15*[(s/R)]^9 *(1/(r - 1)^9 - 1/(r + 1)^9 - 9/(8* r) *(1/(r - 1)^8 - 1/(r + 1)^8)) - [(s/R)]^3 *(1/(r -1)^3 - 1/(r + 1)^3 - 3/(2* r)* (1/(r - 1)^2 - 1/(r + 1)^2))) == 0,r,[-Inf,Inf]);
t
options=optimset('MaxIter',1e3,'TolFun',1e-10);
t1=fminbnd(func,0.1,15,options);
t1
but i dont know if it is correct. I want to plot a graph r x f1 to see the zeros.
Any help how i can find zeros exactly ?
0 个评论
回答(3 个)
dpb
2015-11-23
>> fsolve(f1,sqrt(2))
Equation solved.
fsolve completed because the vector of function values is near zero
as measured by the default value of the function tolerance, and
the problem appears regular as measured by the gradient.
<stopping criteria details>
ans =
1.4121
>> f1(ans)
ans =
-2.4100e-12
>> r=linspace(1.2, 1.6);
>> figure
>> plot(r,f1(r))
>> hold on
>> line([r(1) r(end)],[0 0])
>> ylim([-20 20])
>>
I started by plotting f1(r) for r=linspace(1+eps,100,20); and exploring from there to get a reasonable initial guess for fzero Looking at the algebra first would probably lead to a way to guess that as a closer starting point rather than just brute force, but it works, too, often... :)
Walter Roberson
2015-11-23
You cannot get exact values for r. The exact solutions involve roots of a 10th degree polynomial that does not factor. There is no closed form solutions for such a polynomial.
Star Strider
2015-11-24
It’s not easy, but you can find many of them using fsolve and a loop. I plotted these first to get the approximate range (-10,10):
ss=0.6;
T=0.1;
g=1;
b=0.5;
f1 = @(r) 1 - (b./r).^2 - (g^-2)*((2/15)*(ss)^9 *(1./(r - 1).^9 - 1./(r + 1).^9 - 9./(8*r).*(1./(r - 1).^8 - 1./(r + 1).^8)) -(ss)^3 *(1./(r-1).^3 - 1./(r+1).^3 - 3./(2*r).*(1./(r-1).^2 - 1./(r+1).^2)));
r = linspace(-10, 10, 1000);
r_e = -50:1.2:50;
for k1 = 1:length(r_e)
re_r_sol(k1) = fsolve(f1, r_e(k1));
end
re_r = uniquetol(re_r_sol, 0.05);
for k1 = 1:length(r_e)
cx_r_sol(k1) = fsolve(f1, 1i*r_e(k1));
end
[cx_ru,ia,ic] = uniquetol(abs(cx_r_sol), 0.05);
cx_r = cx_r_sol(ia);
fprintf(1,'\nReal roots:\n')
fprintf(1,'\t%+8.3f\n', re_r)
fprintf(1,'\nComplex roots:\n')
fprintf(1,'\t%+8.3f %+8.3fi\n', [real(cx_r); imag(cx_r)])
Real roots:
-2.528
-1.412
-0.576
+0.576
+1.412
+2.528
Complex roots:
+0.335 +0.377i
+0.000 +79.550i
+0.000 -84.050i
I am not entirely confident of these results, because I would expect the complex roots to be complex conjugates of each other. Using a symbolic solver might be more productive.
0 个评论
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!