Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 11
Текст из файла (страница 11)
д.Алгоритм с отложенными модификациями не столь нужен для тех машинтипа «регистр–регистр», в которых допускается совмещение с арифметикойкак загрузок, так и записей в память. В этом случае операцию (3.7) можнобыло бы удалить, а операцию (3.8) заменить операцией Модифицировать третий столбец матрицы A;загрузить четвертый столбец матрицы A;записать в память второй столбец матрицы A.Таким образом, запись в память второго столбца матрицы A происходитодновременно с загрузкой четвертого столбца.Замечание 3.1. Материал подразд. 3.2, 3.3 и 3.4 из [8] приведен,чтобы объяснить те приложения, для которых создаются алгоритмы с отложенными модификациями: алгоритм на рис.
3.6 и другие, показанные нижев подразд. 3.5. Таким образом, включение таких алгоритмов в лабораторнуюработу (проект) может рассматриваться как задача имитирования алгоритмов векторных или параллельных компьютеров на обычных (скалярных)компьютерах с целью освоения этих современных версий LU -разложения.3.5LU -разложение и его ijk-формыНиже в описании ijk-алгоритмов факторизации (разложения), основанных на методе Гаусса исключения переменных, используем из [8] следующиеобозначения для индексов:k — номер исключаемой переменной,i — номер строки, т. е. модифицируемого уравнения,j — номер столбца, т. е. коэффициента в модифицируемом уравнении.Тогда общую основу всех алгоритмов удобно определить тройкой вложенных циклов видаДляДляДляaij = aij − lik akj623.5 LU-разложение и его ijk-формыЗдесь последняя формула обозначает модификацию j-го элемента i-йстроки матрицы A при исключении k-й переменной вектора неизвестныхx из уравнений системы Ax = f .
Перестановки трех индексов для цикловопределяют 3! = 6 возможных вариантов алгоритмов, образуя так называемые ijk-формы, для каждого вида разложения. Для квадратной матрицыA размера (n × n) возможны четыре вида разложения, а именно:A = LŪ, A = L̄U, A = U L̄, A = Ū L,(3.10)где черта сверху указывает на тот из сомножителей, который имеет единичную главную диагональ.
Поэтому всего возможно построить 24 вариантаijk-алгоритмов разложения матрицы для решения различных задач: решения систем, обращения матрицы и вычисления ее определителя.Рассмотрим все шесть ijk-форм для одного из четырех разложений (3.10),а именно, для L̄U -разложения матрицы A. Для численной иллюстрациивозьмем следующий пример.Пример 3.3.2 4 −4 611 4 2 11 , L̄ = 1/2A=3 8 1 1 3/2112 5 0 51 1/2 2/3 12 4 −4 62 4 −2 , U = .3 −6 4Два алгоритма для L̄U -разложения матрицы Aс немедленными модификациями1) kij-алгоритм, рис. 3.1.Доступ к элементам матрицыA по строкам. Исключение постолбцам. Модификации немедленные.
ГЭ — любая из трехстратегий.Для k = 1 до n − 1Для i = k + 1 до nlik = aik /akkДля j = k + 1 до naij = aij − lik akj2) kji-алгоритм, рис. 3.2.Доступ к элементам матрицыA по столбцам. Исключение постолбцам. Модификации немедленные. ГЭ — любая из трехстратегий.Для k = 1 до n − 1Для s = k + 1 до nlsk = ask /akkДля j = k + 1 до nДля i = k + 1 до naij = aij − lik akj633 Векторно-ориентированные алгоритмы LU-разложенияДва алгоритма для L̄U -разложения матрицы A(столбцово ориентированные с отложенными модификациями)3) jki-алгоритм, рис.
3.6.Доступ к элементам матрицыA по столбцам. Исключение постолбцам. Модификации отложенные. ГЭ — по (j − 1)-мустолбцу.Для j = 2 до nДля s = j до nls,j−1 = as,j−1/aj−1,j−1Для k = 1 до j − 1Для i = k + 1 до naij = aij − lik akj4) jik-алгоритм.Доступ к элементам матрицыA по столбцам. Исключение постолбцам. Модификации отложенные. В цикле по s идет нормировка (j −1)-го столбца. Первый цикл по i вычисляет столбец для U , второй — столбецдля L̄. ГЭ — по (j − 1)-мустолбцу.Для j = 2 до nДля s = j до nls,j−1 = as,j−1/aj−1,j−1Для i = 2 до jДля k = 1 до i − 1aij = aij − lik akj .Для i = j + 1 до nДля k = 1 до j − 1aij = aij − lik akjДва алгоритма ijk-форм для L̄U -разложения матрицы A(строчно ориентированные с отложенными модификациями)5) ikj-алгоритм.Доступ к элементам матрицыA по строкам.
Исключение построкам. Модификации отложенные. ГЭ — по (i − 1)-йстроке.Для i = 2 до nДля k = 1 до i − 1li,k = ai,k /ak,kДля j = k + 1 до naij = aij − lik akj6) ijk-алгоритм.Доступ к элементам матрицыA по строкам. Исключение построкам. Модификации отложенные. Первый цикл по j находит элементы i-й строки L̄.Второй цикл по j — элементыi-й строки U . ГЭ — по (i − 1)-йстроке.Для i = 2 до nДля j = 2 до ili,j−1 = ai,j−1/aj−1,j−1Для k = 1 до j − 1aij = aij − lik akjДля j = i + 1 до nДля k = 1 до i − 1aij = aij − lik akj643.5 LU-разложение и его ijk-формыUzL̄Uk–я строкаzL̄···...k–я строка···Рис. 3.7.
Способ доступа к данным для kij-формы (слева) и для kji-формы (справа)L̄U-разложения. Обозначения: L̄, U — вычисление закончено, обращений большенет; z — главный элемент (ГЭ); — деление на ГЭ (нормировка) [8]UUL̄L̄z...k∅k–я строка∅Рис. 3.8. Способ доступа к данным для jki-формы и для jik-формы (слева) и для ikjформы и для ijk-формы (справа) L̄U-разложения. Обозначения: L̄, U — вычислениезакончено, обращения больше не производятся; z — главный элемент (ГЭ); —деление на ГЭ (нормировка); ∅ — обращений не было [8]Замечание 3.2. В приведенных алгоритмах не содержится процедура выбора главного элемента. Она дословно переносится из описаниялабораторной работы № 1. Аналогичные алгоритмы могут быть написаныдля остальных трех видов разложения матрицы A из списка (3.10). Принаписании программ, соответствующих приведенным алгоритмам, следуетвыполнить требование, согласно которому все вычисления выполняются водном и том же двухмерном массиве, где сначала хранится матрица A.
Впроцессе вычислений матрица A замещается элементами треугольных матриц, составляющих искомое разложение из списка (3.10). Способ доступа кданным для ijk-форм L̄U -разложения показан на рис. 3.7 и рис. 3.8. Расчеты по алгоритмам kij-формы и kji-формы L̄U -разложения достаточноочевидны. Для других четырех форм L̄U -разложения эти вычисления поясняются для примера 3.3 в табл. 3.1–3.4.653 Векторно-ориентированные алгоритмы LU-разложенияТаблица 3.1. Вычисления по алгоритму jki-формы для примера 3.3.
Позиции ГЭ(без их реального выбора) показаны выделенным шрифтом2132A4 −4428150611521/23/22/2j=24 −42221106115нормировка(j − 1)-гостолбца;(j − 1)-кратноеисключениев j-м столбцеисходнаяматрица21/23/22/2j=34 −4242/271/241/2423/22/22/21/22−443261156115нормировка(j − 1)-гостолбца;(j − 1)-кратнаямодификацияj-го столбца21/23/22/221/23/22/221/23/22/2j=44 −4624 −22/23 −81/22/3−14 −4624 −22/23 −61/22/304 −4624 −22/23 −61/22/34Таблица 3.2. Вычисления по алгоритму jik-формы для примера 3.3. Позиции ГЭ(без их реального выбора) показаны выделенным шрифтом2132A4 −4428150исходнаяматрица66611521/23/22/2j=24 −42221106115нормировка(j − 1)-гостолбца;(j − 1)-кратнаямодификацияj-го столбца21/23/22/2j=34 −4242/271/201/2423/22/22/21/22−44326115611521/23/22/2j=44 −4624 −22/23 −61/22/34нормировка(j − 1)-гостолбца;(j − 1)-кратнаямодификацияj-го столбца3.6 Треугольные системыТаблица 3.3.
Вычисления по алгоритму ikj-формы для примера 3.3. Позиции ГЭ(без их реального выбора) показаны выделенным шрифтом21324485A−42106115исходнаяматрицаi=24 −461/224 −23 8112 5052(i − 1)нормировоки вычитанийв i-й строке21/23/2221/23/22i=34 −4624 −227 −85054 −4624 −22/23 −6505(i − 1)нормировоки вычитанийв i-й строке3.621/23/22/2i=44 −4242/23146−2−6−14 −4242/231/226−2−603/24 −4242/232/21/26−2−6421/23/22/221/22/3Треугольные системыПо окончании этапа приведения в гауссовом исключении нам необходиморешить треугольную систему уравнений c1x1u11 u12 · · · u1nu22 · · · u2n x2 c2 .. = ..
...... . . . cnxnunnОбычный алгоритм обратной подстановки описывается формуламиxi = (ci − ui,i+1xi+1 − . . . − uinxn )/uii,i = n, . . . , 1.(3.11)Рассмотрим, как он может быть реализован в векторных операциях. ЕслиU хранится по строкам (так будет, если на этапе приведения A храниласьпо строкам), то формулы (3.11) задают скалярные произведения с длинамивекторов, меняющимися от 1 до n−1, и n скалярных делений (рис. 3.9 слева).Альтернативный алгоритм, полезный, если U хранится по столбцам,представлен в виде псевдокода на рис.
3.9 справа. Он называется столбцовым алгоритмом (или алгоритмом векторных сумм). Как только найдено673 Векторно-ориентированные алгоритмы LU-разложенияТаблица 3.4. Вычисления по алгоритму ijk-формы для примера 3.3. Позиции ГЭ(без их реального выбора) показаны выделенным шрифтом21324485A−4210исходнаяматрица6115i = 2, j = 2, k2 4 −41/2223 812 50=16115i = 3, j = 2, k = 124 −461/224 −23/22112505i = 4, j = 2, k = 124 −461/224 −23/22/23 −62/2105i = 2, j = 3, k2 4 −41/2243 812 50=16115i = 3, j = 3, k = 124 −461/224 −23/22/2712505i = 4, j = 3, k = 124 −461/224 −23/22/23 −62/21/245i = 2, j = 4, k = 12 4 −461/224 −23 8112 505i = 3, j = 3, k = 224 −461/224 −23/22/2312505i = 4, j = 3, k = 224 −461/224 −23/22/23 −62/21/225(i − 1)нормировокв i-й строкеi = 3, j = 4, k = 124 −461/224 −23/22/23 −82505i = 4, j = 4, k = 124 −461/224 −23/22/23 −62/21/22/3−1i = 3, j = 4, k = 2i = 4, j = 4, k = 2(n − i)(i − 1) +i(i − 1)+2вычитанийв i-й строке1/2423/22/2252−44306−2−651/2423/22/2−4432/21/22/326−2−60i = 4, j = 4, k = 3(i − 1)нормировок и(n − i)(i − 1) +i(i − 1)+2вычитанийв i-й строке681/2423/22/2−4432/21/22/326−2−643.7 Задание на лабораторный проект № 2Для j = n до 1 с шагом −1Для j = i + 1 до nci = ci − uij xjxi = ci /uiiДля j = n до 1 с шагом −1xj = cj /ujjДля i = j − 1 до 1 с шагом −1ci = ci − xj uijРис.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
