1611688890-f641c9ec8276824e4686da772eb56520 (826652), страница 54
Текст из файла (страница 54)
д. Столбец с номером j результирующей ортогональной3283. Численные методы линейной алгебрыТаблица 3.3. Ортогонализация Грама-Шмидтасистемы векторов {v1 , v2 , . . . , vn }DO FOR j = 1 TO nqj ← vj ;DO k = 1 TO j − 1αkj ← hqk , vj i ;qj ← qj − αkj qj ;END DOαjj ← kqj k2 ;IF ( αjj = 0 ) THENSTOP, сигнализируя «vj линейно зависитот векторов v1 , v2 , . . . , vj−1 »END IFqj ← qj /αjj ;END DOматрицы равен нормированной линейной комбинации первых j штукстолбцов исходной матрицы. В целом процесс ортогонализации ГрамаШмидта равносилен умножению W справа на верхнюю треугольнуюматрицу, в результате чего должна получиться ортогональная матрица.Фактически, ортогонализацию Грама-Шмидта можно рассматривать как ещё один способ получения QR-разложения матрицы, нарядус методом отражений (§3.7г) или методом вращений (§3.7д).20 Но устойчивость ортогонализации Грама-Шмидта существенно хуже.
Если исходная система векторов близка к линейно зависимой, то, выполняяалгоритм Грама-Шмидта с помощью приближённых вычислений (например, в арифметике с плавающей точкой современных ЭВМ), мы можем получить базис, который существенно отличается от ортогонального: попарные скалярные произведения его векторов будут заметно20 Верно и обратное: любое разложение матрицы на множители, в котором один изсомножителей ортогонален, соответствует некоторому процессу ортогонализации.Вопрос в том, насколько удобны и технологичны соответствующие алгоритмы.3.7.
Методы на основе ортогональных преобразований329отличны от нуля.Причина этого явления довольно прозрачна: при построении QRразложения с помощью матриц отражения или вращения ортогональность соответствующего матричного сомножителя специально контролировалась в процессе работы алгоритма, тогда как в ортогонализацииГрма-Шмидта мы идём обратным путём, от треугольной матрицы, иортогональность получается как побочный продукт вычислительногоалгоритма. Неудивительно, что эта ортогональность и не выполняетсядостаточно строго для результата выполнения реального алгоритма,подверженного влиянию погрешностей.Таблица 3.4. Модифицированный алгоритмортогонализации Грама-ШмидтаDO FOR j = 1 TO nqj ← vj ;DO k = 1 TO j − 1αkj ← hqk , qj i ;qj ← qj − αkj qj ;END DOαjj ← kqj k2 ;IF ( αjj = 0 ) THENSTOP, сигнализируя «vj линейно зависитот векторов v1 , v2 ,.
. . , vj−1 »END IFqj ← qj /αjj ;END DOНедостаточную устойчивость ортогонализации Грама-Шмидта можно до некоторой степени исправить, модифицировав её расчётные формулы так, чтобы вычисление поправочных коэффициентов αkj выполнялось другим способом. Псевдокод модифицированной ортогонализации Грама-Шмидта приведён в Табл. 3.4.В общем случае при ортогонализации Грама-Шмидта построениекаждого следующего вектора требует привлечения всех ранее постро-3303.
Численные методы линейной алгебрыенных векторов. Но если исходная система векторов имеет специальный вид, в определённом смысле согласованный с используемым скалярным произведением, то ситуация упрощается. Важнейший частныйслучай — ортогонализация так называемых подпространств Крылова.Определение 3.7.3 Пусть даны квадратная n × n-матрица A и nвектор r. Подпространствами Крылова Ki (A, r), i = 1, 2, .
. . , n, матрицы A относительно вектора r называются линейные оболочки векторов r, Ar, . . . , Ai−1 r, т. е. Ki (A, r) = lin {r, Ar, . . . , Ai−1 r}.Оказывается, что если A — симметричная положительно определённая матрица, то при ортогонализации подпространств Крылова построение каждого последующего вектора привлекает лишь два предшествующих вектора из строящегося базиса. Более точно, справедливаТеорема 3.7.2 Пусть A — симметричная положительно определённая матрица и векторы r, Ar, A2 r, . . . , An−1 r линейно независимы.Если векторы p0 , p1 , .
. . , pn−1 получены из них с помощью процессаортогонализации, то они выражаются рекуррентными соотношениямиp0 = r,p1 = Ap0 − α0 p0 ,pk+1 = Apk − αk pk − βk pk−1 ,k = 1, 2, . . . , n − 2,где коэффициенты ортогонализации αk и βk вычисляются следующимобразом:αk =βk =hApk , pk i,hpk , pk ik = 0, 1, . . . , n − 2,hpk , pk ihApk , pk−1 i=,hpk−1 , pk−1 ihpk−1 , pk−1 ik = 1, 2, .
. . , n − 2.Этот факт был открыт К. Ланцошем в 1952 году и имеет многочисленные применения в практике вычислений. В частности, он существенно используется в методе сопряжённых градиентов для решенияСЛАУ (см. §3.10г).3.7. Методы на основе ортогональных преобразований331Доказательство. Если векторы p0 , p1 , . . . , pn−1 получены из r, Ar,A2 r, . . . , An−1 r в результате ортогонализации, то из формул (3.89) следуетk−1X (k)(k)ci ∈ R.ci Ai r,pk = Ak r −i=0Как следствие, вектор pk − Apk−1 принадлежит подпространству, являющемуся линейной оболочкой векторов r, Ar, .
. . , Ak r, или, что тоже самое, линейной оболочкой векторов p0 , p1 , . . . , pk . По этой причинеpk+1 выражается через предшествующие векторы как(k)(k)pk+1 = Apk − γ0 p0 − . . . − γk pk(k)(k)с какими-то коэффициентами γ0 , . .
. , γk .Домножая скалярно полученное соотношение на векторы p0 , p1 , . . . ,pk и привлекая условие ортогональности вектора pk+1 всем p0 , p1 , . . . ,pk , получимhApk , pj i(k)γj =,j = 0, 1, . . . , k.hpj , pj iНо при j = 0, 1, . . . , k − 2 справедливо hApk , pj i = 0, так как hApk , pj i =hpk , Apj i, а вектор Apj есть линейная комбинация векторов p0 , p1 , . . . ,pj+1 , каждый из которых ортогонален к pk при j + 1 < k, т.
е. j ≤ k − 2.(k)Итак, из коэффициентов γj ненулевыми остаются лишь два коэффициента(k)αk = γk(k)=βk = γk−1 =hApk , pk i,hpk , pk ihApk , pk−1 i.hpk−1 , pk−1 iДалее,hApk , pk−1 i = hpk , Apk−1 i = hpk , pk + αk−1 pk−1 + βk−1 pk−2 i = hpk , pk i,и поэтомуβk =hpk , pk i.hpk−1 , pk−1 iЭто завершает доказательство теоремы.3323.83. Численные методы линейной алгебрыМетод прогонкиДо сих пор не делалось никаких дополнительных предположений оструктуре нулевых и ненулевых элементов в матрице системы. Но длябольшого числа систем линейных уравнений, встречающихся в практике математического моделирования ненулевые элементы заполняютматрицу не полностью, образуя в ней те или иные правильные структуры — ленты, блоки, их комбинации и т. п. Естественно попытатьсяиспользовать это обстоятельство при конструировании более эффективных численных методов для решения СЛАУ с такими матрицами.Метод прогонки, предложенный в 1952 году И.М.
Гельфандом иО.В. Локуциевским, предназначен для решения линейных систем уравнений с трёхдиагональными матрицами.21 Далее для краткости мынередко будем называть их просто “трёхдиагональными линейными системами”. Это важный в приложениях случай СЛАУ, возникающий, кпримеру, при решении многих краевых задач для дифференциальныхуравнений.
По определению трёхдиагональными называются матрицы,все ненулевые элементы которых сосредоточены на трёх диагоналях —главной и соседних с ней сверху и снизу. Иными словами, для трёхдиагональной матрицы A = ( aij ) неравенство aij 6= 0 имеет место лишьпри i = j и i = j ± 1.00Рис. 3.19. Портрет трёхдиагональной матрицы.В развёрнутом виде система n линейных алгебраических уравненийотносительно неизвестных x1 , x2 , . . . , xn , имеющая трёхдиагональную21 Ванглоязычной литературе этот алгоритм называют также “методом Томаса”.3333.8.
Метод прогонкиматрицу выглядит следующим образом:b1 x1 + c1 x2 = d1 ,ai xi−1 + bi xi + ci xi+1 = di , a x= dn .n n−1 + bn xn2 ≤ i ≤ n − 1,Часто её записывают в следующем специальном каноническом виде,даже без обращения к матрично-векторной форме:ai xi−1 + bi xi + ci xi+1 = di ,1 ≤ i ≤ n,(3.90)где для единообразия полагают a1 = cn = 0 в качестве коэффициентовпри фиктивных переменных x0 и xn+1 .
Подобный вид и обозначенияоправдываются тем, что соответствующие СЛАУ получаются действительно «локально», как дискретизация дифференциальных уравнений,связывающих значения искомых величин также локально, в окрестности какой-либо рассматриваемой точки.Пример 3.8.1 В §2.8 мы могли видеть, что на равномерной сеткес шагом hu(xi−1 ) − 2u(xi ) + u(xi+1 ),u′′ (xi ) ≈h2и правая часть этой формулы помимо самого узла xi , в котором берётся производная, вовлекает ещё только соседние узлы xi−1 и xi+1 .Поэтому решение конечно-разностными методами краевых задач дляразличных дифференциальных уравнений второго порядка приводитк линейным системам уравнений с трёхдиагональными матрицами, укоторых помимо главной диагонали заполнены только две соседние сней.Соотношения вида (3.90)ai xi−1 + bi xi + ci xi+1 = di ,i = 1, 2, .
. . ,называют также трёхточечными разностными уравнениями или разностными уравнениями второго порядка.Пусть для СЛАУ с трёхдиагональной матрицей выполняется прямой ход метода Гаусса без перестановок строк и столбцов матрицы, т. е.без специального выбора ведущего элемента. Если он успешно прорабатывает до конца, то приводит к системе с двухдиагональной матрицей334вида3. Численные методы линейной алгебры×××......00××,××(3.91)в которой ненулевые элементы (обозначенные крестиками) присутствуют лишь на главной диагонали и первой побочной сверху.
Следовательно, формулы обратного хода метода Гаусса вместо (3.60) должны иметьследующий двучленный видxi = ξi+1 xi+1 + ηi+1 ,i = n, n − 1, . . . , 1,(3.92)где, как и в исходных уравнениях, в n-ом соотношении присутствуетвспомогательная фиктивная неизвестная xn+1 . Оказывается, что величины ξi и ηi в соотношениях (3.92) можно несложным образом выразить через элементы исходной системы уравнений.Уменьшим в (3.92) все индексы на единицу —xi−1 = ξi xi + ηi— и подставим полученное соотношение в i-ое уравнение системы, чтодаётai ξi xi + ηi + bi xi + ci xi+1 = di .Отсюдаxi = −cidi − ai ηixi+1 +.ai ξi + biai ξi + biСравнивая эту формулу с двучленными расчётными формулами (3.92),можем заключить, чтоξi+1 = −ηi+1 =ci,ai ξi + bidi − ai ηiai ξi + bi(3.93)(3.94)для i = 1, 2, .
. . , n. Это формулы прямого хода прогонки, целью которого является вычисление величин ξi и ηi , называемых прогоночными3353.8. Метод прогонкикоэффициентами. Вместе с формулами обратного хода (3.92) они определяют метод прогонки для решения систем линейных алгебраическихуравнений с трёхдиагональной матрицей.Для начала расчётов по выведенным формулам требуется знать величины ξ1 и η1 в прямом ходе и xn+1 — в обратном. Формально онинеизвестны, но фактически полностью определяются условием a1 =cn = 0. Действительно, конкретные значения ξ1 и η1 не влияют на результаты решения, потому что в формулах (3.93)–(3.94) прямого ходапрогонки они встречаются с множителем a1 = 0. Кроме того, из формул прямого хода следует, чтоξn+1 = −cn0=−= 0,an ξn + bnan ξn + bnа это коэффициент при xn+1 в обратном ходе прогонки.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
