Самарский А.А., Гулин А.В. Численные методы (1989) (1095856), страница 74
Текст из файла (страница 74)
Прямой ход заключается в нахождении матриц Спч и векторов Р) ' по формулам (8), (9). Обратный ход состоит в нахождении векторов у, из системы (7), начиная с Й=гп. Метод редукции в том виде, как он здесь изложен, не применяется в реальных вычислениях по двум причинам. Во-первых, он не- экономичен из-за того, что на каждом этапе приходится обращать матрицу С("' общей структуры. Во-вторых, вычисление правых частей по формулам (9) неустойчиво.
В следующих пунктах будет показано, как можно устранить указанные недостатки метода редукции. 2. Обращение матриц. Соотношение (8) представляет собой нелинейное разностное уравнение первого порядка для матриц. Его решение можно найти в явном виде. Рассмотрим сначала числовой аналог уравнения (8), а именно разностное уравнение уэ у11 2 А 1 2 уо=х (10) где х — заданное число. Покажем, что решение у„=у„(х) уравнения (10) выражается через многочлен Чебышева, а именно уз (х) = 2Т 4 [ — 1, й = 1, 2,..., (2) (11) где ~ь уь=Ц (х — 2 соз (21 — 1) я1 Приведенные выше выводы имеют место и для матричного уравнения (8).
По индукции легко доказывается, что решением Спч уравнения (8) является многочлен (относительно матрицы С) степени Р со старшим коэффициентом 1. Для этого многочлена справедливо представление, аналогичное (11), т. е. Саа =2Т,э ( — ), и справедливо разложение на линейные множители ~Ф С'м=П (С вЂ” 2соз( ) "Е) . 2ьм !=1 (12) 421 соз(пагссозх), если ~х~ <1, Т„(х) = — [(х + 1/ х' — 1)" + (х — )/х* — 1) "1, если ~ х ~,;э 1.
2 Из выражения для Т„(х) следует, что Т,„(х) =2(Т„(х) )' — 1. Отсюда при л= 2"-' получаем Т,ь( — ) = 2 (Тн-, (~)) — 1, т. е. функция (11) удовлетворяет уравнению (10). Корнями многочлена (11) являются числа х~л — — 2 сов,, 1=1,2, „, 2". (21 — 1) я Поэтому функция у„(х), представляющая собой многочлен относи- тельно х степени 2" со старшим коэффициентом 1, разлагается в произведение Наличие такого разложения позволяет упростить процедуру обращения матриц. В случае разностных аппроксимаций уравнений эллиптического типа матрица С является трехдиагональной, в то время как Сии — матрицы общей структуры, Благодаря разложению (12) обращение матрицы Соч сводится к последовательному обращению трехдиагональных матриц Сна=С вЂ” 2соз „„Е, 1=1, 2,..., 2". (13р Действительно, пусть требуется решить уравнение Скчо=ф, где аа С' ' = Г( Сааь а=а (14) Рассмотрим для наглядности сначала случай, когда й=2.
Тогда (14) принимает вид Са аС, аСааСа ао = ф (15) где оа=ф, о,=о. Точно так же решение системы (14) в общем случае сводится к последовательному решению систем уравнений С,,о,=о, „1=1, 2, ..., 2', (16) где о,=ф, оаа = о. Если матрицы Саа трехдиагональные, то каждую из систем (16) можно решить методом прогонки. Указанный выше способ обращения матрицы Соч предполагает определенный порядок выполнения промежуточных этапов: сначала надо обратить матрицу С„„затем — матрицу С„, и т. д. Однако, поскольку в матрице Сол все сомножители перестановочны, можно использовать и другие способы введения промежуточных значений о„1=1, 2,..., 2' — 1.
Так, систему (15) можно записать в виде Са;Са,аСа аСадо=ар и заменить системой уравнений Сааоа=па, Сааоа=па, Сааоа=оа, Сааоа=он где о,=ф, о,=о. Теоретически при любом порядке выполнения промежуточных этапов мы должны получить одно и то же решение о. Однако, если число промежуточных этапов велико, то при решении систем уравнений (16) на ЭВМ будет происходить накопление погрешностей округления. Рост погрешностей округления зависит от порядка выполнения промежуточных этапов.
Здесь имеет место 422 Обозначим оа = ф, па = Са,аСа,аСа,ао, па= Са,аСазо, оа= Са ао Тогда получим, что решение системы (15) сводится к последовательному решению четырех систем уравнений Садоа = па, Салоа= он Са аоа = па Са аоа= оаа Отсюда, учитывая, что С'"=(С" о)* — 2Е, придем к уравнению (С'й ")' (р';" — р';~ч) + <);" = 2р';ы + <!<й <й, + у,'~~'1 + 1 й С(й-1) ( (й-1) ! (й-1) ! (й-1)) Выберем теперь вектор д< ' таким, чтобы выполнялось соотнои) шение (18) Тогда предыдущее уравнение после умножения на (С"-") ' примет вид Сии(и= ~''+ (о т "о (19) = Р;,й-< + Р<„й-1 т <)< где з <1-1) <и <1-1) < =у< Таким 'образом, в,(числение векторов Р(~) можно заменить нахождением векторов р«, )< из систем,( уравнений (17) — (19).
Сначала, (й) 9) обращая матрицу С" ", находим вектор з<(' ", затем в,(числяем ри)= (й-1) (й-1) (й) = р< + з< и затем по Формуле (!8) в )числяем <)< . Правую часть (й-1) г< уравнения (7) надо заменить согласно (!7) выражением Е(й-1) (й-1) (й-1) (й-1) '< = '- Р< + <)< Тогда придем к уравнению (й-1) <й-1) <1-1) С !< =д, +у<,+ у„„„ (20) где !< =-у< — р< . Из этого уравнения, обращая матрицу С' ", (й-и (й-0 (й-1) находим !и " и затем вычисляем у< = р<й "+ !(<' "'.
4ЗЗ примерно та же ситуация, что и в итерационном методе с чебышевским набором параметров (см. п. 2 $ б гл. 2 ч. П). Поэтому при реальном решении системы (!4) рекомендуется обращать внимание на порядок выполнения промежуточных этапов. Более подробно этот вопрос рассмотрен в книге [35). 3. Вычисление правых частей.
В методе редукции правые части г(1 ') уравнения (7) должны удовлетворять соотношению (9). Будем искать Р(<"' в виде (17) где векторы р(« ', 7<1 ' подлежат определению. Подставляя (17) в (9), получим С )РФ +у! =С Р +<). + Р + <)< = Р< йй-1+ Ц йй-1+ +С' "(С"" "и+ '-' и)+ .' " '," + !,о . 4.
Формулировка н обсуждение алгоритма. Сформулируем более детально алгоритм метода редукции. Вычисления ведутся в цикле по индексу А. Сначала осуществляется прямой ход, т. е. для 1=1, 2, ..., и — 1 решаются уравнения Са-1) а-1) (л-11 (й-13, (1-1) (21) з1 = Р; 1-1+ Р„11- + т1 и находятся векторы (М Й-11, (й-0 Р1 =Р1 „+з1 М~ М1! а-1~ а-1~ 41 =2Р1 + Ч; 11-1 + Ч; 11-1 где 1=2", 2 2', 3.2", ..., 2" — 2'.
Вычисления ведутся, начиная с 1=1, причем при Й=! задаются начальные значения р~м=О, 4~1'=Рь 1=1,2,..., 1Ч вЂ” 1. Решение систем (21) сводится, как это было показано в п. 2, к решению более простых систем уравнений (22) 1=2", 2 2',..., 2т — 2", где (А-а (а-1> й-11 оо,1 = Р1 11-1 + Р1ч11-1 + 41 14-0 О Ь.1;=З1 Системы (22) решаются при фиксированном ! для каждого 1=2', 2 2", ..., 2" — 2'.
Поскольку матрица С„,, не зависит от 1, вычисления целесообразно организовать так, чтобы матрица С„,, обращалась один раз. Подсчитаем число действий, необходимых для решения всех систем вида (22). Предположим, что для решения системы (22) с фиксированными й, 1, ! гоебуется 1! арифметических действий. Тогда при фиксированных й, ! для решения систем (22) с 1=2', 2 2*... 2т ь ..., 2 — 2' требуется д=(2" "— 1)д арифметических дей- 2" ствий. При фиксированном А и для 1=1, 2, ..., 2" ', 1=2", 2 2", ... ..., 2" — 2', требуется 2" '(2"-' — !)д= (2"-' — 21-')д действий. Наконец, для всех 1=1, 2,..., и — 1 необходимо выполнить т 1 д~ (2 ' — 2' ') = (1+ (и — 2) 2" ') д 4=1 арифметических действий.
Самое существенное здесь то, что при больших У число действий является величиной 0(дУ !оа, М). После того как все векторы р1, д1 найдены, осуществляется 1М КЧ обратный ход метода редукции, т. е. начиная с й=и решаются 424 уравнения (й-и %-и (к-а С 1,' =Ф +у,,+у„... (23) н вычисляются искомые векторы (Й-1), (м-1! У~= Р' + где й = т, т — 1,..., 1, 1= 2"-', 3 2'-',..., У вЂ” 2'-', У = 2". Система (23) при фиксированных 1, й решается, как и ранее, путем последовательного решения систем уравнений С...в,з=пь ьь 1=1,2,...,2" ', (24) 1=2" ' 3 2" ', 5 2" ',..., 2" — 2"-', м-и где шох=Ф + У;,к-х+ У;,ь м икь-1,.= 1~ .
Решение системы (23) также требует 0(УУ!оп',У) действий. Рассмотрим применение метода редукции к решению разностного уравнения Пуассона (4), (5) из $5. Предположим, что сетка содержит У,=2" точек по направлению х, и У, точек — по х,.
Согласно (9) из э 5, это разностное уравнение можно записать в векторном виде (1), (2), где У=У, и С=2Е,— й,'Л вЂ” трехдиагональная матрица порядка У,— 1. Поэтому системы вида (22), (24) решаются методом прогонки, что требует У=О(У,) действий. Таким образом, решение разностного уравнения Пуассона осуществляется методом редукции за число действий 0(У,У,1од,У,). На квадратной сетке, когда У,=У,=У, число действий является величиной 0(У*1од,У), т. е. имеет тот же порядок, что и в методе быстрого дискретного преобразования Фурье.
В отличие от последнего, метод редукции не требует знания собственных функций, что позволяет применять его и в более общих случаях, например в случае краевых условий третьего рода. Метод редукции выгодно отличается от метода матричной прогонки не только числом действий, но и требуемой памятью ЭВМ. В то же время следует еще раз подчеркнуть, что метод редукции можно применять только для решения относительно простых систем уравнений, а именно систем, которые можно записать в виде (1) с постоянной матрицей С.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
