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

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

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

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

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

Что можно сделать с помощью любой из программ гл. 1.

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

метод, составлены в соответствии с блок-схемой рис. 3.2. Блоки с номером О, расположены в основной программе, все другие оформлены в виде подпрограмм.

Рис. 3.2. (см. скан) Блок-схема программы интерполяционного метода вычисления собственных значений матрицы

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

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

В подпрограмме блока 1 (строки 100-190) с помощью двойного цикла формируется матрица В путем последовательного ввода ее элементов по строкам.

В подпрограмме блока 2 (строки 200-290) вычисляются узлы и значения функции интерполяционной таблицы. Для уменьшения погрешности вычисления коэффициентов полинома Ньютона рекомендуется [1] размещать узлы интерполяции равномерно между границами поиска собственных значений. Поэтому расстояние между узлами задается переменной Узлы с нулевым номером выбираем в точке В строках 210-280 организован цикл по переменной которой соответствует номер узла С помощью двойного цикла в строках 220-260 по переменным формируется характеристическая матрица к. Для исключения перекрытия переменных для величины к введена переменная После обращения к подпрограмме вычисления определителя (строка 270) результате заносится в массив текущее значение переменной в массив X. Перед завершением цикла по переменной получаем очередное значение переменной для следующего узла путем добавления шага

Блок 3 представляет собой подпрограмму вычисления резделенных разностей, являющихся коэффициентами полинома Ньютона Настоящая подпрограмма (строки 300-390) написана в соответствии с аналогичным блоком программы с учетом исключения перекрытия переменных.

В качестве конкретного метода решения характеристического уравнения (блок 4) выбран метод секущих, подпрограмма которого (строки 400-490) записана на основании программы Для исключения перекрытия переменных начальные приближения соответствуют переменным и XI, для текущего значения аргумента левой части уравнения выбрана переменная Y, значения характеристического полинома присваиваются переменной

Интерполяционный полином Ньютона вычисляется по схеме Горнера (строки 500-510). Найденные собственные значения исключаются путем деления полинома на произведение разностей между собственными значениями и текущим значением аргумента (цикл в строке 520). В результате осуществляется переход от уравнения

к уравнению

где найденные собственные значения; - номер искомого собственного значения.

Для вычисления определителей используется метод Гаусса, подпрограмма которого записана в строках 600-790 из программы с учетом исключения перекрытия переменных.

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

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

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

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

Для контроля программ 3.5 в вещественной области можно использовать пример Для тестирования программы определим собственные значения матрицы [40]

Границы области поиска задаем (-2, 2), (5.1, -7.1), погрешности -

(см. скан)

(см. скан)

(см. скан)

(см. скан)

(см. скан)

(см. скан)