Прямые методы для разреженных матриц. О. Эстербю, З.Златев (984134), страница 8
Текст из файла (страница 8)
Мы не можем утверждать вообще, что 1СуМЯ лучше бМ,"у, т. е. приводит к более устойчивым вычислениям„ хотя имеется немало экспериментальных свидетельств этому. Однако для специальных классов матриц превосходство 1ЙМЯ можно докязять. Теорема 3.11. Если матрииа А имеет диагональное преоблада8ие т и СТРУктУРЙУю симлЯГРию тт ТО гаУссово исключе" иие устойчиво для любой стратегии 1ОМЯ а), Замечание 3.12.
Хорошо известно (см., например, 18О1), что для этого кла а матр ц не нужен Выбор Главного эле мента по устойчивости; однако нам могут понадобиться перестановки для сохранения разреженности. В Доказательство. Пусть 1 ~ 1 е= и — 1. Предположим, что на первых 1 — 1 шагах гауссова исключения на роль главных Выбирались тОлько ДиаГОнальные элементы. При таком Выборе симметрия структуры (ап Ф О.е=: ап Ф О) сохраняется и активная часть матрицы А на 1с-м шаге, т.
е. подматрица Аи, по.прежнему имеет диагональное преобладание (~ аЖ~ > > К ~ф!). Поэтому 1Ф! С(1, 1)=г(1, 1с), 1=1с(1)п, Пусть г (11, 1с) = г (4, Е) =... = г (1„1), 1 ~ з .'. р, Тогда диагональные элементы строк с номерами ~ь 1а, ..., 1, принадлежат Ся, и наибольшим элементом Са ЯВляется Один из них. Он и будет выбран главным элементом 1с-го шага О Т. е. ~ аи ~ >,'»' ~ аи 1, 1= 1 (1) и. — Прим. перев.
1Ф1 У~ Т. е. ап ~ь О:=р- ан Ф О. — Прим. перев. э> Слово ~любая» употреблено авторами потому, что 1ШБ зависит от параметра и и даже при фиксированном и может не определять выбор однозначно (если в Сь несколько элементов с максимальным модулем),— Прим. Иври. любой стратегией 16МЬ независимо от коэффициента устойчивости и. Поскольку это справедливо для произвольного 1~ (1 » » 1» и — 1), то рассуждение по индукции показывает„что главными элементами исключения будут только диагональные элементы. Анализ Уилкинсона 180, 81, 107~ дает оценку Ь„» 2Ь|, свидетельствуя о численной устойчивости. В Пример 3.13. Что выбор главного элемента для сохранения разреженности может быть необходим, демонстрирует матрица с приводимым ниже портретом.
й о ~ о ... о х о о о Если не делать перестановок, то заполнение будет абсолютным, т. е. разреженность будет полностью утрачена. В то же время если матрица имеет диагональное преобладание, то для любой стратегии 1бМЯ никакого заполнения не будет. И Замечание 3.14. Требование диагонального преобладания можно несколько ослабить 184, стр. 37~, Если внимательней просмотреть доказательство теоремы 3,11 и анализ Уилкинсона 180~, обнаруживается, что достаточно, чтобы ~ а,! ~ ,'~'.
~ а,, ~„1 = 1 (1) п, )Ф~ со строгим неравенством хотя бы для одного значения 1, и ~ аи ~ > гпах ~ а.„,. ~, 1 = 1 (1) и. Тестовые матрицы класса Е(п,с) имеют диагональное преобладание в этом более слабом смысле. Ю Теорема 3.15. Если матриу,а А симметрична и положительно определена„то гауссово исключение устойчиво для любой стратегии 16МЗ при и = оо. Доказательство.
Рассуждение следует той же схеме, что и доказательство теоремы 3.11. Поскольку и = оо, все элементы подматрицы порядка з, стоящей на пересечении строк и столбцов Аоо с номерами 1ь 4, ..., 1„принадлежат в то же время Сь Так как эта подматрнца положительно определена, ее наибольший элемент находится на главной диагонали и бу- дет выоран главным элементом 1~-го шага любой стратегией 1бМЯ. й Пример 3.16. Матрица положительно определена.
Если и =. 5, то на первом шаге Гауссова исключениЯ Главным будет назван элемент а13 —— = 900; при этом перестановки нарушат структурную симметрию, Если п 5 и р.==2„то будет выбран элемент а,2= 10 (а при и > 900 и р = 1 — элемент а 1 — — 1), и симметрия сохранится, Мы можем брать сравнительно малые главные элементы, поскольку, как показано Уилкинсоном 180~, для положительно Опред~ленной матрицы это не привОдит к чис- Л~Н~ОЙ НеуСтойчивости.
И Все последовательности Главных элементов являющиеся результатом применения 16МЯ, возможны и в ИМЯ (при тех же р и и), Но если ЙМЯ хотя бы раз выбирает главный элемент вне главной диагонали, то структурная симметрия будет утеряна, и нельзя гарантировать устойчивость.
Как показывает пример 3.8, это действительно случается на практике. Отметим еще, что проблем, с которыми столкнулась бМЗ в примере 3.7, для 16М8 не существует. Можно идентифицировать свойства матрицы, которые делают для нее 1бМЯ более привлекательной стратегией, чем ОМАМ. Если на нескольких шагах исключения имеются кандидаты в Главные элементы, малые приблизительно в такоЙ степени, какую допускает коэффициент устойчивости п, то вероятность выбрать один из них велика, и для ОМАМ часто будет Ь~+1 ~ (1 + п) Ь~ (см. 1861 ) . Пример 3,8 показывает, что последствия этого могут быть катастрофическими. Чтобы провести более полное исследование, мы проделали вычисления для 96 матриц класса Е(п, с); и = 250(50) 1000, с = 4(40)204.
Программа, выбранная, чтобы представлять ИМЯ, — это программа Е01ВКЕ библиотеки ХАб (см. 117„241) с рекомендованным значением п = 10 для исключения; для подстановок использовалась программа ЕО4АХЕ той же библиотеки. Представителем 16МЯ была программа У12М 197, 1021 также при и = 10; правда, для этих матриц при пользовании 16МЯ значение ц не существенно.
Вычисления выполнялись на машине 11КЮАС 1100~82 в КЕСК1! (Вычислительный центр при Университете Копенгагена) с одинарной точностью (е ~ 1.5Š— 8). Ни в одном случае программы, основанные на ИМЯ, не были лучше; при этом приблизительно для 25% матриц разложения, Вы сленные программой РО1ВКЕ, были очень неточны и то же относится к соответствующим решениям, вычисленным посредством РО4АХЕ. Правые части выбирались так, чтобы точным решением всегда был вектор, состоящий из единиц.
Таблана ЗА. Сравнение программ библиотеки ХА6 с программой т"12М для матриц класса Е (800, с). В этом тесте т12М использовалась с Х = 0.01 и итерационным уточнением РО! ВРЕ + РО1АХВ Соинт Погрешноеть Потрешысть Время Время В табл. 3.4 приведены время счета, значение СО13ХТ и погрешность для типичного значения и (и = 8ОО). Видно, что У12М лучше во всех случаях и в любом отношении. В частности, ее преВосходстВО Велико для Фпромежуточных значе ний» с (здесь 44 и 84). Заметим, что типичным для практики значением было бы с= ~п ж 28. Для табл. 3.5 мы взяли с = 44 и сравниваем точность для различных и. В последний столбец помещены результаты, полученные с использованием итерационного уточнения (см. гл. 4).
Таблица 8.б. Сравнение погрешности для программ библиотеки ХА6 и программы т"12М (в режимах 0$ н 1К) для матриц класса Е (и, 44) '1'12М вЂ” 1Н т = О.о1 а Р01ВЙВ + РО4ЛХВ 9420 30424 22868 16778 !3951 12166 1.06 Е-5 1.24 Š— 5 138 Š— 4 1.60 Š— 5 1.79 Š— 5 2.50 Š— 5 2.65 Š— 5 2.98 Š— 5 1,00 Š— 8 5.96 Š— 8 1.12 Š— 6 7,45 Š— 8 1,49 Š— 8 2,98 Š— 8 3.5. Реализаиия стратегии выбора В этом параграфе мы обсудим некоторые практические вопросы, связанные с реализацией стратегии выбора. Будет подробно описана стратегия„воплощенная в программе У12М 197, 1О2~, в основе которой — 1ЙМЯ с малым значением р, Кроме того, мы проведем сопоставление с МА28 1171, использующей 6МЗ при р = и, На 1-м шаге гауссова исключения мы должны просмотреть р (или, точнее, гшп(р, п — Е+ 1)) строк с наименьшим числом ненуле~ы~ эл~ме~~ов, Вместо *ого чтобы Обходить все и — М+ 1 строк в поисках р ~лучших», мы упорядочим строки по возрастанию числа элементов и будем рассматривать р ~первых~.
Для хранения и обра.ботки этой информации нужны три массива длины и: Один для хранения информации и два для ее обновления. В 712М используются три столбца целого массива НА (см. ~ 2,2): 7-й, 8-й и 11-й. Столбец НА(, 7) содержит номера строк, упорядоченных по возрастанию числа элементов. Эту информацию необходимо обновлять после каждого шага исключения, так как мы должны удалить ведущую строку, удалить элементы ведущего столбца и добавить Элементы за~ол~ения.
В НА(1,8) мы храним позицию 1-й строки в упорядоченном списке НА(.,7), а в НА(1, 11) — позицию в НА(,7) первой строки с 1 элементами. Если такой стрОки нет, полагаем НА(1, 11) = О. Пос~о~~~у на од~~~ шаге исключен~я изменяются лишь несколько строк, то посредством этих трех массивов (Зп позиций) мы можем очень эффективно следить за числом ненулевых элементов в строках А~ и поддерживать соответствующее упорядочение строк. Поэтому число Операций 1~ будет величиной порядка 0(пс), где с — среднее число строк, абрабатываемых за шаг„по сравнению с О(п') операциями в случае, когда каждый раз ведется поиск.