Участник:Gukov/Песочница

Материал из MachineLearning.

(Различия между версиями)

Перейти к: навигация, поиск

Текущая версия

Содержание

1 Введение
- 1.1 Метод Ньютона-Рафсона
- 1.2 Случай почти равных корней
2 Изложение метода
3 Числовой пример
4 Рекомендации программисту
5 Заключение
6 Список литературы

Введение

Метод Ньютона-Рафсона

Основная статья: Метод касательных (Ньютона-Рафсона)

Пусть имеется некоторая функция $F(x)$ и необходимо найти такие значения аргумента $x$ , для которых

(1)

$F(x)=0$

Перепишем (1) в виде

(1.1)

$x = f(x)$

и запишем $n+1$ -ое приближения корня (1.1), при этом делая поправку $\alpha$ к очередному значению $x_n$

(2)

$x_{n+1} = x_n + \alpha \Delta x_n,$

где $\Delta x_n = f(x_n) - x_n,$ и положим $\alpha = \frac{1}{1-f'(x_n)}$ .

Тогда (2) перепишется в виде

(3)

$x_{n+1} = \frac{f(x_n) - x_n f'(x_n)}{1 - f'(x_n)}$

Нетрудно видеть, что (3) эквивалентно простому методу последовательных приближений

$x_{n+1} = g(x_n),$

где $g(x) = \frac{f(x) - x f'(x)}{1 - f'(x)}$ .

Вспомним также, что если $|g'(x)| < 1$ , то метод сходится. Имеем

$g'(x) = \frac{f''(x)[f(x) - x]}{[1 - f'(x)]^2}.$

Но так как справедливо соотношение (1.1), то для $x,$ достаточно близких к решению (1.1), выражение в скобках в числителе дроби становится малым. Поэтому итерационный метод, описываемый формулой (3) сходится, если

1. Начальное приближение $x_0$ выбрано достаточно близко к решению $x = f(x)$ .
2. Производная $f''(x)$ не становится слишком большой.
3. Производная $f'(x)$ не слишком близка к 1.

Это и есть метод Ньютона-Рафсона. Обычно его записывают в виде

(4)

$x_{n+1} = x_{n} - \frac{F(x_n)}{F'(x_n)}$ ,

где $F(x) = f(x) - x = 0$

Таким образом, мы вернулись к уравнению в форме (1), и условия сходимости принимают следующий вид

1. Начальное приближение $x_0$ выбрано достаточно близко к корню уравнения $F(x) = 0$ .
2. Производная $F''(x)$ не становится очень большой.
3. Производная $F'(x)$ не слишком близка к 0.

Случай почти равных корней

Рисунок 1

Условие 3 сходимости метода Ньютона-Рафсона означает, что никакие два корня не находятся слишком близко один к другому. Соответствующая ситуация представлена на рисунке 1 (масштаб сильно увеличен). Заметим, что производная $f'(x)$ близка к 1 при x, равном обоим значениям корней, $a_1$ и $a_2$ . Более того, на основании теоремы Лагранжа, можно утверждать, что $f'(x) = 1$ где-то между $a_1$ и $a_2$ .

Рассмотрим, что случится, если принять $x_0$ в качестве исходного значения для корня $a_1$ . Касательная, проведенная через точку $C$ , пересечет прямую $y=x$ в точке $A$ , и следующее приближение будет равно $x_1$ . Касательная, проведенная через точку $B$ , пересекает прямую в точке $D$ , и в качестве следующего приближеня получается снова $x_0$ . Итерационный процесс, таким образом, осциллирует между $x_0$ и $x_1$ до бесконечности, не сходясь ни к одному значению корня. Иначе говоря, не удается отделить эти два корня, потому что они расположены слишком близко.

Поэтому необходимо начальное приближение, достаточно близкое к искомому значению корня. Трудности возникают потому, что вычисление знаменателя в формуле (3) включает в себя вычитание двух почти равных чисел, что приводит к понижению точности.

Изложение метода

Поиск начального приближения

Рисунок 2

Сначала находится значение x, при котором $f'(x)=1$ , то есть решается уравнение

$x=x+f'(x)-1$

Пусть решением этого уравнения будет некоторое $x = x_0$ . Эта точка расположена между двумя корнями, $a_1$ и $a_2$ . Чтобы получить начальное приближение для решения уравнения, предположим, что $x_0$ лежит посредине между $a_1$ и $a_2$ (рисунок 2). Другими словами, мы предполагаем, что $x_0 + d$ и $x_0 - d$ являются корнями уравнения (1.1). Разлагая $f(x)$ в ряд Тейлора в окрестности точки $x_0$ и принимая во внимание, что $f'(x)=1$ , получаем

$f(x) = f(x_0) + (x - x_0) + \frac{1}{2} f''(x_0) (x - x_0)^2 + \ldots$

Ограничим ряд тремя членами. Подставляя $x + d$ вместо $x$ , имеем

$f(x_0+d) = f(x_0) + d + \frac{1}{2} f''(x_0)(d^2)$

Но по условию

$f(x_0+d) = x_0+d,$

поэтому, решая эти уравнения относительно $d$ , получаем

$d=\sqrt{\frac{2(x_0 - f(x_0))}{f''(x_0)}}$

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

$x = x + f'(x) - 1$

и определить значение $x_0$ . Для решения этого уравнения можно применить, например, метод Ньютона-Рафсона. Найдя значение $x_0$ , можно определить значение $d$ . И, наконец, значения $x_0-d$ и $x_0+d$ используются в качестве начальных приближений для определения соответственно $a_1$ и $a_2$ .

Метод Рунге

Пусть существует асимптотическое разложение вида:

$D_N(f) = I_N(f) - I_0 = \alpha _2 h^2 + \alpha _4 h^4 + \alpha _6 h^6 + \ldots$ ,

где $f(x)$ - достаточно гладкая функция, а $hN = b - a$ .

Проведем расчеты на двух равномерных сетках с шагами $h_1$ и $h_2$ соответственно и найдем выражения $I^{h_1}[f] = I_{N_1}[f]$ и $I^{h_2}[f] = I_{N_2}[f]$ , $h_1 N_1 = h_2 N_2 = b - a$ . Потребуем, чтобы погрешность для их линейной комбинации:

${\tilde D}^{h}(f) = \sigma D^{h_1}(f) + (1 - \sigma)D^{h_2}(f)$

была величиной более высокого порядка по сравнению с $D^{h_1}$ и $D^{h_2}$ . Если для $D^h = D_N$ имеет место формула вида

$D^{h} = I^{h}[f] - I[f] = \alpha _p h^p + \alpha _q h^q + \ldots,$ $q > p$ ,

то для

${\tilde D}^{h} = (\sigma I^{h_1}[f] + (1 - \sigma)I^{h_2}[f]) - I[f]$

получим

${\tilde D}^{h}(f) = \alpha _p(\sigma h_1^{p} + (1 - \sigma)h_2^{p}) + \alpha _q (\sigma h_1^{q} + (1 - \sigma)h_2^{q}) + \ldots$

Выберем параметр \sigma из условия $\sigma h_1^{p} + (1 - \sigma)h_2^{p}$ :

$\sigma = \frac{h_2^p}{h_2^p - h_1^p}$ .

Имеем

${\tilde D}^h(f) = \alpha _q(\sigma h_1^q + (1 - \sigma)h_2^q) + \ldots = O(h^q)$ , $h = max(h_1, h_2)$ ,

причем $\sigma h_1^q + (1 - \sigma)h^q_2 < 0$ . Так, если $p = 2, q = 4$ , то ${\tilde D}^{h}(f) = - \alpha _4 h_1^2 h_2^2 + \ldots = O(h^4)$ . Таким образом, проведя вычисления на двух сетках с шагами $h_1$ и $h_2 \ne h_1$ , мы повысили порядок точности на $2$ (на $(q - p)$ ) для ${\tilde I} = \sigma I^{h_1} + (1 - \sigma)I^{h_2}$ .

Процесс Эйткена

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

$D^{h}(f) = \alpha _p h^{p} + \alpha _q h^{q} + \ldots$ , $q > p$ ,

так что

$I^{h}[f] = I[f] + \alpha _p h^p + \alpha _q h^{q} + \ldots$

Проведем вычисления на трех сетках: $h_1 = h$ , $h_2 = \rho h$ , $h_3 = \rho^{2} h$ ( $0 < \rho < 1$ ). Определим $p$ . При этом пренебрегаем членом $O(h^q)$ . Образуем отношение

$A = \frac{I^{h_1}[f] - I^{h_2}[f]}{I^{h_2}[f] - I^{h_3}[f]} \approx \frac{h^p_1 - h^p_2}{h^p_2 - h^p_3} = \frac{1 - \rho^p}{\rho^p(1-\rho^p)} = (\frac{1}{\rho})^p$

и найдем

$p \approx \frac{\ln A}{\ln(\frac{1}{\rho})}$ .

Зная приближенное значение $p$ , можно методом методом Рунге повысить порядок точности. Для этого образуем комбинацию ${\tilde I}^h = \sigma I^{h_1} + (1 - \sigma)I^{h_2}$ и выберем $\sigma$ так, чтобы $\sigma h_1^p+(1 - \sigma)h_2^p = (\sigma + (1 - \sigma)\rho^p)h^p = 0$ , то есть $\sigma = \frac{\rho^p}{\rho^p - 1} = \frac{1}{1 - A}$ . Тогда для погрешности ${\tilde D}^h = {\tilde I}^h - I$ получаем

${\tilde D}^h(f) = O(h^q)$ .

Числовой пример

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

$\int^{5}_{1} \ln x\, dx = x \ln x - x |_1^{5} = 5\ln 5 - 4$

В нижеследующей таблице представлены результаты работы программы:

r	Исходная формула	Экстраполированная формула	Точное значение	Погрешность вычислений	Погрешность формулы
2	3.98277278	4.04665506	4.04718956	0.0005345	0.00275556
4	4.03068449	4.04714980	4.04718956	0.00003976	0.00017222
8	4.04303347	4.04718692	4.04718956	0.00000264	0.00001076
16	4.04614856	4.04718939	4.04718956	0.00000017	0.00000067
32	4.04692918	4.04718955	4.04718956	0.00000001	0.00000004
64	4.04712446	4.04718956	4.04718956	0	0
20384	4.04718956

Здесь $r$ - коэффициент измельчения шага $h$ . Исходная величина шага $h = 2$ .

На иллюстрации черная сплошная линия - исходная формула, красная пунктирная - экстраполированная.

Как мы видим, разница между экстраполированными и неэкстраполированными результатами значительна. Уже при величине шага в $h = \frac{1}{64}$ мы можем найти значение интеграла с точностью $10^{-8}$ , тогда как в исходной формуле нам для достижения такой точности пришлось бы задать величину шага $h = \frac{1}{10192}$ .

Заключение

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

Список литературы

А.А.Самарский, А.В.Гулин. Численные методы М.: Наука, 1989.
Fundamental Methods of Numerical Extrapolation With Applications, mit.edu

Это незавершённая статья. Вы поможете проекту, исправив и дополнив её.

Источник — «http://machinelearning.ru/wiki/index.php?title=%D0%A3%D1%87%D0%B0%D1%81%D1%82%D0%BD%D0%B8%D0%BA:Gukov/%D0%9F%D0%B5%D1%81%D0%BE%D1%87%D0%BD%D0%B8%D1%86%D0%B0»

Категория: Незавершённые статьи

@@ Строка 45: / Строка 45: @@
 === Случай почти равных корней ===
-[[Изображение:macon.png|thumb|180px|Рисунок 1]]
+[[Изображение:macon.png|thumb|200px|Рисунок 1]]
 Условие 3 сходимости метода Ньютона-Рафсона означает, что никакие два корня не находятся слишком близко один к другому. Соответствующая ситуация представлена на рисунке 1 (масштаб сильно увеличен). Заметим, что производная <tex>f'(x)</tex> близка
 к 1 при x, равном обоим значениям корней, <tex>a_1</tex> и <tex>a_2</tex>. Более того, на основании теоремы Лагранжа, можно утверждать, что <tex>f'(x) = 1</tex> где-то между <tex>a_1</tex> и <tex>a_2</tex>.
@@ Строка 61: / Строка 61: @@
 === Поиск начального приближения ===
-[[Изображение:macon2.PNG|thumb|300px|Рисунок 2]]
+[[Изображение:macon2.PNG|thumb|200px|Рисунок 2]]
 Сначала находится значение x, при котором <tex>f'(x)=1</tex>, то есть решается уравнение
 :<tex>x=x+f'(x)-1</tex>