Численные методы решения нелинейных уравнений

Геометрический смысл метода касательных  состоит в замене дуги y = f (x) касательной, одной к одной из крайних точек . Начальное приближение x ₀=a или    x₀ = b брать таким, чтобы вся последовательность приближения х _k принадлежала интервалу (a;b). В случае существования производных f ’, f ”, сохраняющих свои знаки в интервале, за х₀ берется тот конец отрезка [a;b], для которого выполняется условие  f ’(х₀) * f (х₀) > 0. Для оценки приближения используется общая формула :

|c-x _k-1 | £ | f (x _k+1)/m| , где m = _min f ’(x) на отрезке [a;b] .

На практике проще пользоваться другим правилом :

Если на отрезке [a;b] выполняется условие 0 < m <  | f (x)|  и e - заданная точность решения, то неравенство | x _k+1-x _k| £ e влечет выполнение неравенства       |c-x _k-1| £ e .

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

|c-x _k-1| £ e .

Пример 1. Определим корни уравнения х³ + 0,1х² + 0,4х – 1,2 = 0 аналитически. Находим: f  (x) =  х³ + 0,1х² + 0,4х – 1,2

f ‘ (x) =  3х² + 0,1х + 0,4

f (–1)   = –2,5 < 0  f (0)   = –1,2 < 0  f (+1)   = 0,3 > 0

x	- ¥	-1	0	+1	+ ¥
sign f (x)	-	-	-	+	+

Следовательно, уравнение имеет  действительный корень, лежащий в  промежутке [0; +1].

Приведем уравнение к виду x =  j (x), так, чтобы |j ‘ (x)|<1 при 0 £ x £ +1.

Так как max | f ’(x) | = f ’(+1) = 3 + 0,1 + 0,4 = 3,5 то можно взять R = 2.

Тогда j (x) = x – ( f (x) / R) = x – 0,5 х³ – 0,05 х² – 0,2 х + 0,6 = – 0,5 х³ – 0,05 х² + 0,8 х + 0,6.

Пусть х₀ = 0 , тогда х _n+1 = j (х _n).

Вычисления расположим в таблице.     

n	хn	х2n	х3n	j (хn).	f (x)
1	1	1	1	0,85	-0,17363
2	0,85	0,7225	0,614125	0,9368125	0,08465
3	0,9368125	0,87761766	0,822163194	0,89448752	-0,04651
4	0,89448752	0,800107923	0,715686552	0,917741344	0,024288
5	0,917741344	0,842249174	0,772966889	0,905597172	-0,01306
6	0,905597172	0,820106238	0,74268589	0,912129481	0,006923
7	0,912129481	0,83198019	0,758873659	0,908667746	-0,0037
8	0,908667746	0,825677072	0,750266124	0,910517281	0,001968
9	0,910517281	0,829041719	0,754856812	0,909533333	-0,00105
10	0,909533333	0,827250884	0,752412253	0,910057995	0,000559
11	0,910057995	0,828205555	0,753715087	0,909778575	-0,0003
12	0,909778575	0,827697055	0,753021048	0,909927483	0,000159
13	0,909927483	0,827968025	0,753390861	0,909848155	-8,5E-05
14	0,909848155	0,827823665	0,753193834	0,909890424	4,5E-05
15	0,909890424	0,827900583	0,753298812	0,909867904	-2,4E-05
16	0,909867904	0,827859602	0,753242881	0,909879902	1,28E-05
17	0,909879902	0,827881437	0,753272681	0,90987351	-6,8E-06
18	0,90987351	0,827869803	0,753256804	0,909876916	3,63E-06
19	0,909876916	0,827876002	0,753265263	0,909875101	-1,9E-06
20	0,909875101	0,827872699	0,753260756	0,909876068	1,03E-06

 

График функции  y =  х³ + 0,1х² + 0,4х – 1,2

 

 

 

 

Блок-схема  программы

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Программа на языке Paskal

program metod_kasatel;

var  xn,xn1,a,b,c,mx,y0,x0:real;

function f1(x1:Real): Real;

begin

  f1 := x1*x1*x1*(-0.5)-0.05*x1*x1+0.8*x1+0.6;

 end;

function f2(x4:Real): Real; {Производная от основной функции}

 begin

  f2 := x4*x4*x4+0.5*x4*x4+0.1*x4*x4+0.4*x4–1.2;

 end;

begin  

a:=0;b:=1;c:=0.00000001;

 Writeln(' От A=',a,' до B=',b);

Writeln(' Погрешность с=',c);

Readln;

xn:=b;

xn1:= f1(xn);

 y0:=f2(b);

while ABS(y0)>c do {Проверка по точности вычисления корня}

 begin   xn:=xn1;

  xn1:=f1(xn);

  y0:= f2(xn1);

   {Печать  промежуточного результата}

   Writeln('xn=',xn,' xn+1=',xn1,' f(xn+1)=',y0);

  Readln;

end;

Writeln('Конечные значения'); {Печать полученного результата}

Writeln(' xn+1=',xn1,' f(xn+1)=',y0);

end.

 

Результат выполнения программы

От A= 0.0000000000E+00 до B= 1.0000000000E+00

 Погрешность с= 1.0000000000E-08

От A= 0.0000000000E+00 до B= 1.0000000000E+00

 Погрешность с= 1.0000000000E-08

 

xn= 8.5000000000E-01 xn+1= 9.3681250000E-01 f(xn+1)= 8.4649960270E-02

xn= 9.3681250000E-01 xn+1= 8.9448751986E-01 f(xn+1)=-4.6507647892E-02

xn= 8.9448751986E-01 xn+1= 9.1774134381E-01 f(xn+1)= 2.4288343840E-02

xn= 9.1774134381E-01 xn+1= 9.0559717189E-01 f(xn+1)=-1.3064617920E-02

xn= 9.0559717189E-01 xn+1= 9.1212948085E-01 f(xn+1)= 6.9234699658E-03

xn= 9.1212948085E-01 xn+1= 9.0866774587E-01 f(xn+1)=-3.6990702320E-03

xn= 9.0866774587E-01 xn+1= 9.1051728099E-01 f(xn+1)= 1.9678960780E-03

xn= 9.1051728099E-01 xn+1= 9.0953333295E-01 f(xn+1)=-1.0493249720E-03

xn= 9.0953333295E-01 xn+1= 9.1005799543E-01 f(xn+1)= 5.5884091853E-04

xn= 9.1005799543E-01 xn+1= 9.0977857497E-01 f(xn+1)=-2.9781681224E-04

xn= 9.0977857497E-01 xn+1= 9.0992748338E-01 f(xn+1)= 1.5865717614E-04

xn= 9.0992748338E-01 xn+1= 9.0984815480E-01 f(xn+1)=-8.4537703515E-05

xn= 9.0984815480E-01 xn+1= 9.0989042365E-01 f(xn+1)= 4.5040009354E-05

xn= 9.0989042365E-01 xn+1= 9.0986790364E-01 f(xn+1)=-2.3997676180E-05

xn= 9.0986790364E-01 xn+1= 9.0987990248E-01 f(xn+1)= 1.2785800209E-05

xn= 9.0987990248E-01 xn+1= 9.0987350958E-01 f(xn+1)=-6.8122881203E-06

xn= 9.0987350958E-01 xn+1= 9.0987691573E-01 f(xn+1)= 3.6295678001E-06

xn= 9.0987691573E-01 xn+1= 9.0987510095E-01 f(xn+1)=-1.9338276616E-06

xn= 9.0987510095E-01 xn+1= 9.0987606786E-01 f(xn+1)= 1.0303429008E-06

xn= 9.0987606786E-01 xn+1= 9.0987555269E-01 f(xn+1)=-5.4896190704E-07

xn= 9.0987555269E-01 xn+1= 9.0987582717E-01 f(xn+1)= 2.9248803912E-07

xn= 9.0987582717E-01 xn+1= 9.0987568093E-01 f(xn+1)=-1.5583464119E-07

xn= 9.0987568093E-01 xn+1= 9.0987575885E-01 f(xn+1)= 8.3031409304E-08

xn= 9.0987575885E-01 xn+1= 9.0987571733E-01 f(xn+1)=-4.4236003305E-08

xn= 9.0987571733E-01 xn+1= 9.0987573945E-01 f(xn+1)= 2.3572283681E-08

xn= 9.0987573945E-01 xn+1= 9.0987572766E-01 f(xn+1)=-1.2558302842E-08

xn= 9.0987572766E-01 xn+1= 9.0987573394E-01 f(xn+1)= 6.6920620156E-09

Конечные значения

xn+1= 9.0987573394E-01 f(xn+1)= 6.6920620156E-09

Метод хорд.

Пусть на отрезке [a;b] функция непрерывна, принимает на концах отрезка значение разных знаков, а производная f '(x) сохраняет знак. В зависимости от знака второй производной возможны следующие случаи расположения кривых (рис. 1).

 
 

Рис. 1. Возможные случаи расположения кривых.

Алгоритм приближенного вычисления корня методом хорд. 
Исходные данные: f (x) – функция; ε – требуемая точность; x₀ – начальное приближение. 
Результат: xпр – приближенный корень уравнения f (x) = 0.

Метод решения: 

Рис. 2. Геометрическая интерпретация  метода хорд для случая f '(x) f ''(x)>0.

Рассмотрим случай, когда f '(x)  и f ''(x)  имеют одинаковые знаки (рис. 2).

График функции проходит через  точки A₀(a,f(a)) и B₀(b,f(b)). Искомый корень уравнения (точка x*) нам неизвестен, вместо него возьмет точку х₁ пересечения хорды А₀В₀ с осью абсцисс. Это и будет приближенное значение корня. 
В аналитической геометрии выводится формула, задающая уравнение прямой, проходящей через две точки с координатами (х1; у1) и (х2; у2): . 
Тогда уравнение хорды А₀В₀ запишется в виде: . 
Найдем значение х = х₁, для которого у = 0:  . Теперь корень находится на отрезке [x₁;b]. Применим метод хорд к этому отрезку. Проведем хорду, соединяющую точки A₁(x₁,f(x₁)) и B₀(b,f(b)), и найдем х₂ - точку пересечения хорды А₁В₀ с осью Ох: x₂=x₁ . 
Продолжая этот процесс, находим: x₃=x₂. Получаем рекуррентную формулу вычисления приближений к корню x_n+1=x_n. 
В этом случае конец b отрезка [a;b]  остается неподвижным, а конец a перемещается. 
Таким образом, получаем расчетные формулы метода хорд:

x_n+1=x_n;         x₀=a.                                     (4)

Вычисления очередных приближений  к точному корню уравнения  продолжается до тех пор, пока не достигнем  заданной точности, т.е. должно выполняться  условие: |x_n+1-x_n|<, где - заданная точность. 
Теперь рассмотрим случай, когда первая и вторая производные имеют разные знаки, т.е. f '(x) f ''(x)<0. (рис. 3).

 
Рис. 3. Геометрическая интерпретация метода хорд для случая f '(x) f ''(x)<0.

Соединим точки A₀(a,f(a)) и B₀(b,f(b))  хордой А₀В₀. Точку пересечения хорды с осью Ох будем считать первым приближение корня. В этом случае неподвижным концом отрезка будет являться конец а. 
Уравнение хорды А₀В₀:. Отсюда найдем x₁, полагая y = 0: x₁=b. Теперь корень уравнения x[a;x₁]. Применяя метод хорд к этому отрезку, получим x₂=x₁. Продолжая и т.д., получим x_n+1=x_n. 
Расчетные формулы метода:

x_n+1=x_n, x₀=0 .                                                           (5) 
Условие окончания вычислений: |x_n+1-x_n|<. Тогда хпр = xn+1 с точностью Итак, если f '(x) f ''(x)>0 приближенное значение корня находят по формуле (4), если f '(x) f ''(x)<0, то по формуле (5).

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

Пример. Проиллюстрировать действие этого правила на уравнении

(x-1)ln(x)-1=0, если отрезок изоляции корня [2;3].

Решение. Здесь f(x)=(x-1)ln(x)-1.

f '(x)=ln(x)+;

 f ''(x)=.

Вторая производная в этом примере  положительна на отрезке изоляции корня [2;3]: f ''(x)>0, f(3)>0, т.е. f(b) f''(x)>0. Таким образом, при решении данного уравнения методом хорд для уточнения корня выбираем формулы (4).

 

 

program horda; 
var e,c,a,b,y,ya,yb,yn,x,x1,x2,xn,f1,f2:real; 
begin e:=0.0001; 
writeln('vvedi nachalo otrezka'); 
readln(a); 
writeln('vvedi konec otrezka'); 
readln(b); 
x:=a; 
y:=((x-1)*ln(x))-1; 
ya:=y; x:=b; 
y:=((x-1)*ln(x))-1; 
yb:=y; c:=(a+b)/2; x:=c; 
y:=((x-1)*ln(x))-1; 
f1:=ln(x) + (x-1)/x ; 
f2:= 1/x + 1/(x*x);

if (ya*yb < 0) and (f1*f2 > 0) 
then begin x1:=a; while abs(x2 - x) > e do 
begin x:=x1; y:=((x-1)*ln(x))-1; yn:=y; 
x2:=x1 - (yn*(b-x1))/(yb - yn); 
x1:=x2; end; 
writeln('koren uravneniya xn = ', x2) 
end elsebegin x1:=b; 
while abs(x2 - x) > e do 
begin x:=x1; y:=((x-1)*ln(x))-1; yn:=y; 
x2:=x1 - (yn*(x1- a))/(yn - ya); 
x1:=x2; end; 
writeln('koren uravneniya xn = ', x2); 
end; 
end. 

 

 

Метод простых  итераций.

Рассмотрим  уравнение f(x)=0  (1) с отделенным корнем X[a, b]. Для решения уравнения (1) методом простой итерации приведем его к равносильному виду: x=φ(x).  (2)

Это всегда можно  сделать, причем многими способами. Например:

x=g(x) · f(x) + x ≡ φ(x), где g(x) - произвольная непрерывная функция, не имеющая корней на отрезке [a,b].

Пусть x⁽⁰⁾ - полученное каким-либо способом приближение к корню x (в простейшем случае x⁽⁰⁾=(a+b)/2). Метод простой итерации заключается в последовательном вычислении членов итерационной последовательности:

x^(k+1)=φ(x^(k)), k=0, 1, 2, ... (3)

начиная с приближения x⁽⁰⁾.

УТВЕРЖДЕНИЕ: 1 Если последовательность {x^(k)} метода простой итерации сходится и функция φ непрерывна, то предел последовательности является корнем уравнения x=φ(x)

ДОКАЗАТЕЛЬСТВО: Пусть . (4)

Перейдем  к пределу в равенстве x^(k+1)=φ(x^(k)) Получим с одной стороны по (4), что а с другой стороны в силу непрерывности функции φ и (4) .

В результате получаем x^*=φ(x^*). Следовательно, x^* - корень уравнения (2), т.е. X=x^*.

Чтобы пользоваться этим утверждением нужна сходимость последовательности {x^(k)}. Достаточное условие сходимости дает:

ТЕОРЕМА 1: (о сходимости) Пусть уравнение x=φ(x) имеет единственный корень на отрезке [a,b] и выполнены условия:

1) φ(x) C¹[a,b];

2) φ(x) [a,b] " x [a,b];

 

3) существует константа q > 0: | φ '(x) | ≤ q < 1 x [a,b]. Tогда итерационная последовательность {x^(k)}, заданная формулой x^(k+1) = φ(x^(k)), k=0, 1, ... сходится при любом начальном приближении x⁽⁰⁾ [a,b].

ДОКАЗАТЕЛЬСТВО: Рассмотрим два соседних члена последовательности {x^(k)}: x^(k) = φ(x^(k-1)) и x^(k+1) = φ(x^(k)) Tак как по условию 2) x^(k) и x^(k+1) лежат внутри отрезка [a,b], то используя теорему Лагранжа о средних значениях получаем:

x ^(k+1) - x ^(k) = φ(x ^(k)) - φ(x ^(k-1)) = φ '(c _k )(x ^(k) - x ^(k-1)), где c _k (x ^(k-1), x ^(k)).

Отсюда получаем:

| x ^(k+1) - x ^(k) | = | φ '(c _k ) | · | x ^(k) - x ^(k-1) | ≤ q | x ^(k) - x ^(k-1)| ≤

≤ q ( q | x ^(k-1) - x ^(k-2) | ) = q ² | x ^(k-1) - x ^(k-2) | ≤ ... ≤ q ^k | x ⁽¹⁾ - x ⁽⁰⁾ |. (5)

Рассмотрим  ряд S_∞ = x ⁽⁰⁾ + ( x ⁽¹⁾ - x ⁽⁰⁾ ) + ... + ( x ^(k+1) - x ^(k) ) + ... . (6)

Если мы докажем, что этот ряд сходится, то значит сходится и последовательность его  частичных сумм

S_k = x ⁽⁰⁾ + ( x ⁽¹⁾ - x ⁽⁰⁾ ) + ... + ( x ^(k) - x ^(k-1) ).

Но нетрудно вычислить, что

S_k = x ^(k)). (7)

Следовательно, мы тем самым докажем и сходимость итерационной последовательности {x^(k)}.

Для доказательства сходимости pяда (6) сравним его почленно (без первого слагаемого x⁽⁰⁾) с рядом

q ⁰ | x ⁽¹⁾ - x ⁽⁰⁾ | + q ¹ |x ⁽¹⁾ - x ⁽⁰⁾| + ... + |x ⁽¹⁾ - x ⁽⁰⁾| + ...,  (8)

который сходится как бесконечно убывающая геометрическая прогрессия (так как по условию q < 1). В силу неравенства (5) абсолютные величины ряда (6) не превосходят соответствующих членов сходящегося ряда (8) (то есть ряд (8) мажорирует ряд (6). Следовательно ряд (6) также сходится. Tем самым сходится последовательность {x⁽⁰⁾}.

Получим формулу, дающую способ оценки погрешности |X - x ^(k+1)|

метода простой  итерации.

Имеем

X - x^(k+1) = X - S_k+1 = S_∞ - S_k+1 = (x^(k+2) - ^(k+1) ) + (x^(k+3) - x^(k+2) ) + ... .

Следовательно

|X - x^(k+1)| ≤ |x^(k+2) - ^(k+1) | + |x^(k+3) - x^(k+2) | + ... ≤ q^k+1 |x⁽¹⁾ - x⁽⁰⁾ | + q^k+2 |x⁽¹⁾ - x⁽⁰⁾ | + ... = q^k+1|x⁽¹⁾ - x⁽⁰⁾ | / (1-q).

В результате получаем формулу

|X - x^(k+1)| ≤ q^k+1|x⁽¹⁾ - x⁽⁰⁾ | / (1-q).   (9)

Взяв за x⁽⁰⁾ значение x^(k), за x⁽¹⁾ - значение x^(k+1) (так как при выполнении условий теоремы такой выбор возможен) и учитывая, что при имеет место неравенство q^k+1 ≤ q выводим:

|X - x^(k+1)| ≤ q^k+1|x^(k+1) - x^(k) | / (1-q) ≤ q|x^(k+1) - x^(k) | / (1-q).

Итак, окончательно получаем:

|X - x^(k+1)| ≤ q|x^(k+1) - x^(k) | / (1-q).   (10)

Используем  эту формулу для вывода критерия окончания итерационной последовательности. Пусть уравнение x=φ(x) решается методом простой итерации, причем ответ должен быть найден с точностью ε, то есть

|X - x^(k+1)| ≤ ε.

С учетом (10) получаем, что точность ε будет достигнута, если выполнено неравенство