2 семестр / ПОСОБИЕ_ВычМат

21

Readln(a,b); sgna:=sgn(f(a)); if sgna*sgn(f(b))<0 then begin

repeat

x:=(a+b)/2; sgnx:=sgn(f(x));

if sgnx=0 then break else if sgna=sgnx then a:=x else b:=x; until (b-a)<e;

writeln('x=',x:12:8); end

else writeln('Метод не применим'); readkey;

end.

1.2.8.2 Метод итераций

program Pitr; uses crt;

const e=0.0000001; {точность приближения} function f(x:real):real; {функция f(x)}

begin f:=exp(x*ln(5))-6*x-3; end;

function f1(x:real):real; {производная функции f(x)} begin f1:=ln(5)*exp(x*ln(5))-6; end;

function sgn(z:real):integer; {функция вычисления знака}

begin if z=0 then sgn:=0 else if z>0 then sgn:=1 else sgn:=-1; end; {основная программа}

var a,b,x0,x1,x,L,dx,M,d,p:real; begin

clrscr;

writeln('Корни отделены на отрезках [-1;0], [1;2]. Введите значения a,b'); Readln(a,b);

{Найти наибольшее значение модуля производной } M:=(f1(a)); x:=a; dx:=0.01;

repeat x:=x+dx; p:=f1(x); if abs (p) > abs (M) then M:=p; until x>b; L := sgn (M)*2/abs(M); {Определить параметр}

{Начать процесс итераций} x0:= a;

repeat x1:=x0 - L*f(x0); d:=abs(x1-x0); x0:=x1; until d<e; writeln('Корень x = ', x1:12:7); readkey;

end.

1.2.8.3 Метод Ньютона

program PNtn; uses crt;

const e=0.0000001; {точность приближения}

22

function f(x:real):real; {функция f(x)} begin f:=exp(x*ln(5))-6*x-3; end;

function f1(x:real):real; {производная функции f(x)} begin f1:=ln(5)*exp(x*ln(5))-6; end;

function f2(x:real):real; {вторая производная функции f(x)} begin f2:=sqr(ln(5))*exp(x*ln(5)); end;

function sgn(z:real):integer; {функция вычисления знака}

begin if z=0 then sgn:=0 else if z>0 then sgn:=1 else sgn:=-1; end; {основная программа}

var a,b,x0,x1,x,L,dx,M2,m1,d,p:real; begin

clrscr;

writeln('Корни отделены на отрезках [-1;0], [1;2]. Введите значения a,b'); Readln(a,b);

{Найти наибольшее значение модуля второй производной} M2:=f2(a); x:=a; dx:=0.01;

repeat x:=x+dx; p:=f2(x); if abs (p) > abs (M2) then M2:=p; until x>b; {Найти наименьшее значение модуля производной} m1:=abs(f1(a)); x:=a; dx:=0.01;

repeat x:=x+dx; p:=abs(f1(x)); if p < m1 then m1:=p; until x>b; {Определить начальное приближение корня}

if sgn (M2) = sgn (f(a)) then x0 := a else x0 := b; {Начать процесс итераций}

repeat x1:=x0 - f(x0)/f1(x0); d:=sqr(x1-x0); x0:=x1; until d<2*m1*e/abs(M2); writeln('Корень x = ', x1:12:7); readkey;

end.

1.2.8.4 Метод хорд

program Phorda; uses crt;

const e = 0.0000001;

function f(x:real):real; {функция, для которой вычисляется нуль}

begin	f:=exp(x*ln(5))-6*x-3; end;
function f1(x:real):real; {производная функции f(x)}
begin	f1:=ln(5)*exp(x*ln(5))-6; end;
function f2(x:real):real; {вторая производная функции f(x)}
begin	f2:=sqr(ln(5))*exp(x*ln(5)); end;
function sgn(z:real):integer; {функция вычисления знака}
begin	if z=0 then sgn:=0 else if z>0 then sgn:=1 else sgn:=-1; end;

var {основная программа}

a, b, dx, fx, x, M, m1, p, v, d : real; begin

clrscr;

writeln('Корни отделены на отрезках [-1;0], [1;2]. Введите значения a,b');

23

Readln(a,b);

if sgn(f(a)) <> sgn(f(b)) then begin

{Установить вариант сужения интервала изоляции корня} if sgn(f1(a)) = sgn (f2(a)) then v := 1 else v := 2;

{Найти наибольшее значение модуля производной} M:=abs(f1(a)); x:=a; dx:=0.01;

repeat x:=x+dx; p:=abs(f1(x)); if p > M then M:=p; until x>b; {Найти наименьшее значение модуля производной} m1:=abs(f1(a)); x:=a; dx:=0.01;

repeat x:=x+dx; p:=abs(f1(x)); if p < m1 then m1:=p; until x>b; repeat {цикл уточнения корня} x:=a-f(a)*(b-a)/(f(b)-f(a)); {вычислить оценку приближения и выполнить сужение интервала изоляции корня}

if v =1 then begin d := x - a; a:=x end else begin d := b - x; b:=x; end; until d < m1/(M - m1)*e;

writeln('Корень x = ', x:12:9); end

else writeln('Метод хорд не применим'); readkey;

end.

1.2.8.5 Комбинированный метод

program PHordiNuton; uses crt;

const e = 0.00000001;

function f(x:real):real;{вычисление функции f(x)} begin f:=exp(x*ln(5))-6*x-3; end;

function f1(x:real):real; {вычисление производной функции f'(x)} begin f1:=exp(x*ln(5))*ln(5)-6; end;

function f2(x:real):real; {вторая производная функции f(x)} begin f2:=sqr(ln(5))*exp(x*ln(5)); end;

function sgn(z:real):integer; {функция вычисления знака}

begin if z=0 then sgn:=0 else if z>0 then sgn:=1 else sgn:=-1; end; function fHord(c:real):real; {вычисление приближения по методу хорд} begin fHord:=c-f(c)*(b-a)/(f(b)-f(a)); end;

function fNuton(c:real):real; {вычисление приближения по методу Ньютона} begin fNuton:=c-f(c)/f1(c); end;

var a, b : real; v:byte; begin

clrscr;

writeln('Корни отделены на отрезках [-1;0], [1;2]. Введите значения a,b'); Readln(a, b);

if sgn(f(a)) <> sgn(f(b)) then

24

begin

{Установить вариант сужения интервала изоляции корня} if sgn(f1(a)) = sgn (f2(a)) then v := 1 else v := 2;

repeat {цикл уточнения значения корня}

if v = 1 then begin a:=fHord(a); b:=fNuton(b); end

else	begin a:=fNuton(a);b:=fHord(b); end;
until b-a<e;
writeln('Корень x= ',(a+b)/2:12:9);
end
else	writeln('Комбинированный метод не применим');

readkey;

end.

1.2.9 Уточнение корней уравнений средствами Excel

Для организации циклических вычислений в Excel имеются два механизма –

циклическая ссылка и подбор параметра.

При организации циклической ссылки точность вычислений зависит от числа выполненных итераций. Это число определяет пользователь, руководствуясь формулами оценки точности вычислений приближенного значения корня, в зависимости от применяемого метода, например, в методе половинного деления число итераций, определяемого формулой (6), при этом погрешность вычислений определяется по формуле (5).

Вычислительный процесс должен быть запущен после набора всех команд, реализующих алгоритм. В связи с этим предлагается следующая схема действий:

1.В меню Сервис выбрать команду Параметры. В появившемся окне Па- раметры выбрать вкладку Вычисления, на которой снять флажок Итерации; в поле Предельное число итераций ввести число, обеспечивающее точность вычислений.

2.Выполнить ввод команд, реализующих алгоритм.

3.Выполнить команду Сервис/Параметры. На вкладке Вычисления окна диалога Параметры установить флажок Итерации и нажать кнопку Ок. Выполняется расчет.

Замечания:

1.Если Excel регистрирует циклическую ссылку, следует в появившемся окне нажать кнопку Отмена.

2.Если повторный расчет требует восстановления исходного состояния, то следует на вкладке Вычисления снять флажок Итерации, а затем выполнить двойной щелчок мышкой по соответствующей ячейке. Ячейка открывается для корректировки, формулу не корректируем, а нажимаем на клавишу Enter.

3.Расчеты с использованием циклической ссылки могут быть выполнены пошагово:

1) В меню Сервис выбираем команду Параметры. В появившемся окне Па- раметры выбираем вкладку Вычисления, на которой устанавливаем флажок Ите-

25

рации, в поле Предельное число итераций вводим число 1 и нажимаем кнопку Ok. Выполняется одна итерация.

2) При каждом нажатии клавиши F9 будет выполняться одна итерация. Рассмотрим решение примера средствами Excel.

1.2.9.1Метод половинного деления

Втаблице Excel выполним следующие действия.

1.В ячейку D4 вводим значение a = –1 – начало интервала изоляции нуля функции f(x). В ячейку E4 вводим значение b = 0 – конец интервала изоляции нуля функции f(x).

2.В ячейку G4 вводим константу ε = 10–7, определяющую точность вычис-

лений.

3.В ячейку H4 вводим формулу (5), определяющую число итераций для достижения заданной точности:

=ОКРУГЛ((LN($E$4-$D$4)-LN($G$4))/LN(2);0).

На вкладке Вычисления в поле Предельное число итераций вводим полученное число.

4. В ячейку G9 вводим формулу (4), определяющую достигнутую точность вычислений:

=($E$4-$D$4)/2^($H$7+1).

6. В ячейку H7 вводим формулу, определяющую число выполненных итера-

ций

=$H$7+1.

7. Распределяем ячейки и вводим соответствующие формулы, реализующие метод половинного деления. Данные см. в таблице.

Имя ячейки	Обозначение	Формула	Комментарий
E15	an	=ЕСЛИ($H$7=1;D4;E21)	Установить начало интер-
			вала изоляции корня
			функции f(x)
E16	bn	=ЕСЛИ($H$7=1;E4;E22)	Установить конец интер-
			вала изоляции корня
			функции f(x)
E17	x	=(E15+E16)/2	Вычислить срединную
			точку интервала – при-
			ближенное значение нуля
			функции f(x)
E18	f(x)	= 5^E17-6*E17-3	Вычислить значение
			функции f(x) в точке x.
E19	sign f(x)	=ЕСЛИ(E18>0;1;-1)	Определить знак функции
			f(x) в точке x
E20	sign f(a)	=ЕСЛИ(5^E15-6*E15>0;1;-1)	Определить знак функции
			f(x) в точке x = a
E21	an+1	=ЕСЛИ(E20=E19;E17;E15)	Вычислить новое значение

26

			начало интервала изоля-
			ции нуля функции f(x)
E22	bn+1	=ЕСЛИ(E20=E19;E16;E17)	Вычислить новое значение
			конца интервала изоляции
			нуля функции f(x)

8.Выполняем расчет по алгоритму п. 1.2.9.

1.2.9.2Метод итераций

	Метод итераций реализуется следующим алгоритмом.
1.	Выполните распределение ячеек и запись расчетных формул:
–	диапазон A5:A15 – значения аргумента: –1 ≤ x ≤ 0; x = 0,1;

–диапазон B5:B15 – значения производной f '(x). В ячейку B5 введите формулу =LN(5)*5^A5-6 и распространите ее на указанный диапазон;

–С5 – максимальное значение модуля производной max |f '(x)|. Введите формулу =ЕСЛИ(B5<0;ABS(МИН(B5:B15));МАКС(B5:B15));

–D5 – параметр λ = sgn f '(η) M2 . Введите формулу =ЗНАК(B5)*2/C5;

–E5 – предыдущее приближение x0. Введите формулу =ЕСЛИ(G5=0;A5;F5). Если число итераций G5=0, то значение для x0 устанавливается равным началу отрезка отделения корня, иначе – значению, вычисленному при предыдущей итерации;

–F5 – последующее приближение x1. Введите формулу =E5-$D$5*(5^E5- 6*E5-3);

–G5 – счетчик числа итераций. Введите формулу =G5+1;

–H5 – приближенное значение корня x. Введите формулу =(E5+F5)/2;

–I5 – значение функции в точке x. Введите формулу =5^H5-6*H5-3.

2.Выполните расчет.

3.Конец алгоритма.

1.2.9.3Метод Ньютона

Втаблице Excel выполним следующие действия.

1. Выполните распределение ячеек и запись расчетных формул: − диапазон A5:A15 – значения аргумента: –1 ≤ x ≤ 0; x = 0,1;

−диапазон B5:B15 – значения производной f '(x). В ячейку B5 введите формулу =LN(5)*5^A5-6 и распространите ее на указанный диапазон;

−диапазон C5:C15 – значения второй производной f ''(x). В ячейку C5 введите формулу =LN(5)*5^A5*LN(5) и распространите ее на указанный диапазон;

−D5 – максимальное значение модуля второй производной max |f ''(x)|. Введи-

те формулу =ЕСЛИ(C5<0;ABS(МИН(C5:C15));МАКС(C5:C15));

−E5 – минимальное значение модуля производной min |f '(x)|. Введите форму-

лу=ЕСЛИ(B5<0;ABS(МАКС(B5:B15));МИН(B5:B15));

27

−F5 – предыдущее приближение x0. Введите формулу

=ЕСЛИ(H5=0;ЕСЛИ(ЗНАК(B5)=ЗНАК(C5);A15;A5);G5). {Если число итераций

H5=0, то значение для x0 устанавливается равным началу отрезка отделения корня, если не совпадают знаки первой и второй производной функции и концу отрезка отделения корня, если эти знаки совпадают, иначе – значению, вычисленному при предыдущей итерации};

−G5 – последующее приближение x1. Введите формулу

=F5-(5^F5-6*F5-3)/(LN(5)*5^F5-6);

−H5 – счетчик числа выполненных итераций. Введите формулу =H5+1;

−I5 – значение функции в точке x. Введите формулу =5^G5-6*G5-3;

−J5 – значение d. Введите формулу =(G5-F5)^2;

−K5 – погрешность вычислений. Введите формулу =2*E5*0,0000001/D5.

2.Выполните расчет.

3.Конец алгоритма.

1.2.9.4Метод хорд

Втаблице Excel выполним следующие действия.

1. Выполните распределение ячеек и запись расчетных формул: − диапазон A5:A15 – значения аргумента: –1 ≤ x ≤ 0; x = 0,1;

−диапазон B5:B15 – значения производной f '(x). В ячейку B5 введите формулу =LN(5)*5^A5-6 и распространите ее на указанный диапазон;

−диапазон C5:C15 – значения второй производной f ''(x). В ячейку C5 введите формулу =LN(5)*5^A5*LN(5) и распространите ее на указанный диапазон;

−D5 – максимальное значение модуля второй производной max |f ''(x)|. Введи-

те формулу =ЕСЛИ(C5<0;ABS(МИН(C5:C15));МАКС(C5:C15));

−E5 – минимальное значение модуля производной min |f '(x)|. Введите форму-

лу=ЕСЛИ(B5<0;ABS(МАКС(B5:B15));МИН(B5:B15));

−F5 – определить вариант v расчета. Введите формулу =ЕСЛИ(ЗНАК(B5)=ЗНАК(C5);1;2). Если совпадают знаки первой и второй производной функции f(x), то принимаем v =1 – расчет приблизительного значения корня с недостатком, иначе v = 2 – расчет приблизительного значения корня с избытком.

−G5 –начало интервала изоляции корня, a. Введите формулу

=ЕСЛИ(J5=0;A5;ЕСЛИ(F5=1; I5;G5)).

Если число итераций J5=0, то значение для a устанавливается равным началу отрезка отделения корня, иначе, если вариант расчета v = 1 – значению x, в других случаях значение a не изменяется;

−H5 – конец интервала изоляции корня, b. Введите формулу

=ЕСЛИ(J5=0;A15;ЕСЛИ(F5=2; I5;H5)).

Если число итераций J5=0, то значение для a устанавливается равным началу отрезка отделения корня, иначе, если вариант расчета v = 2 – значению x, в других случаях значение b не изменяется;

−I5 – приближенное значение корня. Введите формулу

28

=G5-(5^G5-6*G5-3)*(H5-G5)/(5^H5-6*H5-3-(5^G5-6*G5-3));

−J5 – счетчик числа выполненных итераций. Введите формулу =J5+1;

−K5 – значение функции в точке x. Введите формулу =5^I5-6*I5-3;

−L5 – значение d. Введите формулу =ЕСЛИ(F5=1;I5-G5;H5-I5);

−M5 – погрешность вычислений. Введите формулу

=E5*0,0000001/ABS((D5-E5)).

2.Выполните расчет.

3.Конец алгоритма.

1.2.9.5Комбинированный метод

Втаблице Excel выполним следующие действия.

1.Выполните распределение ячеек и запись расчетных формул:

−диапазон A5:A6 – значения аргумента: a = –1; b = 0;

−диапазон B5:B6 – значения производной f '(x). В ячейку B5 введите формулу =LN(5)*5^A5-6 и распространите ее на указанный диапазон;

−диапазон C5:C6 – значения второй производной f ''(x). В ячейку C5 введите формулу =LN(5)*5^A5*LN(5) и распространите ее на указанный диапазон;

−D5 – определить вариант расчета v. Введите формулу =ЕСЛИ(ЗНАК(B5)=ЗНАК(C5);1;2). Если совпадают знаки первой и второй производной функции f(x), то принимаем v =1, иначе v = 2;

−E5 – начало интервала изоляции корня, a. Введите формулу

=ЕСЛИ(H5=0;A5;ЕСЛИ(D5=1; E5-(5^E5-6*E5-3)*(F5-E5)/(5^F5-6*F5-3-(5^E5- 6*E5-3));E5-(5^E5-6*E5-3)/(LN(5)*5^E5-6))).

Если число итераций H5=0, то значение для a устанавливается равным началу отрезка отделения корня, иначе, если вариант расчета v = 1 – значение a рассчитывается по формуле хорд, в других случаях значение a рассчитывается по формуле Ньютона;

−F5 – b, конец интервала изоляции корня. Введите формулу

=ЕСЛИ(H5=0;A6;ЕСЛИ(D5=2; F5-(5^F5-6*F5-3)*(F5-E5)/(5^F5-6*F5-3-(5^E5- 6*E5-3));F5-(5^F5-6*F5-3)/(LN(5)*5^F5-6))).

Если число итераций H5=0, то значение для b устанавливается равным концу отрезка отделения корня, иначе, если вариант расчета v = 2 – значение b рассчитывается по формуле хорд, в других случаях значение b рассчитывается по формуле Ньютона;

−G5 – приближенное значение корня. Введите формулу =(E5+F5)/2;

−H5 – счетчик числа выполненных итераций. Введите формулу =H5+1;

−I5 – значение функции в точке x. Введите формулу =5^G5-6*G5-3;

−J5 – значение d – погрешность вычислений. Введите формулу

=F5-E5.

2.Выполните расчет.

3.Конец алгоритма.

29

1.2.9.6 Использование команды "Подбор параметра"

Рассмотрим решение примера с использованием команды подбор параметра из меню сервис.

1. В ячейку А5 вводим формулу оценки погрешности

=ABS(СТЕПЕНЬ(5;А6)-6*А6-3)/A4, где A4 = min |f ′(x)|.

x [α ; β ]

2.В ячейку А6 вводим приближенное значение корня из интервала его изоляции, например, –0,5 для отрезка [–1; 0].

3.Делаем ячейку А5 активной. В меню сервис выбираем команду подбор па-

раметра. В окне диалога подбор параметра в поле установить в ячейке уста-

навливаем значение $A$5. В поле значение вводим число, определяющее точность вычислений, например 0,00000001. В поле изменяя значение ячейки вводим ссылку на ячейку А6. Нажимаем кнопку ОК. EXCEL выводит окно результат подбора параметра. Ответ в смотрите в ячейке А6. Аналогично получаем приближенное значение корня на отрезке [1; 2].

Замечание. Команда подбор параметра не всегда обеспечивает требуемую точность вычислений. Можно приближать также значение функции к нулю – иногда это повышает точность вычислений.

1.2.10Решение уравнений средствами MathCAD

Для численного решения нелинейного уравнения f(x)=0 MathCAD имеет встроенную функцию root, которая имеет вид

root(f(х),x,a,b),

где f (х) – левая часть уравнение, х — скалярная переменная, относительно которой решается уравнение, а, b — границы интервала, на котором находится корень уравнения (необязательные параметры).

Поиск корня уравнения осуществляется итерационным методом с заданной точностью (точность по умолчанию 10-3 ; системная переменная TOL отвечает за точность).

В системе MathCAD выполните следующие действия:

1.Введите уравнение f(x):= 5x − 6 x − 3 ;

2.Введите заданную точность TOL:= 10-7 ;

3.С помощью функции root уточните корень на отрезке [-1;0] root( f (x) , x ,−1 ,0) = −0.4144636 ;

4.Уточните корень на отрезке [1;2]

root( f (x) , x ,1 ,2) = 1.5633283 .

Задания к данной теме приведены в приложении А.

30

2.РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ

2.1Постановка задачи

Пусть дана система n линейных уравнений с n неизвестными.

a x + a x

2

+ ... + a

x

n

= b , ü

11 1

12

1n

1

ï

a21x1 + a22 x2

+ ... + a2n xn

= b2

,ï

(1)

........................................

ý

ï

a

x + a

n2

x

2

+ ... + a

nn

x

n

= b .ï

n1 1

n þ

Совокупность чисел x1, x2, …,xn, обращающих систему (1) в тождество, называется решением этой системы.

Методы решения систем линейных уравнений можно разделить на две груп-

пы:

1)точные методы, представляющие собой конечные алгоритмы для вычисления корней системы. Таковы, например, правило Крамера, метод Гаусса, метод главных элементов, метод квадратных корней и другие;

2)итерационные методы, позволяющие получить корни системы с заданной точностью путем сходящихся бесконечных процессов. К числу их относятся метод итерации, метод Зейделя, метод релаксации и другие.

Вследствие неизбежных округлений результаты даже точных методов являются приближенными, причем оценка погрешностей корней в общем случае затруднительна. При использовании итерационных методов, сверх того, добавляется погрешность метода.

2.2Метод Гаусса

Метод Гаусса распадается на 2 этапа.

Этап 1 носит название прямого хода и представляет собой процесс исключения неизвестных, в результате которого заданная система уравнений (1) приводится к эквивалентной системе (4) треугольного вида за n – 1 шагов. Пусть k – номер шага, k = 1, 2, ..., n–1. На шаге k выполняется исключение неизвестной xk из уравнений системы (1), начиная с k+1 и последующих уравнений. Таким образом, если i – номер уравнения, из которого исключается неизвестная xk, то i = k+1, .., n.

Принято k-е уравнение называть разрешающим, а элемент ak k – разрешаю- щим элементом. Необходимое требование метода akk¹0. В противном случае можно поменять местами k-е уравнение с i-м уравнением (i > k), в котором коэффициент при неизвестной xk отличен от нуля.

Обозначим через ai(kj) преобразованный на k шаге коэффициент в i уравнении при xj неизвестной; bi(k ) – преобразованный на k шаге свободный член i уравнения.

Исключение неизвестной xk из i уравнения системы (1) выполняется путем последовательного вычитания из каждого i уравнения (i = k + 1, ..., n) разрешающего уравнения, умноженного на соответствующий коэффициент ai k и деленного