Вычислительные методы алгебры и оценивания. И.В. Семушкин (2011) (1185350), страница 18
Текст из файла (страница 18)
В массивеa хранятся построчно элементы матрицы от первого ненулевого до диагонального включительно. В массиве b на i-м месте стоит положение i-го диагонального элемента матрицы P в массиве a. Для определения положения(i, j)-гo элемента матрицы P в массиве a надо воспользоваться следующималгоритмом. Сначала вычисляем k = b(i) − (i − j). Затем, если k > b(i − 1),1016 Разложения Холесскогото этот элемент стоит на k-м месте. В противном случае (i, j)-й элементстоит левее первого ненулевого элемента i-й строки, поэтому он равен нулюи в массиве a не хранится. Способ хранения по столбцам строится аналогичным образом, но в этом случае хранятся все элементы от диагональногодо последнего ненулевого элемента столбца включительно.Таким образом, существуют 4 варианта хранения разреженной ленточнойматрицы P , и выбор конкретного варианта должен соответствовать заданному варианту разложения Холесского и разновидности ijk-форм.Замечание 6.7.
С учетом положительной определенности матрицэтой лабораторной работы, процедура выбора главного элемента, а такжепроцедуры перестановки строк и столбцов матрицы P отсутствуют как длязаполненных, так и для разреженных матриц.6.7Задание на лабораторный проект № 5Написать и отладить программу, реализующую заданный вариант алгоритма разложения Холесского, для решения системы P x = f , где P —симметричная положительно определенная матрица (ПО-матрица P , иликратко, матрица P > 0). Отделить основные части программы:а) подпрограмму генерации ПО-матрицы P ;б) подпрограмму разложения Холесского;в) подпрограмму решения систем линейных алгебраических уравнений;г) подпрограмму решения систем линейных алгебраических уравнений сленточной матрицей P ;д) сервисные подпрограммы, включая демонстрацию разложения наэкране, подпрограмму контроля правильности разложения и др.Уделить особое внимание эффективности программы (в смысле экономииоперативной памяти).
Для этого в одномерном массиве достаточно хранитьтолько нижнюю (или верхнюю) треугольную часть и диагональ матрицыP . Результат разложения замещает исходную матрицу. Предусмотретьпошаговое выполнение алгоритма исключения с выводом результата наэкран.Выполнить следующие пункты задания:1. Дать формулировку и доказательство теоремы о заданном вариантеразложения Холесского. Теорема может быть сформулирована как утверждение о том, что предлагаемый студентом алгоритм действительно дает1026.7 Задание на лабораторный проект № 5единственное разложение Холесского. Для иллюстрации дать численныйпример работы алгоритма по шагам для матрицы P размера (4 × 4).2.
Провести подсчет количества операций:а) извлечения квадратного корня;б) умножения и деления.Подсчет выполнить тремя способами:а) фактически — при конкретном расчете разложения;б) теоретически точно в зависимости от размерности матрицы n;в) теоретически приближенно в зависимости от n при n → ∞.3. Определить скорость решения задач (решение систем линейных алгебраических уравнений), для чего спроектировать и провести эксперимент,который охватывает сгенерированные случайным образом ПО-матрицы Pпорядка от 5 до 100 (через 5 порядков). Результаты представить в видетаблицы и графика зависимости времени выполнения от порядка матриц.Сравнить со своим вариантом из лабораторного проекта № 1.4. Оценить точность решения систем линейных алгебраических уравнений с матрицами из п. 3.
Для этого выбрать точное решение x∗ и образовать правые части f = P x∗ . В качестве точного решения взять векторx∗ = (1, 2, . . . , n). Для оценки точности использовать норму вектора (2.18) излабораторного проекта № 1. Провести анализ точности вычисленного решения x в зависимости от порядка матрицы. Результаты представить в видетаблицы и графика. Сравнить со своим вариантом из проекта № 1.5. Для заполнения матрицы P использовать случайные числа из диапазона от −100 до 100. Сначала заполнить треугольную часть матрицы P , т. е.элементы pij , где i > j Затем заполнить диагональ. В качестве диагонального элемента pii, i = 1, 2, . .
. , n, выбрать случайное число из интервалаXX|pij | + 101 ,(6.6)|pij | + 1,j6=ij6=iчтобы обеспечить выполнение условияX|pij | + 1,pii ≥j6=iгарантирующего положительную определенность матрицы P .1036 Разложения Холесского6. Определить скорость и точность решения систем линейных алгебраических уравнений с разреженными ленточными матрицами. Для этогоспроектировать и провести эксперимент для систем порядка от 100 до 200(через 5).
Результаты представить в виде таблиц и графиков зависимостискорости и точности решения от порядка матриц. Для этих же систем найтианалогичные зависимости для обычного метода Холесского. Pезультатысравнить.7. Для случайного заполнения разреженной ленточной матрицы Pиспользовать следующий алгоритм:а) случайным образом заполнить половину матрицы (верхнюю или нижнюю), включая диагональ;б) в случае заполнения нижней треугольной части матрицы P в i-йстроке, i = 1, 2, . .
. , n, случайным образом определить количество ненулевыхэлементов (от 1 до 10), их местоположение (номер столбца oт max{1, i − 50}до i − 1) и значение (ненулевые целые числа, лежащие в интервале от −100до 100);в) при заполнении верхней треугольной части матрицы P применять тотже алгоритм, что и в п. б), с той лишь разницей, что номер столбца лежитв интервале от i + 1 до min{i + 50, n};г) диагональный элемент в i-й строке, i = 1, 2, . .
. , n, определить случайным образом на интервале (6.6).В качестве точного решения взять вектор x∗ = (1, 2, . . . , n), n — порядокматрицы P . В случае, если при решении системы P x = f выяснится,что матрица P вырождена (плохо обусловлена), то сгенерировать новуюматрицу того же порядка и решить систему линейных алгебраических уравнений с новой матрицей P и новой правой частью. Для оценки точностирешения использовать норму вектора (2.18) из лабораторного проекта № 1.Замечание 6.8.
По ходу проведения всех численных экспериментовна экран дисплея должны выводиться следующие таблицы.Число вычислительных операций:ПорядокматрицыКвадратные корниабУмножение и делениесабгде а, б, в означает способ вычисления числа действий (см. п. 2).104в6.7 Задание на лабораторный проект № 5Решение систем линейных алгебраических уравненийс заполненной матрицей P :ПорядокматрицыВремяметодГауссаметодХолесскогоТочностьметодГауссаметодХолесскогоТаким образом, в данный проект следует включить работу, выполненную ранее в проекте № 1, чтобы иметь возможность сравнения методаХолесского с методом Гаусса как по времени, так и по точности вычислений.Решение систем линейных алгебраических уравненийс разреженной матрицей P :ПорядокматрицыВремяТочностьЗаполненная Разреженная Заполненная РазреженнаяматрицаматрицаматрицаматрицаЭто означает, что для каждого текущего значения n порядка матрицы Pнеобходимо решать две системы: одну — с заполненной матрицей P (см.
п. 5задания), другую — с разреженной матрицей P (см. п. 7 задания).Замечание 6.9. Необходимо вывести на экран следующие графики:Графики решения систем с заполненной матрицей P :••зависимость времени решения от порядка матриц для методов Гаусса иХолесского;зависимость точности решения от порядка матриц для методов Гауссаи Холесского.Графики решения систем с разреженной матрицей P :••зависимость времени решения от порядка матриц для обычного методаХолесского и с учетом разреженности матрицы P ;зависимость точности решения от порядка матриц для обычного методаХолесского и с учетом разреженности матрицы P .
При построении графиков использовать данные из соответствующей таблицы. Для этого ихнеобходимо записать в текстовый файл.1056 Разложения Холесского6.8Варианты задания на лабораторный проект № 5Как отмечалось в конце подразд. 6.5, всего по данной теме предлагается 40 различных вычислительных схем разложений Холесского, которыеи составляют весь набор вариантов задания на лабораторный проект № 5(см. подразд. 6.7)В табл.
6.1 каждому номеру варианта соответствует своя разновидностьразложения Холесского и свой способ организации вычислений.Таблица 6.1. Варианты задания на лабораторный проект № 5Окаймлениеijk-формыВидразложенияcjkijikikjijkизвестнойчастинеизвестнойчастиabcb12345678910P = LLT111213141516171819202122232425262728293031323334353637363940P = UU TbkjiP = L̄D L̄TP = ŪD ŪakijT— строчный алгоритм;— алгоритм скалярных произведений;— алгоритм линейных комбинаций.Если нет других указаний преподавателя, выбирайте ваш вариант по вашему номеру в журнале студенческой группы.1067Ортогональные преобразования7.1Ортогональные матрицы и приложенияВ этом разделе напомним определение и некоторые свойства ортогональных матриц, полезные для дальнейшего.Определение 7.1.
Матрица T , имеюшая размер (n × n), т. е. T (n, n),есть ортогональная матрица, когда T T T = I.Свойство A. Если T1 и T1 суть две ортогональные матрицы, то ихпроизведение T1 T2 есть тоже ортогональная матрица.Свойство B. T −1 = T T и T T T = I.Свойство C. Ортогональное преобразование сохраняет скалярное произведение векторов, т. е. ∀x, y ∈ Rn : y T x , (x, y) = (T x, T y), в частности,оно сохраняет (евклидову) норму вектора: kT yk = kyk.Свойство D. Если v есть вектор случайныхс математиче Tпеременныхским ожиданием E {v} = 0 и ковариацией E vv = I, то теми же характеристиками обладает вектор v̄ = T v, т. е.E {v̄} = 0,E v̄v̄ T = I.Хотя это свойство легко проверяется, немного удивительно, что компоненты преобразованного вектора остаются взаимно некоррелированными.Свойства C и D играют существенную роль в квадратно-корневых алгоритмах решения прикладных задач оптимального моделирования и оптимального оценивания методом наименьших квадратов (рис.
7.1).7 Ортогональные преобразованияzA••СистемаМодель(a)v = z − AxAx−xx̂kvk2 → minxНаблюдениеxAx•ОценкаvAx̃z+•x̃x̂E {kek2 } → min(b)x̃e = x̃ − xНеизвестный вектор−Рис. 7.1. Алгебраически эквивалентные задачи, решаемые методом наименьшихквадратов значений невязки v или среднего квадрата погрешности e: (a) — оптимальное моделирование неизвестной системы по экспериментальным условиям A иданным z; (b) — оптимальное оценивание неизвестного вектора по наблюдениямAx Tв присутствии случайных помех v с характеристиками E {v} = 0 и E vv = I1087.2 Линейная задача наименьших квадратов7.2Линейная задача наименьших квадратовЛинейная задача наименьших квадратов (см.
рис. 7.1) ставится следующим образом (см. также подразд. 11.2) [13, 15].Дано линейное уравнениеz = Ax + v,(7.1)в котором известны вектор z ∈ Rm и (m × n)-матрица A ∈ Rm×n , т. е.A = A(m, n). Разностный вектор v , z − Ax, называемый невязкой, зависитот переменного вектора x ∈ Rn . Требуется найти значение x̂ вектора x,минимизирующее квадратический критерий качестваJ(x) = (z − Ax)T (z − Ax) = kvk2 → min .(7.2)Если ни при каком x невязка v не может быть обращена в 0 — нулевойвектор, то система Ax = z — несовместная, в противном случае совместная,т.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
