Итерационное уточнение решения линейной системы.

Лабораторная работа №5.

Решение систем линейных уравнений.

Часть 1.

Содержание.

Введение 3.

5.1 Метод Гаусса 3.

5.2 Погрешность решения систем. Нормы матриц. Обусловленность матриц. 5.

5.3 Итерационное уточнение решения линейной системы 8.

Задачи для самостоятельного решения 10.

Задания к лабораторной работе 12.

Приложение 1 13.

Литература 15.

Методы численного решения систем линейных уравнений принято разбивать на две группы, К первой группе относят так называемые точные или прямые методы - это алгоритмы позволяющие получить решение за конечное число арифметических операций. Сюда относится метод Крамера, метод Гаусса и метод прогонки[1]. Другую группу состав­ляют приближенные методы, в частности итерационные методы решения систем линейных алгебраических уравнений.

Правило Крамера, при решении задач на ЭВМ не применяется, так как требует большого числа арифметических операций, и, следовательно, дают большую вычислительную погрешность. Метод Гаусса, как правило, используется для решения систем до порядка , а итерационные методы – до порядка .

Мы уже рассматривали системы уравнений, но эти системы имели низкий порядок и системы задавались точно аналитически. Для систем линейных уравнений эта ситуация не имеет места и поэтому крайне важной становится задача оценки погрешности решения системы, по погрешности поставленной задачи и проблема оптимизации алгоритма, с целью уменьшения вычислительной погрешности.

Метод Гаусса.

Этот метод наиболее часто применяется для решения систем т.к. он имеет два важных достоинства: выбор ведущего элемента, позволяющий минимизировать погрешность и корректный учет влияния ошибок округления.

В общем случае гауссово исключение состоит из двух этапов: прямого хода и обратного хода. Прямой ход состоит из шагов. На -том шаге кратные -того уравнения вычитаются из оставшихся уравнений с целью исключить -тое неизвестное. Обратный ход состоит в решении последнего уравнения относительно , затем предпоследнего уравнения относительно и т.д., пока из первого уравнения не будет вычислено .

Мы всегда будем использовать для системы следующие обозначения:

(1_a)

или в матричном виде:

, (1_b)

где , , .

Алгоритм в целом может быть компактно записан в матричных обозначениях. Для системы (1) положим:

,

и умножая равенство (1) слева на матрицу получим систему:

(2).

Это преобразование выполнимо, если элемент , называемый ведущим, не равен нулю.

В этом состоит первый шаг прямого хода метода Гаусса. Второй шаг состоит в умножении равенства (2) на матрицу

.

Выполняя шагов прямого хода[2], получим систему с матрицей верхнего треугольного вида, у которой на главной диагонали стоят единицы

.

Замечание. При выполнении каждого из шагов прямого хода можно уменьшить вычислительную погрешность, для этого в качестве ведущего элемента следует выбирать максимальный по модулю элемент разрешающего столбца, находящийся ниже элемента стоящего на главной диагонали[3]. Говоря другими словами, при каждом шаге на диагональ переставляется уравнение с максимальным по модулю коэффициентом при исключаемой неизвестной.

Далее выполняя обратный ход, получим решение системы. Т.е. для реализации этого алгоритма на ЭВМ достаточно иметь процедуры: перемножения матриц, перестановки строк матрицы и построения матрицы. Для законченности алгоритма необходимо уметь подсчитывать погрешность решения системы, чем мы и займемся в следующем параграфе.

5.2 Погрешность решения системы. Нормы матриц. Обусловленность систем.

Поскольку мы договорились оценивать погрешность решения как расстояние между точным и приближенным решениями системы, то нам потребуется определить норму матрицы, и ее определим, согла­сованно с определением нормы вектора:

(1).

Отсюда непосредственно следует, что (1) удовлетворяет аксиомам нормы

, для [4].

Таким образом, множество матриц может быть интерпретировано как нормированное векторное пространство размерности . Также из (1) следует, что

. (2)

Более того, неравенство (2) является “неулучшаемым”, в том смысле, что для каждой матрицы найдется вектор , такой, что

.

Квадратная матрица задает линейное преобразование каждого вектора в вектор , можно считать, что формула (1) определяет как наибольшую длину образа вектора единичной сферы при преобразовании .

Если выбрать евклидову норму в качестве нормы вектора, то можно интерпретировать норму матрицы следующим образом: При выполнении линейного преобразования единичная сфера деформируется в эллипсоид[5], тогда норма матрицы , представляет собой длину максимальной полуоси эллипсоида. Или, что норма диагональной матрицы равна максимальному по модулю элементу диагонали. Также следует, что умножение матрицы на любую ортогональную матрицу не изменяет ее нормы.

Мы интерпретировали только одну из возможных норм векторов. Две другие, обычно используемые в численном анализе нормы таковы:

.

Эти две нормы порождают соответствующие матричные нормы:

.

Предлагается, в качестве упражнения, интерпретировать их в двух и трехмерном пространствах.

Для расчета погрешности решений системы будем вначале полагать, что матрица системы задана точно, а вектор правых частей приближенно и выясним, как велико может быть . Имеем следующую цепочку соотношений:

,

так как , то . Перемножая, обе части последних неравенств, имеем:

.

Предполагая, что , получим

.

Для любой невырожденной матрицы определим сейчас ее число обусловленности, обозначаемое , как произведение [6] и предыдущая формула перепишется в виде

. (3)

Теперь это соотношение можно интерпретировать следующим образом: есть мера относительной неопределенности в задании вектора . Аналогично можно интерпретировать как меру относительной неопределенности вектора решения системы. Поскольку неравенство является не улучшаемым, то правая часть (3) есть предельная относительная погрешность.

До сих пор мы предполагали, что матрица A известна точно, вектор правых частей имеет неопределенность. Предположим теперь, что матрица A известна приближенно

,

тогда , после этого, воспользовавшись тождеством

,

найдем . Переходя в этой формуле к нормам: получаем

.

Объединяя полученные результаты, мы будем иметь соотношение для погрешности решения системы вследствие представления данных в арифметике с плавающей запятой[7]

,

где -вычисленное решение, -основание системы счисления, -длина разрядной сетки мантиссы. Следует заметить, что этот вывод справедлив, если мы пренебрегаем погрешностью возникающей при выполнении арифметических операций. Это предположение оправдано тем обстоятельством, что если эта погрешность и возникает, то она минимизируется выбором ведущего элемента при прямом ходе метода Гаусса.

Все результаты рассмотрения можно сформулировать в окончательном виде: Предположим, что все элементы и точно представляются в арифметике с плавающей запятой, и пусть вектор , получен на выходе подпрограммы решения системы линейных уравнений. Предположим также, что точная вырожденность матрицы, если она имеет место, не была обнаружена, и что не было ни машинных нулей, ни переполнений. Тогда имеют место следующие неравенства:

,

.

Здесь - основание системы счисления арифметики, - длина разрядной сетки мантиссы.

Второе неравенство утверждает, что, как правило, можно рассчитывать на то, что относительная невязка будет иметь величину, сравнимую с ошибкой округления, независимо от того, насколько плохо обусловлена матрица.

Итерационное уточнение решения линейной системы.

Несмотря на то, что ошибки округления, как правило, невелики, иногда возникает необходимость получить решение с более высокой точностью. Ключом для улучшения точности первого решения является вычисление с двойной точностью его невязки . Зная , мы затем решаем систему

.

Если мы знаем точно, то будет точным решением системы, потому что

.

Если мы не знаем точно, то можно вычислить и продолжить процесс уточнения решения: зная , мы затем решаем систему

.

Заметим, что каждая система имеет ту же матрицу и поэтому прямой ход метода Гаусса можно выполнить лишь один раз, если мы представим его в виде умножения на некоторую матрицу . Если не производить выбор ведущего элемента, то

.

Как нетрудно видеть перестановку -той и -той строк можно свести к умножению на матрицу у которой по диагонали, за исключением этих строк, стоят единицы, а в -той строке на -том месте стоит единица и в -той строке на -том месте стоит единица, другие элементы этой матрицы равны нулю. С учетом этого обстоятельства матрицу можно записать в виде:

.

Исключительно важно, чтобы вычисление невязок выполнялось с более высокой точностью, чем остальные операции. Это основной принцип для любых способов решения уравнения: вычисление невязки является критической операцией, и должно выполнятся точнее других операций.

Детальный анализ ошибок округления дан в предыдущем параграфе, сходимость процесса итерационного уточнения решения мы рассмотрим чуть позже, а сейчас дадим некоторое правдоподобное представление о природе алгоритма. В частности, мы покажем, что величина приближенно определяет быстроту сходимости к своему предельному вектору. Обратно, последовательность может быть использована для оценки величины .

Замечание. Ошибка решения системы связана с простым соотношением , таким образом

.

По результатам предыдущего параграфа имеем

.

Используя предыдущее равенство и замечание, имеем:

.

Следовательно, при

.

Если , то можно положить , в другом случае

,

т. е. у приблизительно первых разрядов верны, а остальные в пределах погрешности.

После вычисления , мы можем предположить, что невязка имеет - верных разрядов. Когда мы решаем систему уравнений для , то получим еще - верных разрядов. Следовательно, имеет верных разрядов.

Задачи для самостоятельного решения.

№1.

Доказать, имеет место неравенство .

№2.

Доказать, что имеют место равенства:

а) ,

б) ,

где нормы матриц согласованы с соответствующими нормами векторов.

№3.

Пусть и имеют следующий вид:

, .

Собственные значения матрицы приближенно равны 0.0558, 0.2007, 84.74.

а) Описать множество , т.е. образ единичной сферы при преобразовании .

б) Найти , , .

в) Рассмотрим систему уравнений . Предположим, что мы имеем векторы и , относительно которых известно только, что

.

1) Какова наименьшая верхняя граница для абсолютной ошибки .

2) Какова наименьшая верхняя граница для относительной ошибки .

№4.

Пусть

.

Доказать, что данная матрица имеет наибольшее число обусловленности [8] из всех невырожденных матриц второго порядка, элементами которых являются положительные числа меньшие или равные 100.

№ 5.

Пусть симметричная положительно определенная матрица. Доказать, что есть монотонно убывающая функция , при .

№ 6.

Найти евклидову норму диагональной матрицы.

№ 7.

Найти евклидову норму ортогональной матрицы.

№ 8.

Показать, что на основании следующего алгоритма можно получить оценку числа обусловленности матрицы, по норме, согласованной с векторной . Для расчета , находят вектора и , для этого решаются две системы уравнений

,

,

где - вектор с компонентами , выбираемый так, чтобы максимизировать в процессе обратной подстановки для , причем .

Задания к лабораторной работе.

Решить системы уравнений, с максимально возможной точностью, для данных одинарной точности :

.

Основная матрица системы порядка одинакова для всех вариантов и ее элементы равны .

Вектор правых частей, в каждом варианте получается расчетом значений функции , для , с шагом .

Порядок выполнения работы:

1. Составить и отладить процедуры, реализующие прямой ход метода Гаусса.

2. На основании полученных процедур составить подпрограмму расчета матрицы, приводящей матрицу к верхнему треугольному виду.

3. Составить и отладить подпрограмму расчета числа обусловленности матрицы и провести расчет.

4. Выполнить обратный ход метода Гаусса.

5. Провести расчет невязки и уточнить полученное решение.

Приложение 1.

Вычисление нормы матрицы согласованной с евклидовой нормой вектора.

Мы уже отмечали, что величина нормы матрицы и число обусловленности зависят от выбора типа нормы. Если для норм и их расчет производится по явным формулам достаточно просто, то вычисление числа обусловленности требует знания обратной матрицы[9], что, вообще говоря, эквивалентно решению системы и поэтому простота расчета указанных норм не гарантирует простоту расчета числа обусловленности матрицы.

Для подсчета евклидовой нормы матрицы , с отличным то нуля определителем введем , заметим, что данная матрица является симметричной. Более того, для любого вектора

,

следовательно, матрица положительно определена. Это означает, что все ее собственных значений положительны, так что мы можем их обозначить через , причем можно положить

,

эти числа называют сингулярными числами матрицы . Из теории вещественных симметрических матриц следует, что можно найти ортогональную матрицу , такую, что

,

где - диагональная матрица с элементами . Теперь определим матрицу как

.

Тогда , так что матрица является диагональной. Более того, последнее равенство показывает, что различные строки матрицы ортогональны друг другу . Построим ортонормированную систему вектор-строк , для . Пусть - матрица строками которой являются вектора

,

тогда и следовательно . Тем самым доказана:

Теорема. Для любой невырожденной матрицы существуют две вещественные ортогональные матрицы и , такие, что диагональная матрица . Более того, можно выбрать и так, чтобы диагональные элементы имели вид:

.

Из этой теоремы следует, что и норма , следовательно . Т.е., для вычисления числа обусловленности матрицы достаточно найти максимальное и минимальное собственные значения матрицы и окончательно . Укажем итерационный алгоритм, принадлежащий Л. А. Люстернику, который позволяет найти наибольшее собственное значение симметричной положительно определенной матрицы . Рассмотрим итерационный процесс

, причем - произвольный ненулевой вектор. Покажем, что

.

Для этого представим в виде разложения по полной системе собственных векторов матрицы ,

.

После -той итерации будем иметь:

.

Перепишем последнее равенство в виде:

,

где в первой сумме выписаны собственные векторы, принадлежащие собственному значению . На основании этого равенства имеем

.

Переходя к пределу, при получим:

.

Переходя в этом соотношении к норме, получим искомый результат:

. (1)

Минимальное собственное значение матрицы можно получить, если найти максимальное собственное значение матрицы , которое связано с соотношением

.

Замечание. При реализации итерационного процесса основанного на формуле (1), для сокращения числа арифметических операций дробь в правой части формулы можно вычислять не для всех значений , а с некоторым разумным интервалом, например через каждые десять шагов итерации.

Литература.

1. Н.С. Бахвалов, Н.П. Жидков, Г.М. Кобельков. Численные методы. М:, Наука, 1987.

2. Дж. Форсайт, М. Малькольм, К. Моулер Машинные методы математических вычислений. М., “Мир” 1980.

3. В.В. Воеводин Вычислительные основы линейной алгебры. М., “Наука”, 1977.

4. Е.А. Волков Численные методы. М., “Наука”, 1987.

6. Дж. Форсайт, К. Моулер Численное решение систем линейных алгебраических уравнений. М. “Мир”, 1969.

7. В. В. Смелов. Основы методов вычислительной математики. НГУ 1986.

[1] Вариант метода Гаусса, для систем специального вида.

[2] Если при выполнении -того шага =0, то перестановкой строк с номером большим можно добиться чтобы 0, если система однозначно разрешима.

[3] Поскольку при делении на малое число происходит рост погрешности.

[4] -матрица состоящая из нулей.

[5] Мы предполагаем, что преобразование невырождено.

[6] Таким образом, число обусловленности зависит от используемой нормы.

[7] Естественно мы предполагаем, что ошибка возникает только вследствие приближенного представления системы.

[8] Это обозначение показывает, что число обусловленности подсчитывается по норме согласованной с евклидовой.

[9] Обычно вместо вычисления числа обусловленности получают только его оценку.

Не нашли, что искали? Воспользуйтесь поиском по сайту: