Г.М. Кобельков - Курс лекций по численным методам (1160467), страница 11

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

Текст из файла (страница 11)

Что это значит? Это значит, что мы будем искатьрешение x, для которого норма kAx − bk минимальна. То есть мы ищем перпендикуляр (то есть кратчайшее расстояние) к гиперплоскости, задающей решение системы. Сам по себе метод состоит в «покоординатном спуске»к этому перпендикуляру.А далее идёт «Метод Поспелова» в вольном изложении А. Воронцова.Все написанное ниже — описание того, как я понимаю, что такое метод Поспелова. Не исключено, что я сильно заблуждаюсьи речь в этом вопросе должна идти о чем-то совсем другом.

Те рассуждения, которые я додумал сам, я буду выделять вот такимшрифтом.Итак, мы хотим научиться решать уравнение Ax = b, не предполагая ничего о матрице A. Ясно, что еслиматрица вырождена, то система может не иметь решений или иметь не единственное решение.Та же проблема возникает, если мы решаем задачу с невырожденной матрицей, которая близка к вырожденной, например: A = C + εI, det C = 0. Такая матрица из-за вычислительных погрешностей может вести себякак вырожденная.Заменим исходную задачу на задачу минимизации функционала Φ = kAx − bk2 . Оказывается, что эта задачапоставлена корректно. Действительно, рассмотрим вариацию функционалаΦ(x + δ) = kAx − b + Aδk2 = kAx − bk2 + 2(AT (Ax − b), δ) + kAδk2 .40(131)Это значит, условие экстремума — равенство AT Ax = AT b.

Оказывается, что эта задача всегда имеет решение.Действительно, по теореме Фредгольма для того чтобы система имела решение необходимо и достаточно, чтобыAT b было ортогонально ядру AT A.Утверждение 4.5. Ядра операторов AT A и A совпадают. Ax = 0 ⇒ AT Ax = 0, AT Ax = 0, поэтому 0 = (AT Ax, x) = (Ax, Ax). Значит, Ax = 0. Утверждение 4.6. AT b ⊥ Ker AT A. Пусть x ∈ Ker AT A ⇒ Ax = 0 (по предыдущему утверждению).

Тогда (AT b, x) = (b, Ax) = 0. Замена уравнения Ax = b уравнением AT Ax = AT b называется первой трансформацией Гаусса.Осталось понять, как решать новую систему. Попытаемся минимизировать функционал Φ методом оптимального покоординатного спуска. Для этого сначала зафиксируем базис, в котором будем искать решение:выберем w1 , . . . , wq такие, что kAwi k 6= 0 и Awi линейно независимы.Собственно, насколько я понимаю, весь фокус в том, как выбирается этот базис.

Я думаю, что надо передтем как добавить очередной вектор wn проверить, что kAwn k > ε и получившаяся система достаточно линейно независима (то есть если представить Awn = en + r, где r ∈ hAw1 , . . . , Awn−1 i, а en ортогонально этомупространству ken k должна быть больше ε).В частности, если рассмотреть матрицу, в которой первые q собственных значений имеют порядок 1, аостальные λi ≪ 1, то мы выберем базис e1 , . . .

, eq .Далее организуем процесс вычислений следующим образом. Пусть у нас есть значение xn . Индекс у x вверхуозначает номер итерации. Для каждого k определим ck , которое минимизирует функционал Φ(xn + ck wk ). Этозначение определяется из уравненийΦ(xn + ck wk ) = kAxn − bk2 + c2k kAwk k2 + 2ck (Awk , Axn − b),(132)(Awk , Ax − b).(133)kAwk k2Из всех векторов wi выбираем тот, спуск вдоль которого дает наибольший выигрыш (в смысле уменьшениязначения функционала Φ).Остается показать, что описанный процесс сходится.

Для этого нам понадобитсяЛемма 4.7. Пусть g1 , . . . , gk — линейно независимые единичные векторы. Обозначим Lk = hg1 , . . . , gk i.Тогда существует γ такое, что для любого x ∈ Lkck = −kx − (x, gi )gi k 6 γkxk.(134)Здесь i выбирается так, чтобы значение |(x, gi )| было максимальным. Предположим противное. Пусть существует последовательность xk , таких что kxk k = 1 и1.(135)kНа единичной сфере можно выбрать подпоследовательность, которая сходится к некоторому элементу x∗ .

Дляэтого элемента имеемkx∗ − (x∗ , gi )gi k = kx∗ k2 − 2(x∗ , gi )2 + (x∗ , gi )2 > 1.(136)kxk − (xk , gi )gi k > 1 −Поскольку kx∗ k = 1 получаем, что (x∗ , gi ) = 0 для всех i. А поскольку gi образуют базис в Lk , это означает,что x∗ = 0. Эта лемма означает, что если рассмотреть проекцию невязки на подпространство, натянутое наAw1 , . . . , Awq , ее норма будет убывать со скоростью геометрической прогрессии.

Что и требовалось. От метода Тихонова метод Поспелова, насколько я понимаю, отличается тем, что здесь мы просто забилина тот кусок матрицы, который плохо обусловлен, и считаем, что решение имеет вид (x1 , . . . , xq , 0, . . . , 0).5. Нелинейные и дифференциальные уравнения5.1. Нелинейные уравненияПусть нам задано нелинейное уравнение f (x) = 0.5.1.1. Метод половинного деленияЕсли известно, что функция хорошая, и известны две точки x0 и x1 , для которых f (x0 )f (x1 ) < 0, то можнозапускать метод половинного деления: делим отрезок пополам, смотрим, на котором из двух новых отрезковзначения разных знаков, выбираем его и так далее.Он, разумеется, плох тем, что совсем не обобщается на многомерный случай.

Кроме того, на функциюнакладываются существенные ограничения.415.1.2. Метод простой итерацииМетод простой итерации заключается в том, что мы переходим от уравнения f (x) = 0 к эквивалентномууравнению x = g(x), а затем рассматриваем итерационный процесс xn+1 = g(xn ). Если отображение g сжимающее, то он сойдётся к некоторому решению этого уравнения.В одномерном случае достаточным условием для сжимающего отображения будет условие 0 < |g ′ (x)| < 1 втом интервале, куда попадают итерации отображения. В самом деле, пусть x — точное решение. Тогдаxn+1 − x = g x + (xn − x) − g(x).(1)Разложим правую часть в ряд Тейлора:xn+1 − x = a1 (xn − x) + a2 (xn − x)2 + .

. .(2)Если g ′ (x) 6= 0, то a1 6= 0. А если |a1 | < 1, то, понятное дело, метод будет сходиться со скоростью геометрическойпрогрессии.5.1.3. Метод НьютонаРешаем уравнение f (x) = 0. Будем считать, что наша функция достаточно хорошая (что под этим понимается, уточним чуть позже). Напишем формулу Тейлора:f (x) = f (x) + f ′ (ξ)(x − x),(3)ξ = ξ(x).Имеем f (x) = 0, поэтому возникает идея написать такой итерационный метод:f (xn ) + f ′ (xn )(xn+1 − xn ) = 0,то естьxn+1 = xn −(4)f (xn ).f ′ (xn )(5)Это и есть метод Ньютона.Будем считать, что функция f строго дифференцируема, то есть2kf (x) − f (y) − f ′ (y)(x − y)k 6 M1 kx − yk .(6)′−1Чтобы при решении уравнения на xn+1 не возникало проблем, естественно требовать существования (f (x)) .Более того, нам потребуется, чтобы в некоторой окрестности она была ограничена некоторой константой M2 .Пусть kxn − xk 6 b.

Покажем, что при этих условиях отображение будет сжимающим.В самом деле, имеем2kf (x) − f (xn ) − f ′ (xn )(x − xn )k 6 M1 kxn − xk .(7)Снова замечая, что f (x) = 0 и подставляя в эту формулу значение f (xn ), получими после сокращений получаем′kf ′ (xn )(xn+1 − xn − x + xn )k 6 M1 kxn − xk2 ,(8)kf ′ (xn )(xn+1 − x)k 6 M1 kxn − xk2 ,(9)′−1Обозначим z := f (xn )(xn+1 − x), тогда xn+1 − x = (f (xn ))z, откуда2kxn+1 − xk 6 M2 kzk 6 M1 M2 kxn − xk 6 M1 M2 b2 .(10)Значит, если мы хотим, чтобы всё это сходилось, нам нужно выбрать такую окрестность, чтобы выражениесправа было меньше b, значит, b < M11M2 .Имеем2δn+1 := M1 M2 kxn+1 − xk 6 M12 M22 kxn − xk .(11)Значит, имеется неравенство δn+1 6 δn2 . Значит,nδn 6 δ02 ,(12)то есть скорость сходимости просто бешеная.Замечание.

Метод Ньютона в одномерном случае имеет довольно наглядный геометрический смысл, поэтому его часто называют методом касательных.√Пример 1.1. Вычислим a с помощью метода Ньютона. Рассмотрим уравнение x2 − a = 0. Итерационныйпроцесс имеет видx2 − ax2 + axn+1 = xn − n= n.(13)2xn2xnВся проблема метода Ньютона — это попасть в ту окрестность, в которой он начинает сходиться.425.2. Дифференциальные уравненияРассмотрим простейшую задачу Коши:(y ′ = f (x, y),y(0) = y0 .(14)5.2.1.

Метод Эйлера и его модификацииПусть нам известно значение y(x). Будем строить значение y(x + h). Конечно, можно разложить функциюв ряд Тейлора. Но это плохо. Гораздо лучше (с вычислительной точки зрения) поступать так: проинтегрируемнаше уравнение от x до x + h. По формуле Ньютона – Лейбница получимy(x + h) = y(x) +Zh0f x + t, y(x + t) dt.(15)Конечно, мы ничего не знаем про то, чему равна функция y на отрезке [x, x + h]. Поэтому мы заменим интегралсамой простой квадратурной формулой. В итоге получим такое соотношение:yn+1 = yn + hf (xn , yn ).(16)При этом значение интеграла будет вычислено с точностью до O(h2 ).

Характеристики

Список файлов лекций