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

Размещено на http://www.allbest.ru/

МИНИСТЕРСТВО НАРОДНОГО ОБРАЗОВАНИЯ РЕСПУБЛИКИ УЗБЕКИСТАН

НАВОИЙСКИЙ ГОСУДАРСТВЕННЫЙ ПЕДАГОГИЧЕСКИЙ ИНСТИТУТ

РЕФЕРАТ

Тема: Алгебраическая проблема собственных значений для матриц специального вида и ее программное обеспечение

НАВОИЙ - 2011

ОГЛАВЛЕНИЕ

ВВЕДЕНИЕ

1. Алгебраическая проблема собственных значений и собственных векторов

1.1 Общая постановка

1.2 Характeристическое уравнение

1.3 Алгебраическая кратность собственного значения

2. Классификация задач на собственные значения

2.1 Полная проблема собственных значений

2.2 Частичная проблема собственных значений

3. Вычислительные методы собственных значений и собственных векторов

3.1 Вычислительные методы полной проблемы собственных значений

3.2 Вычислительные методы частичной проблемы собственных значений

4. Программное обеспечение некоторых алгоритмов нахождения собственных значений и собственных векторов

4.1 Программы на языке С++

4.2 М - файлы для системы MatLab

ЗАКЛЮЧЕНИЕ

СПИСОК ЛИТЕРАТУРЫ

ВВЕДЕНИЕ

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

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

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

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

Отметим, что все предлагаемые ниже методы, кроме метода Леверье (1840 г.) и метода Якоби (1846 г.), появились в тридцатых годах нашего столетия или позднее.

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

1. Алгебраическая проблема собственных значений и собственных векторов

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

1.1 Общая постановка

Пусть - матрица порядка и - обратимая матрица со столбцами . Легко видеть, что равенство эквивалентно системе равенств

, .

Эти равенства подводят нас к важным понятиям собственного значения матрицы и собственного вектора.

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

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

Доказательство. Пусть - линейно независимая система собственных векторов матрицы , соответствующих собственным значениям :

. .

Матрица обратима как матрица с линейно независимыми столбцами.

Пример недиагонализуемой матрицы: . Допустим, что

.

Отсюда

.

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

Теорема. Собственные векторы, соответствующие попарно различным собственным значениям матрицы, являются линейно независимыми.

Пусть - собственные векторы для попарно различных собственных значений матрицы . Пусть . Умножим обе части слева на матрицу :

.

Из данного равенства вычтем предыдущее, умноженное на :

.

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

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

1.2 Характeристическое уравнение

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

.

Число является собственным значением матрицы данная система имеет нетривиальное решение .

Определение. Уравнение относительно называется характеристическим уравнением матрицы . Левая часть этого уравнения есть многочлен степени от , называемый характеристическим многочленом матрицы .

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

,

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

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

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

В частности, - величина, называемая следом матрицы . Обозначение: . В силу формул Виета, след равен сумме всех собственных значений с учетом кратностей. Заметим также, что .

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

1.3 Алгебраическая кратность собственного значения

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

Теорема. Любая комплексная матрица А порядка имеет комплексных собственных значений с учетом алгебраических кратностей совпадают

Теорема. Характеристические многочлены подобных матриц совпадают.

Доказательство. Пусть где - обратимая матрица. Тогда

.

Следствие. Собственные значения и их алгебраические кратности для подобных матриц совпадают.

2. Классификация задач на собственные значения

Перечислим некоторые возможные постановки проблемы собственных значений.

I. Полная проблема собственных значений. Заключается в задаче поиска всех собственных чисел заданной матрицы. В некоторых случаях вычисление также и собственных векторов.

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

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

2.1 Полная проблема собственных значений

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

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

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

Выяснению алгебраической сущности преобразования А.Н.Крылова посвящены работы Н.Н.Лузина [1], [2], И.Н.Хлодовского [1], Ф.Р. Гантмахера [1], Д.К.Фаддеева [1].

2.2 Частичная проблема собственных значений

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

Для построения этих методов используются две, по существу различные, основные идеи.

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

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

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

3. Вычислительные методы собственных значений и собственных векторов

Методы определения собственных значений и собственных векторов матрицы естественным образом делятся на две основные группы: прямые методы и итерационные.

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

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

Пусть

есть характеристический многочлен матрицы А. Рассмотрим последовательность векторов , где

, .

Тогда коэффициенты характеристического многочлена будут удовлетворять системе линейных алгебраических уравнений

.

Оказывается [31], что многие методы (Крылова, Данилевского, Самуэльсона и т.д.) по существу являются методами решения этой системы специальном способом.

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

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

3.1 Вычислительные методы полной проблемы собственных значений

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

Метод непосредственного развертывания

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

 (2.5)

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

Нахождение собственных значений матрицы методом А.M. Данилевского.

Сущность метода Данилевского заключается в приведении векового определителя к так называемому нормальному виду Фробениуса

Если нам удалось записать вековой определитель в такой форме, то, разлагая его по элементам первой строки будем иметь:

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

b_{i j}=a_{i j} - a_{i n-1}a_{n i} / a_{n n-1} при i=1,..,n; j не равно n-1;

b_{i n-1}=a_{i n-1} / a_{n n-1}, при i=1,..,n

Последняя строка построенной матрицы B будет удовлетворять нашим условиям, но не будет подобна матрице A, поэтому проведем еще одно преобразование и получим матрицу С подобную A и сохраняющую последнюю строку:

с_{i j}=b_{i j} при i=1,..,n-2

c_{n-1 j}=a_{n 1}b_{1 j}+a_{n 2}b_{2 j}+ . . . +a_{n n}b_{n j} при j=1,...,n

Таким образом получили матрицу С подобную A и с последней строкой как в матрице привеленного вида Фробениуса. Продолжая, преобразуем аналогично n-1 строку матрицы С и т.д. Допустим, что при преобразовании матрицы A в матрицу Фробенниуса P мы пришли к матрице вида:

причем оказалось, что d_{k k-1}=0. Тогда преобразования методом Данилевского нельзя продолжить. Здесь возможны два случая: 1. Существует элемент d_{k l} отличный от нуля, где l < k-1 , тогда переставляем местами (k-1) и l столбцы и (k-1) и l строки получаем матрицу подобную D для которой возможны дальнейшие преобразования по методу Данилевского. 2. Все элементы k строки равны нулю, тогда полученная матрица имеет вид:

Причем D₂ имеет нормальный вид Фробениуса, а матрицу D₁ можно привести к нему методом Данилевского. Полином же порождающий собственные значения матрицы А, есть произведение аналогичных полиномов для D₁ и D₂, при этом коэффициенты полинома для D₂ определены. Процедура на вход получает матрицу A, а на выходе выдает матрицу C подобную A, при этом С имеет нормальный вид Фробениуса. Флаговая переменная S используется для определения случая когда дальнейшее преобразование не возможно, точнее: S=0 когда алгоритм успешно отработал, S > 0 когда на некотором шаге возник случай описанный в пункте 2, при этом S колличество строк в матрице D₁.

Нахождение собственных значений матрицы методом Леверрье.

Процедура определяет коэффициенты характерестического полинома

матрицы A. Положим:

тогда справедливы формулы Ньютона:

Отсюда получаем линеиную алгебраическую систему:

Из которой шаг за шагом определяются коэффицитенты p₁,...,p_n. Следует заметить, что s_k равен следу матрицы A^k, которая находятся непосредственным перемножением используя алгоритм умножения матриц.

Методика решения задачи

1. Для заданной матрицы А составить характеристическое уравнение (2.5): .

Для развертывания детерминанта можно использовать различные методы, например метод Крылова, метод Данилевского или другие методы [2,9,14].

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

3. Для каждого собственные значения составить систему (2.4):

, ,

и найти собственные векторы .

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

Пример. Найдите собственные числа и собственные векторы матрицы

Решение. Составляем характеристическую матрицу :

Находим характеристический многочлен

Решим характеристическое уравнение

Подбором находим, что один корень уравнения равен . Есть теорема, которая говорит, что если число является корнем многочлена , то многочлен делится на разность , то есть , где- многочлен. В соответствии с этой теоремой многочлен должен делиться на . Выделим в характеристическом многочлене этот множитель :

Находим корни трехчлена . Они равны и 3. Таким образом,

- корень кратности 2 17.7 b, - простой корень. Итак, собственные числа матрицы равны , . Найдем соответствующие им собственные векторы. Пусть, тогда для собственного вектора получаем матричное уравнение

что соответствует системе уравнений

Решаем ее методом Гаусса (раздел "Алгоритм нахождения решений произвольной системы линейных уравнений (метод Гаусса)"). Выписываем расширенную матрицу системы

Первую строку, умноженную на числа и прибавляем соответственно ко второй и третьей строкам

Меняем местами вторую и третью строки

Возвращаемся к системе уравнений

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

Полагаем , находим , . Итак, собственному числу соответствует собственный вектор .

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

что соответствует системе уравнений

Решаем ее методом Гаусса. Выписываем расширенную матрицу

Первую строку умножаем на числа 2 и 3 и прибавляем соответственно ко второй и третьей строкам

Вторую строку умножаем на и прибавляем к третьей

Возвращаемся к системе уравнений

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

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

Ответ: Собственные числа: ,, соответствующие собственные векторы: , .

3.2 Вычислительные методы частичной проблемы собственных значений

Метод итераций

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

Пусть матрица А имеет линейно независимых собственных векторов , и собственные значения матрицы А таковы, что

.

Степенной метод решения частичной проблемы собственных значений.

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

 , (10.1)

Теорема 10.1. Пусть матрица A имеет полную систему из ортонормированных собственных векторов e_i , которым соответсвуют собственные значения ⁽ⁱ⁾ , причем

|⁽¹⁾| > |⁽²⁾| ... |⁽ⁿ⁾| (т.е. вектора занумерованы в порядке невозрастания модуля собственного значения, причем собственное значение ⁽¹⁾ - не кратное).

Тогда итерационный процесс (10.1) сходится к собственному значению ⁽¹⁾, причем

.(10.2).

При этом величины сходятся к собственному вектору e₁ (c точностью до постоянного сомножителя, по модулю равного 1):

.(10.3)

4. Программное обеспечение некоторых алгоритмов нахождения собственных значений и собственных векторов

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

Матрицы широко используются при решении обыкновенных дифференциальных уравнений (ОДУ) и уравнений в частных производных, решении оптимальных задач и т. п.

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

В вычислительном плане раздел линейной алгебры поддержан пакетами прикладных программ LINPACK, EISPACK, разработанными в 60-70-е годы ведущими специалистами, к числу которых принадлежит и основатель фирмы The MathWorks, Inc. Моулер (C. Moler). Изначальное назначение системы MATLAB состояло именно в том, чтобы создать диалоговую среду для работы с пакетами программ линейной алгебры.

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

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

Вычисление собственных значений и сингулярных чисел

· EIG, CDF2RDF - собственные значения и собственные векторы матрицы (Matritsaning xos son va xos vektorlari)

· BALANCE - масштабирование матрицы (matritsa ko'lmamining kengligi)

· HESS - приведение к форме Хессенберга (Xessenberg formasiga keltirish)

· SCHUR, RSF2CSF - приведение к форме Шура (Shura formasiga keltirish

· CPLXPAIR - сортировка комплексно сопряженных пар (kompleks bog'langanlikni saralash)

· QZ - приведение пары матриц к обобщенной форме Шура

· POLYEIG - вычисление собственных значений матричного полинома (xos sonlarning matritsa polinomligini hisoblash)

· SVD - сингулярное разложение матрицы

алгебраический программный алгоритм mathlab

4.1 Программы на языке С++

1. Поиск собственных чисел и векторов симметричной матрицы

/******************************************************************

Алгоритм ищет собственные пары симметричной матрицы, приводя её к трехдиагональной и используя QL/QR алгоритм.

Входные параметры:

A - симметричная матрица, заданная верхним или нижним треугольником. Массив с нумерацией элементов [0..N-1, 0..N-1]

N - размер матрицы

IsUpper- формат хранения

ZNeeded- флаг, сообщающий, требуются ли собственные векторы.

Если ZNeeded:

* равно 0, то собственные векторы не возвращаются

* равно 1, то собственные векторы возвращаются

Выходные параметры:

D - собственные числа в порядке возрастания.

Массив с нумерацией элементов [0..N-1].

Z - если ZNeeded:

* равно 0, не модифицируется

* равно 1, содержит собственные векторы

Массив с нумерацией элементов [0..N-1, 0..N-1], при этом

собственные векторы хранятся в столбцах матрицы.

Результат:

True, если алгоритм сошелся

False, если алгоритм не сошелся (редчайший случай)

******************************************************************/

#include <stdafx.h>

#include "sevd.h"

bool smatrixevd(ap::real_2d_array a,

int n,

int zneeded,

bool isupper,

ap::real_1d_array& d,

ap::real_2d_array& z)

{

bool result;

ap::real_1d_array tau;

ap::real_1d_array e;

ap::ap_error::make_assertion(zneeded==0||zneeded==1, "SMatrixEVD: incorrect ZNeeded");

smatrixtd(a, n, isupper, tau, d, e);

if( zneeded==1 )

{

smatrixtdunpackq(a, n, isupper, tau, z);

}

result = smatrixtdevd(d, e, n, zneeded, z);

return result;

}

/******************************************************************

Obsolete 1-based subroutine

******************************************************************/

bool symmetricevd(ap::real_2d_array a,

int n,

int zneeded,

bool isupper,

ap::real_1d_array& d,

ap::real_2d_array& z)

{

bool result;

ap::real_1d_array tau;

ap::real_1d_array e;

ap::ap_error::make_assertion(zneeded==0||zneeded==1, "SymmetricEVD: incorrect ZNeeded");

totridiagonal(a, n, isupper, tau, d, e);

if( zneeded==1 )

{

unpackqfromtridiagonal(a, n, isupper, tau, z);

}

result = tridiagonalevd(d, e, n, zneeded, z);

return result;

}

2. Разложение Холецкого

/******************************************************************

Алгоритм осуществляет разложение Холецкого симметричной положительно определенной матрицы

Результатом работы алгоритма является представление матрицы A в виде A = U'*U или A = L*L'.

Входные параметры:

A - верхний или нижний треугольник факторизуемой матрицы.

Массив с нумерацией элементов [0..N-1, 0..N-1]

N - размер матрицы

IsUpper-если IsUpper=True, A содержит верхний треугольник симметричной матрицы, иначе A содержит нижний треугольник.

Выходные параметры:

A - результат факторизации. Если IsUpper=True, в верхнем треугольнике находится матрица U, такая, что A = U'*U, а элементы, лежащие ниже главной диагонали, не модифицируются. Аналогично, если IsUpper=False.

Результат:

Если матрица положительно определена, функция возвращает True.

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

******************************************************************/

#include <stdafx.h>

#include "cholesky.h"

bool spdmatrixcholesky(ap::real_2d_array& a, int n, bool isupper)

{

bool result;

int i;

int j;

double ajj;

double v;

//

// Test the input parameters.

//

ap::ap_error::make_assertion(n>=0, "Error in SMatrixCholesky: incorrect function arguments");

// Quick return if possible

result = true;

if( n<=0 )

{

return result;

}

if( isupper )

{

// Compute the Cholesky factorization A = U'*U.

for(j = 0; j <= n-1; j++)

{

// Compute U(J,J) and test for non-positive-definiteness.

v = ap::vdotproduct(a.getcolumn(j, 0, j-1), a.getcolumn(j, 0, j-1));

ajj = a(j,j)-v;

if( ajj<=0 )

{

result = false;

return result;

}

ajj = sqrt(ajj);

a(j,j) = ajj;

// Compute elements J+1:N of row J.

if( j<n-1 )

{

for(i = 0; i <= j-1; i++)

{

v = a(i,j);

ap::vsub(&a(j, j+1), &a(i, j+1), ap::vlen(j+1,n-1), v);

}

v = 1/ajj;

ap::vmul(&a(j, j+1), ap::vlen(j+1,n-1), v);

}

}

}

else

{

// Compute the Cholesky factorization A = L*L'.

for(j = 0; j <= n-1; j++)

{

// Compute L(J,J) and test for non-positive-definiteness.

v = ap::vdotproduct(&a(j, 0), &a(j, 0), ap::vlen(0,j-1));

ajj = a(j,j)-v;

if( ajj<=0 )

{

result = false;

return result;

}

ajj = sqrt(ajj);

a(j,j) = ajj;

// Compute elements J+1:N of column J.

if( j<n-1 )

{

for(i = j+1; i <= n-1; i++)

{

v = ap::vdotproduct(&a(i, 0), &a(j, 0), ap::vlen(0,j-1));

a(i,j) = a(i,j)-v;

}

v = 1/ajj;

ap::vmul(a.getcolumn(j, j+1, n-1), v);

}

}

}

return result;

}

/******************************************************************

Obsolete 1-based subroutine.

******************************************************************/

bool choleskydecomposition(ap::real_2d_array& a, int n, bool isupper)

{

bool result;

int i;

int j;

double ajj;

double v;

int jm1;

int jp1;

//

// Test the input parameters.

//

ap::ap_error::make_assertion(n>=0, "Error in CholeskyDecomposition: incorrect function arguments");

//

// Quick return if possible

//

result = true;

if( n==0 )

{

return result;

}

if( isupper )

{

//

// Compute the Cholesky factorization A = U'*U.

//

for(j = 1; j <= n; j++)

{

//

// Compute U(J,J) and test for non-positive-definiteness.

//

jm1 = j-1;

v = ap::vdotproduct(a.getcolumn(j, 1, jm1), a.getcolumn(j, 1, jm1));

ajj = a(j,j)-v;

if( ajj<=0 )

{

result = false;

return result;

}

ajj = sqrt(ajj);

a(j,j) = ajj;

//

// Compute elements J+1:N of row J.

//

if( j<n )

{

for(i = j+1; i <= n; i++)

{

jm1 = j-1;

v = ap::vdotproduct(a.getcolumn(i, 1, jm1), a.getcolumn(j, 1, jm1));

a(j,i) = a(j,i)-v;

}

v = 1/ajj;

jp1 = j+1;

ap::vmul(&a(j, jp1), ap::vlen(jp1,n), v);

}

}

}

else

{

//

// Compute the Cholesky factorization A = L*L'.

//

for(j = 1; j <= n; j++)

{

//

// Compute L(J,J) and test for non-positive-definiteness.

//

jm1 = j-1;

v = ap::vdotproduct(&a(j, 1), &a(j, 1), ap::vlen(1,jm1));

ajj = a(j,j)-v;

if( ajj<=0 )

{

result = false;

return result;

}

ajj = sqrt(ajj);

a(j,j) = ajj;

//

// Compute elements J+1:N of column J.

//

if( j<n )

{

for(i = j+1; i <= n; i++)

{

jm1 = j-1;

v = ap::vdotproduct(&a(i, 1), &a(j, 1), ap::vlen(1,jm1));

a(i,j) = a(i,j)-v;

}

v = 1/ajj;

jp1 = j+1;

ap::vmul(a.getcolumn(j, jp1, n), v);

}

}

}

return result;

}

4.2 М - файлы для системы MatLab

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

function [L1,x1]=iter_eg(A,n,ep)

for i=1:n

x1(i)=1;

end;

L0=1;

while 1

for i=1:n

x0(i)=x1(i);

end;

for i=1:n

S=0;

for j=1:n-1

S=S+A(i,j)*x0(j);

end;

x1(i)=(S+A(i,n))/L0;

end;

S=0;

for j=1:n-1

S=S+A(n,j)*x1(j);

end;

L1=S+A(n,n);

if abs(L0-L1)<ep

break;

end;

L0=L1;

end;

end

Результат 1:

>> A=[4,2,2; 2,5,1; 2,1,6];

>> [L1,x1]=iter_eg(A,3,0.001)

L1 =

8.3884 % наибольшего собственного значения

x1 =

0.8080 0.7725 1.0000 % соответствующего ему собственного вектора

>> (A-L1*eye(3))*x1'

ans =

-0.0008

-0.0015

0

>>

Результат 2:

>> A=[2,1,4,1;3,0,1,1;-1,2,3,4;3,1,1,1]

A =

2 1 4 1

3 0 1 1

-1 2 3 4

3 1 1 1

>> [L1,x1]=iter_eg(A,4,0.001)

L1 =

6.9176

x1 =

1.3033 0.8737 1.1341 1.0000

>> (A-L1*eye(4))*x1'

ans =

0.0012

0.0000

0.0013

0

2. function [x]=holecki(A,b,n)

for i=1:n

D(i,1)=A(i,1);

C(1,i)=A(1,i)/D(1,1);

C(i,i)=1;

end;

for i=2:n

for j=2:n

if i>=j

S=0;

for k=1:j-1

S=S+D(i,k)*C(k,j);

end;

D(i,j)=A(i,j)-S;

end;

if i<j

S=0;

for k=1:i-1

S=S+D(i,k)*C(k,j);

end;

C(i,j)=(A(i,j)-S)/D(i,i);

end;

end;

end;

y(1)=b(1)/D(1,1);

for i=1:n

S=0;

for k=1:i-1

S=S+D(i,k)*y(k);

end;

y(i)=(b(i)-S)/D(i,i);

end;

x(n)=y(n);

i=n-1;

while i>=1

S=0;

for k=i+1:n

S=S+C(i,k)*x(k);

end;

x(i)=y(i)-S;

i=i-1;

end;

end

Результат:

>> A=[4,2,2; 2,5,1; 2,1,6];

>> b=[6,10,-3];

>> [x]=holecki(A,b,3)

x =

1.2250 1.7500 -1.2000

>> A*x'

ans =

6.0000

10.0000

-3.0000

ЗАКЛЮЧЕНИЕ

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

Среди итерационных методов особое место занимает метод вращений (метод Якоби) решения полной проблемы собственных значений вещественной симметричной матрицы [21]. Основан он на построении последовательности матриц, ортогонально подобных исходной матрице и имеющих монотонно убывающие до нуля суммы квадратов всех внедиагональных элементов.

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

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

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

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

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

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

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

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

СПИСОК ЛИТЕРАТУРЫ

1. Ф. Р. Гантмахер. Теория матриц. -- М.: Наука, 1967. -- 576 с.

2. Уилкинсон Д.Х. Алгебраическая проблема собственных значений. - М.: Наука, -1970. -С.564с.

3. Д.К.Фаддеев, В.Н.Фаддеева. Вычислительные методы линейной алгебры. М.: Наука, -1963.

4. В.В.Воеводин. Численные методы алгебры (теория и алгорифмы). М.: Наука, -1966.

5. Л.Коллатц. Задачи на собственные значения. М.: Наука, -1968.

6. К.Ю.Богачев. Методы решения линейных систем и нахождения собственных значений. М. -1998.

7. Белов С.А., Золотых Н.Ю. Численные методы линейное алгебры. Лабораторный практикум. - нижний Новгород: Изд-во Нижегородского госуниверситета им Н.И.Лобачевского, 2005. - 264 с.

8. Воеводин В.В., Решение полной проблемы собственных значений обобщенным методом вращений, Вычислительные методы и программирование, 3, 1965, 89 - 105.

9. Воеводин В.В., Решение полной проблемы собственных значений степенными методами, Вычислительные методы и программирование

10. Тыртышников. Мачричный анализ и линейная алгебра. М. 2005.

11. В.И.Киреев, А.В.Пантелеев. Численные методы в примерах и задачах. М.:Высшая школа, - 2006.

12. О.В.Мантуров, Ю.К.Солнцев, Ю.И.Соркин, Н.Г.Федин. Математика терминлари изо?ли лу?ати. Т.:Ў?итувчи, 1974 й.

13. Потемкин В.Г. Система MATLAB. Справочное пособие.- М., ДИАЛОГМИФИ, 1997. -370с.

14. Н.Б.Культин. С/C++ в задачах и примерах. - СПб.:БХВ - Петербург

Размещено на Allbest.ru