Джордж, Лю - Численное решение больших разреженных систем уравнений (947498), страница 14
Текст из файла (страница 14)
Назначение подпрограммы — вычислить компоненты введенного и разделе 4.4.1 массква ХЕХЧ, который используется при хранении множителя Е матрицы РАРТ. Выдается также значение ЕХУЫЕ, представляющее собой размер оболочки Ь и равное аа Гм е. Ленточные и нрофинвнвве методы ХЕНЧ (ХЕЯХЗ + 1) — 1. Здесь, как и раньше, ЫЕЯ)в)З— число уравнений или узлов. Подпрограмма проста и почти не нуждается в пояснениях,. Цикл 00 200 1= ... обрабатывает каждую строку; индекс первого ненулевого элемента рй строки (1Р1КЗТ) матрицы РАР' определяется циклом ОО!00 3 = .... В конце каждого исполнения этого цикла соответствуюшим образом изменяется значение переменной Е)ч!ЧЗЕЕ. Использование РЕЯМ и 1ИЧР связано с тем, что пара массивов (ХА03, АО,)НСУ) хранит структуру А, а структура Ь, которую нужно найти, соответствует матрице РАР'.
й 4.5. Подпрограммы численного решения ЕЗГСТ, Е).З).Ч и Е()З).Ч В этом разделе мы опишем подпрограммы, осуществляющие численное разложение и решение треугольных систем. Они используют профильную схему хранения, о которой говорилось в разделе 4.4.!. Подпрограммы для решения треугольных систем — ЕЕЗЕЧ (Епче!оре-(.очеег-Зо(.Че) и Е~/З(.Ч (Епче!ореУррег-Зо(.Че) — будут описаны до подпрограммы разложения ЕЗРСТ (Епче!оре-Зупнпе1Пс-РаСТог(хаПоп), поскольку последйяя использует подпрограмму Е(.ЗЬЧ, 4.5.1. Подпрограммы для решения треугольных систем Е(З(.Ч и ЕЫЗ1Ч Эти подпрограммы выполняют численное решение нижней и верхней треугольных систем ~у-й и тт соответственно, где Т,— нижняя треугольная матрица, хранимая в соответствии с описанием в разделе 4.4.1.
В подпрограмме Е(.З1Ч имеются несколько важных особенностей, которые заслуживают разъяснения. Прежде всего определяется номер (1Р1КЗТ) первого ненулевого элемента в правой части (КНЗ). Затем программа выполняет цикл (00500 1 = ...) по строкам Т, с номерами 1Р1КЗТ, 1Р!КЗТ -1- 1, ..., ЫЕЯЫЗ; при этом используется схема скалярных произведений (см. раздел 2.2.1). Однако программа пытается использовать последовательности нулевых компонент решения.
Переменная ЕАЗТ хранит индекс последней (на данный момент) вычисленной ненулевой компоненты решения. (Компо- б 4.б. Подпрограммы ЕЯГСТ, ЕьБ1Ч и ЕУоьЧ Вр ненты решения записываются на места соответствующих компонент правой части в массив КНЗ.) Настоятельно рекомендуем читателю смоделировать работу подпрограммы для примера, изображенного иа рис. 4.5.1. Проверьте, что ЕЕЯЕЧ реально использует лишь ненулевые элементы, обозначенные звездочкой в квадрате.
:в Щ ЮБЯ БЯКК ЕЯ яа Решение Ф~юУазр еаплтз <ХЕ)ЧЧ, ЕМЧ) <инз) <янз) рнс. 4ЛЛ. Элементы ь, реально используемые подпрограммой Е).о1Ч, обозна. чены звездочкой в квадрате, Заметим, что 1АЬТ позволяет опустить операции с некоторыми строками; мы все же выполняем умножения с нулевыми операндами внутри цикла РО 300 К = ... по той причине, что для большинства машин проверка, помогающая избежать такого умножения, стоит дороже, чем его выполнение, Проверка и подгонка переменной 1ВАНР, предшествующие циклу РО 300 '..., также требуют некоторых объяснений.
В некоторых обстоятельствах ЕЕБЕЧ используется для решения нижней треугольной системы, матрица козффициентов ко. торой является лишь нодматрицей матрицы 1., переданной подпрограмме парой массивов (ХЕНЧ, ЕМЧ), Ситуация изображена на рисунке. Некоторые строки оболочки простираются вовне используемой подматрицы, и 1ВАИР подгоняется так, ВО т"л.
4. Ленточные и нрофильные методы чтобы учесть это обстоятельство. В примере на рисунке имеет порядок 16. Если нужно решить систему с указанной че еуррррозиеетрр чрррть йагпрре4ет + ~ нозтрфициенто5 4З Г чь зь ррир ттррррр чржь 4т зт 4т тте зт 4а 4е зь ! подматрнцей порядка 11 и правой частью ВНЗ, следует произвести вызов СА1Л. Е1.Я Ч(11, ХЕ)ч!Ч(5), Е5!Ч, О(АО(5), КНЗ), В подпрограмме ХЕ5!Ч(5) интерпретируется как ХЕ)ч!Ч(1), ХЕг!Ч(6) становится ХЕ5!Ч(2) и т. д. Этот «трюк» используется в подпрограмме ЕЗРСТ, обращающейся к Е! 3!.Ч, Перейдем теперь к описанию подпрограммы Е!)З(.Ч, решающей задачу Етх = у, причем Ь хранится по той же схеме, что и в Е(.З(.Ч. Это значит, что мы имеем удобный доступ к столбцам Ат и разреженность можно использовать полностью, применяя схему внешних произведений (см.
раздел 2.2.1). Столбец с номером ! матрицы 5т используется в вычислениях только в том случае, если Ря компонента решения не равна нулю. Как и ЕЕЗЕЧ, подпрограмму ЕБЗ(.Ч можно использовать для решения верхних треугольных систем, матрицы которых являются подматрицами для Ь, заданной парой массивов (ХЕ!ч!Ч, Е!ч!Ч). Техника здесь аналогична описанной выше. Значение переменной !ВА!ч!О подгоняется соответствующим образом для тех столбцов Ет, которые простираются вовне реально используемой части е,. Все подпрограммы, выполняющие численное решение, еодержат помеченный общий блок ЗРКОРЗ, состоящий из единственной переменной ОРЗ.
Каждая. подпрограмма подсчитывает число проделанных ею операций (умножений и делений) Б 4.5 Подпрограммы ЕБРСТ, Е15(.Ч и Е()5(.Ч 91 С С ПОДПРОГРАННА РЕШАЕТ НИЖНЮЮ ТРЕУГОЛЬНУЮ СНСТЕНУ С СХ «ННЗ. МНОЖИТЕЛЬ С ХРАНИТСЯ В ПРОФИЛЬНОМ С ФОРНАТЕ ° С С С С С С С ИЗМЕНЯЕМЫЕ ПАРАМЕТРЫС ВНБ - НЛ ВХОДЕ СОДЕРЖИТ ВЕКТОР ПРАВЫХ ЧАСТЕЙ. С НА ВЫХОДŠ— ВЕКТОР РЕШЕННЯ С ОРБ - ПЕРЕМЕННАЯ.С ДВОЙНОЙ ТОЧНОСТЬЮ> СОДЕРЖАЩАЯСЯ В ОБЩЕМ БЛОКЕ ЗРКОРЗ. ЕЕ ЗНАЧЕНИЕ С БУДЕТ УВЕЛИЧЕНО НА ЧИСЛО ОПЕРАЦИЯ> с ВЫЖ)ЛНЯЕНЫХ ОТОЙ ПОДПРОГРАММОЙ.
С ° ° ° ° ° > > ° ° ° « ° ° ° ° ° БОВВООТ1МЕ ЕЬЯЛ ( МЕ((МБ, ХЕМУ, ЕМУ, 01АС, ВНБ ) С ° ° «« ° *«« С ВХОДНЫЕ ПАРАМЕТРЫ " ме(гм5 - чнслО УРАВнений. (ХЕМУ, ЕМУ) - НАССИВЬ! ДЛЯ ОБОЛОЧКН !. 01АС - МАССИВ ДЛЯ ДИАГОНАЛИ !.* ОООВЬЕ РВЕС1510М СОСМТ, ОР5 СОМНОМ /БРКОРБ/ ОР5 ВЕАЬ 01АС(!), ЕМУ(!). ВНБ(!), 5 1МТЕСЕВ ХЕНЧ(!), 1, 1ВА>(О, 1РТВБТ, К, К5ТОР, 1 К5ТВТ, С, САБУ, МЕ((МБ с« « ° ° ° > ««« ° >««« ° ««« ° > ° ° « ° Ф«« ° « ° ° ° ° С С С С С НАЙТИ НОМЕР ПЕРВОГО НЕНУЛЕВОГО ЗЛЕНЕНТА В ННЗ И ПОМЕСТИТЬ ЕГО В ХРХИЗТ. 1Р1ВБТ 0 1Р1ВБТ 1Р1ВБТ « ! 1Р ( ВНБ( 1Р1В5Т) .МЕ. О.ОЕО ) 00 ТО 200 1Р ( 1Р1ВБТ .1.Т. МЕОМБ ) 00 ТО 100 ВЕТОВМ !.А5Т 0 ЬАЗТ СОДЕРЖИТ ЛОНЕР ПОСЛЕДНЕЙ ВЫЧИСЛЕННОЙ НЕНУЛЕВОЙ КОМПОНЕНТЫ РЕИЕНИЯ . 100 200 С С С С 00 500 1 1Р1В5Т, МЕ()МБ 1ВАМО ХЕ>(У(1«!) - ХЕНЧ(1) 1Р ( 1ВАМО .СЕ.
1 ) 1ВАМО Б ВН5(1) «1 - 1ВАБ(О ВНБ(1) О.ОЕО 1 - ! С е СТРОКА ОБОЛОЧКИ ПУСТА ИЛИ ВСЕ СООТВЕТСТВУЮВЬИЕ КОМПОНЕНТЫ РЕЮЕНИЯ " ИУЛИ ° с С « ° '>' ° Е«БЬУ ... РЕВЕНИЕ НИЖНЕЙ ТРЕУГОЛЬНОЙ СИСТЕМЫ« С ° «« ° ° ««««ПРОФИЛЬНА(Н МЕТОДОМ ° «« '« ° ° ° ««« 99 Гл е. Ленточные и нрофиленые методы 11' ( !В*МО Е(). 0 .Ой 1.АЯТ 1.Т Ь ) 00 ТО Е00 кзтйт ХЕМА(1+1) — !БАМО кзтОР хенч(1+!) " 1 ВО ЗОО К КБТЕТ, КЯТОР Я Я - ЕМУ(К)"ЕЯЕ(Ь) 1, 1.+1 СОМТ1МОЕ СООМТ ! ВАМО ОРБ ОРБ е СООМТ 1Р ( Б .ЕО.
О.ОЕО ) 00 ТО БОО йнБ(1) Б/Р(Аа(1) Я ОРБ+! 000 ьйт.) СОМТ(МОЕ БЕТОКМ ЕМО 500 4.б.й. Подпрограмма разложения ЕЕРСТ В этом разделе мы обсудим некоторые детали подпрограммы численного разложения ЕОРСТ. Подпрограмма вычисляет разложение Холесского е.(т данной матрицы А, хранимой по профильной схеме (см. раздел 4.4.1). Используется вариант метода Холесского, называемый методом окаймления (см. раздел 2.1.2).
Напомним, что если представить А в виде и соответственным образом увеличивает значение ОРЗ. Поэтому если пользователь желает контролировать число выполняемых операций, он может описать такой же общий блок в своей вызывающей программе и следить за значением ОРО. Переменная ОРЯ была описана как переменная с двойной точностью, чтобы избежать ошибок округления при подсчете числа операций. Наши подпрограммы могут использоваться для решения очень больших систем, так что ОРО вполне может принимать значения, достигающие 1О' или 10', хотя отдельная подпрограмма может увеличивать ОРЗ на сравнительно небольшую величину. На многих машинах прн использовании обычной точности сложение в арифметике с плавающей точкой малого числа (скажем, меньшего чем 10) с числом 1Ое снова дает 105.
(Проверьте это, моделируя шестизначну(о арифметику с плавающей точкой) Использование двойной точности для ОРБ делает ошибку округления при подсчете числа операций чрезвычайно маловероятной, С С' С ° г г г * С 50ВВООТ1МЕ ЕОБ!.Ч ( МЩЧ5, ХЕМУ, ЕМУ, 01АС, ВН5 ) С г ° г и г ° гг ° г 1 МЕ()М5 С С г ° г 100 200 С С С С С С С С С С С С С С С С ф 4 й Подпрограммы ЕБЕСТ, ЕЕБ1Ч и Е(/БЕЧ 93 ЕО5ЬУ * ° ° .
° РЕШЕНИЕ ВЕРХНЕЙ ТРЕУГОЛЬНОЙ СИСТЕМЫ ПРОФИЛЬНЫМ МЕТОДОМ "° ° ° ° ПОДПРОГРАММА РЕШАЕТ ВЕРХНЮЮ ТРЕУГОЛЬНУЮ СИСТЕМУ ОХ = ВНВ. МНОЖИТЕЛЬ Н КРЛНИТСЯ В ПРОФИЛЬНОМ ФОРНЛТЕ . В)(ОДНЫЕ ПАРАНЕТРЫМЕ()МВ - ЧИСЛО УРАВНЕНИЙ. (ХЕМУ, ЕМУ) - МАССИВЫ ДЛЯ ОБОЛОЧКИ О. 01АС - МАССИВ ДЛЯ ДИАГОНАЛИ О, ИЗМЕНЯЕМЫЕ ПАРАМЕТРЫВНВ - МА входе содержит некто~ ЛРАВых члстей. НА ВЫХОДЕ - ВЕКТОР РЕШЕНИЯ . ОР5 - ПЕРЕМЕННАЯ С ДВОйНОИ ТОЧНОСТЬЮг СОДЕРЖАЩАЯСЯ В ОБЩЕМ БЛОКЕ ВРКОРВ ЕЕ ЗНЛЧЕНИЕ БУДЕТ УВЕЛИЧЕНО НА ЧИСЛО ОПЕРЛЦИйг ВЫПОЛНЯЕМЫХ ЭТОЙ ПОДПРОГРЛММОЙ .
ВОЛОВЬЕ РВЕС1510М ОХЗМТ, ОР5 СОННОМ /5РКОРВ/ ОР5 ВЕЛЬ 01АС(1), ЕМУ(!), ВН5(!). 5 1мтесен хемг(!), 1, 1ВАНО, к. к5тОР, к5тнт, !., 1 МЕОМ5 + 1 1 1 — 1 1Р ( 1 ЕО 0 ) ВЕТОК)Ч 1Р ( ВН5(1) ЕО 0 ОЕО ) 60 ТО 100 5 ВН5(1)/ОГАС(1) ВН5(1) . 5 ОР5 ОРВ + ! 000 1ВАМО ХЕ))У( 1+1 ) - ХЕМУ( 1 ) 1Р(!ВАНО Сйт) 1ВЛМО-1-1 1Р ( 1ВАМО Е() 0 ) 60 ТО !00 КВТВТ 1 - 1ВА% КВТОР 1 — 1 Ь ХЕМ'(1+!) — 1ВАМО ОО 200 К КВТВТ, КВТОР ВН5(К) ВН5(К) 5 Е)(У(3.) Ь ! + 1 СОМТПЧОЕ СООМТ 1ВА!ЧО ОРВ ОР5 + СООМТ СО ТО !00 С ° ч ° чч ° ° ° ° ь С Счь ч ЕБРСТ ...