Page 9 - 6736
P. 9
% N - максимальне число ітерацій
% eps - граничне відхилення f(x) у нулі
%------------------------------------------------------------------------
N=50;
epselon=1e-7;
eps=1e-7;
%------------------------------------------------------------------------
%Вихід:с - корінь рівяння f(x)=0
% Yc - значення функції f(x) при x=c
% err - похибка у знаходженні кореня
%====================================
%Функцію f(x) задати у підпрограмі fun_Solve
%Вхід: х0 - початкова точка інтервалу
% xk - кінцева точка інтервалу
% delta -приріст аргументу функції f(x)
%Вихід: Графік функції f(x)
%====================================
N=100;
x0=0.51;
xk=3;
delta=(x0-xk)/N;
%------------------------------------------------------------------------
X=x0:delta:xk;
q=2*c*(sqrt(1-(b./X).^2)+sqrt(1-m^2*(b./X).^2))./(pi*X);
Y=X.*sqrt(q/(k*a1).*(a0*sqrt(k*q/a1-a2));
plot(X,Y,'k-')
%====================================
%Із графіка функц. знаходимо інтервал
%локального нул. Їх знач. вводимо
%у командному вікні
%====================================
disp('Введіть а');
a=input('a=');
disp('Введіть b');
b=input('b=');
%------------------------------------------------
Ya=fun_Solve(a,B,S);
Yb=fun_Solve(b,B,S);
if Ya*Yb>0
error('Неправильний інтервал локального нул.');
break
end
%====================================
%Реаліз. процедури методу хорд
%====================================
for k=1:N
Dx=Yb*(b-a)/(Yb-Ya);
c=b-Dx;
ac=c-a;
Yc=fun_Solve(c,B,S);
if Yc==0,break
elseif Yb*Yc>0
b=c;
Yb=Yc;
else
a=c;
Ya=Yc;
end
Dx=min(abs(Dx),ac);
if abs(Dx)<epselon,break,end
if abs(Yc)<eps,break,end
end
omega=sqrt(a2/a0);
%xsol=(b+a)/2;
Yc=fun_Solve(c,B,S);
9