Численные методы. Ионкин (2012) (неоффициальные) (косяки есть) (1160437), страница 7
Текст из файла (страница 7)
· H2T H1T = QTИтак, Q−1 = QT ортогональная матрица, а A представима в виде:A = QRПолучили, что любую матрицу можно разложить в виде QR, где Q — ортогональная матрица, а R — матрица верхнетреугольной формы.Если A — полная матрица, то разложение (1) пропорционально m3 действий.Если матрица A имеет верхнюю почти треугольную форму, то разложение требуетпорядка m2 действий. А если A — трёхдиагональная, то необходимо всего порядкаm действий.На первом этапе в QR алгоритме обозначимA = A0 .Вначале делается QR разложение этой матрицы:A0 = Q0 R0 ,(2)A1 = R0 Q0 .(3)Ts где Q−10 = Q0 , a R - ВТФ.Далее находим A1 :Утверждается, что она сохраняет спектр A0 , а значит и спектр A.Тогда из (2) получаем, чтоR0 = Q−10 A0 .47Затем, подставляя в (3), получимA1 = Q−10 AQ0 .Откуда видно, что матрицы A1 и A0 подобны, следовательноA01λAk = λk ,k = 1, m.В результате этих двух шагов получаем матрицу A1 , которая сохраняет спектри симметрию(если она была) матрицы A0 , так как A1 получается из неё c помощьюпреобразования, использующего ортогональную матрицу.
Следующим шагом разложим A1 в видеA1 = Q1 R1 .Тогда на втором шаге получим матрицы A2 :A2 = R1 Q1 .Ясно, что мы записали частный случай рекуррентного итерационного метода:Ak = Qk Rk ,(4)Tгде Q−1k = Qk , Rk — ВТФ, аAk+1 = Rk Qk ,kλnAk+1 = λAn ,k = 0, 1, 2, . . .(5)n = 1, mКроме тогоAk+1 = Q−1k Ak Qk .Мы построили итерационный процесс (4), (5).В подавляющем большинстве случаев, если λk — вещественные, Ak сходится кматрице верхней треугольной формы:× × ... × × 0 × . . . × ×lim Ak = 0 0 .
. . × × .k→∞ .. ... . . . . . . . . .. 0 0 ... 0 ×Заметим, что под главной диагональю могут и не получаться нули в математическом смысле. Достаточно, чтобы значения под главной диагональю были по модулюменьше некоторого числа (так называемый машинный ноль), определяющего точность вычислений. А на главной диагонали будут находиться собственный значенияматрицы.48У вещественной матрицы могут быть и комплексные собственные значения.
Тогда, в случае с комплексными с/з матрица примет вид:× ................... 0 × . . . . . . . . . . . . . . . .0 0 xy . . . . . .lim Ak = 0 0 −y x . . . . . . .k→∞ .. ...... ... .. .. ...0 0 0 ... 0 ×Теперь на диагонали появились клетки второго порядка такие, что на главнойдиагонали стоят вещественные части собственных значений, а под ней — мнимые.Таким образом QR алгоритм (4)-(5) сходится в случае вещественных собственных значений к предельной верхнетреугольной матрице, а в случае комплексныхсобственных значений появятся на главной диагонали клетки второго порядка.Подчеркнем тот факт, что QR алгоритм применим к любой матрице и не ухудшает верхнюю почти треугольную форму матрицы.
Таким образом, для уменьшениячисла итераций, матрицу сначала можно свести к верхней почти треугольной форме.К минусам QR алгоритма можно отнести следующее:1. надо помнить промежуточные матрицы2. большое число действий, необходимых для получения предельных матриц.§12 Предварительное преобразование матрицы кверхней почти треугольной формеAk = Qk RkAk+1 = Rk Qk ,k = 0, 1, 2, . . .(1)(2)Для того, чтобы доказать, что QR алгоритм не будет портить верхнюю почтитреугольную форму, докажем 2 леммы.Лемма 1.
Пусть матрицы A, B, C — матрицы одного порядка, такие что B имеет верхнюю треугольную форму, A – верхнюю почти треугольную форму, тогдаC = BAимеет верхнюю почти треугольную форму.Доказательство: Элементы матрицы C имеют видcij =mXα=1biα aαj(3)49Теперь будем учитывать, что B имеет ВТФ, тогда biα = 0 ∀i > α. Откудаполучим, чтоmXcij =biα aαj .α=iТеперь учтем, что A — ВПТФ ⇒ aαj = 0 ∀α > j + 1, откудаcij =j+1Xbiα aαj .α=iРассмотрим cij : i > j + 1. Такие элементы, очевидно, будут равны нулю, а это иозначает, что матрица C - ВПТФ.Лемма 2.
Пусть матрицы A, B, C — матрицы одного порядка, такие что B имеет верхнюю треугольную форму, A – верхнюю почти треугольную форму, тогдаC = AB(3)имеет верхнюю почти треугольную форму.Доказательство: Элементы матрицы C имеют видcij =mXaiα bαjα=1Теперь будем учитывать, что B имеет ВТФ, тогда bαj = 0 ∀α > j. Откудаполучим, чтоjXcij =aiα bαj .α=1Теперь учтем, что A — ВПТФ ⇒ aiα = 0 ∀i > α + 1, откудаcij =jXaiα bαj .α=i−1Рассмотрим cij : i > j + 1. Такие элементы, очевидно, будут равны нулю, а это иозначает, что матрица C - ВПТФ.Вернемся к QR алгоритму.
Rk имеет верхнюю треугольную форму. Из (1)Qk = Ak Rk−1 . Изначально мы A0 привели к верхней почти треугольной форме ⇒Ak имеет верхнюю почти треугольную форму, тогда мы получаем, что матрица Qkимеет верхнюю почти треугольную форму. Из (2) получаем, что Ak+1 тоже имеетверхнюю почти треугольную форму.Глава 2Интерполирование и приближениефункций§1 Постановка задачи интерполированияПусть есть функция f (x) и известны ее значения в узлах x ∈ [a, b],a 6 x0 < x1 < .
. . < xn 6 b{xi }n0 — узлы интерполирования. В этих узлах заданы значения функцииf (xi ) = fi ,i = 0, n(1)Надо построить такую непрерывную функцию, которая в узлах совпадает с fi .Эта задача имеет бесконечное число решений. Мы рассмотрим случай, когда функция f (x) будет приближаться алгебраически полиномами.
Пусть есть полином n-ойстепениPn (x) = a0 + a1 x + · · · + an xn , ai – вещественные числа, an 6= 0(2)Нужно построить такой полином, который в узлах будет совпадать с интерполированной функцией. Он будет называться интерполяционным полиномом, если:Pn (xi ) = fi ,i = 0, n(3)Эта задача имеет всегда единственное решение. Для нахождения коэффициентов aiиз равенства (2) получаем из (3) систему:a0 + a1 x0 + a2 x20 + · · · + an xn0 = f0a + a x + a x2 + · · · + a xn = f01 12 1n 11·································a0 + a1 xn + a2 x21 + · · · + an xnn = fn50(4)51Выпишем определитель матрицы из коэффициентом этой системы(определительВандермонда): 1 x0 x2 · · · xn 00 1 x2 x2 · · · xn Y22(xi − xj ) .. ......
=. ..···. n>i>j>01 xn x2n · · · xnn Если узлы разные, то определитель отличен от нуля и, согласно критерию, система (4) имеет единственное решение.Замечание. Поскольку показано существование и единственность интерполирующего полинома, то при его поиске в любой форме он будет тождественно равенвсем своим представлениям в иных формах, полученных с помощью других методов.§2 Интерполяционная форма Лагранжа Ln(x)Пусть существуют (n+1) несовпадающий узел, функция интерполирования f (x),такая чтоf (xi ) = fi , i = 1, n.Мы хотим построить полином Ln (x) так, чтоLn (xi ) = fi .(1)Тогда говорим, что мы построили интерполяционный полином для функции f (x)по узлам xi .Если мы его ищем в видеLn (x) =nXcn (x)f (xk ),(2)k=0где cn (x) — полином n-ой степени, тогда говорят, что строят полином Лагранжа.Введем полином (n + 1)-ой степени:ω(x) = (x − x0 )(x − x1 ) · · · (x − xn ).Обозначим через [·] все члены произведения в ω, за исключением (x − xk ).
Тогдаω 0 (x) = [·] + (x − xk )[·]0 .(3)Следовательно,ω 0 (xk ) =nY(xk − xi ).k=0,k6=iТогда из (1)ω(x),ck (x) =(x − xk )ω 0 (xk )Ln (x) =nXk=0ω(x)f (xk ).(x − xk )ω 0 (xk )(4)52Легко можно получить оценку погрешности Ln (x):ψLn (x) = f (x) − Ln (x).Введем функцию:g(s) = f (s) − Ln (s) − kω(s),k = const.Выберем постоянную k из условия g(x) = 0 в данной точке x 6= xi . Константаf (x) − Ln (x)k(x) получается равной k =. Затем, полагая f (x) ∈ C n+1 [a, b] и примеω(x)няя (n+1) раз теорему Ролля, получаем, что ∃ξ ∈ [a, b], в которой g n+1 (ξ) = 0.
Тогдаочевидно что:f (n+1) (ξ)ψLn (x) =ω(x)(n + 1)!Возьмем это выражение по модулю, тогда|ψLn (x) | 6Mn+1|ω(x)|,(n + 1)!где Mn+1 = sup |f n+1 (x)|xЗамечание. Если интерполируемая функция f (x) полином степени не выше n, тополином Лагранжа точно её приближает, так как Mn+1 = 0.Замечание. Полином Лагранжа, вообще говоря, не сходится к f (x).§3 Разделенные разностиПусть у нас есть отрезок [a, b], на котором выбраны узлы {xi }n0 , в которых заданызначения некоторой функции f (xi ) = fi , i = 0, n:a 6 x0 < x1 < · · · < xn 6 bОпределение. Разделенной разностью 1-го порядка, построенной по узлам xi , xjназывается отношениеf (xj ) − f (xi )f (xi , xj ) =.xj − xiУзлы xi и xj могут быть и не соседними.
Далее для простоты будем рассматривать разделенные разности, построенные по соседним узлам. По разделенной разности первого порядка, можно построить разделенную разность 2-го порядка, например:f (x0 , x1 , x2 ) =f (x1 , x2 ) − f (x0 , x1 )x2 − x0f (x1 , x2 , x3 ) =f (x2 , x3 ) − f (x1 , x2 )x3 − x153Можно ввести разделенную разность любого порядка.
Пусть задана разделеннаяразность k-го порядка f1 (xj , xj+1 , . . . , xj+k ) по узлам (xj , xj+1 , . . . , xj+k ) и разделеннаяразность k-го порядка f2 (xj+1 , xj+2 , . . . , xj+k+1 ) по узлам (xj+1 , xj+2 , . . . , xj+k+1 ), тогдаразделенная разность k + 1 порядка равна:f (xj , xj+1 , xj+2 , . . . , xj+k+1 ) =f2 − f1xj+k+1 − xjУтверждение. Разделенная разность k − го порядка представляется в видеkXf (xi ),f (x0 , x1 , .
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
