Файл: М.А. Тынкевич Система Matlab Справочное пособие к курсу Численные методы анализа.pdf

20

[Q,R,P]=qr(A) отличается от предыдущего упорядочением по убыванию модулей диагональных элементов R и наличием соответствующей матрицы перестановок Р (A*P’=Q*R);

[Q,R]=qr(A,0) при m>n отличается тем, что вычисляются лишь n столбцов матрицы Q.

Если после выполнения QR-разложения выполнить команду [Q1,R1]=qrdelete(Q,R,k), то будет выполнен пересчет матриц для варианта, когда в матрице А удален k-й столбец. Если после QRразложения выполнить команду [Q1,R1]=qrdelete(Q,R,k,X), то будут пересчитаны матрицы для варианта, когда в матрице А перед столбцом k вставлен столбец Х.

X=nnis(A,B), X=nnis(A,B,t) позволяют искать решение систе-

мы АХ=В методом наименьших квадратов, где отыскиваются не-

отрицательные решения Х, минимизирующие norm(A*X-B) или га-

рантирующие точность ε при задании t=max(m,n)*norm(A,1)* ε. Особое место в библиотеке занимают средства для вычисления

собственных чисел и векторов.

В простейшем варианте отыскиваются ненулевые решения системы АХ=λX командами d=eig(A) или [X,d]=eig(A) (d- диагональная матрица собственных чисел, Х –матрица из нормированных собственных векторов):

Для проверки качества поиска (при значительных размерностях и многочисленных особых случаях такая проверка весьма желательна) достаточно проверить на близость к нулю значений A*X=R*D:

» a*R-R*d

ans = 1.0e-014 * 0.0583 0.0666 -0.7549 0.0166 -0.0444 -0.5329 0.0902 -0.3886 -0.7105

Функции d=eig(A,B) и [V,D]= eig(A,B) позволяют решать полную (обобщенную) проблему собственных значений АХ=λBX.

Решение задачи осуществляется на основе QR-алгоритма и его модификаций и при числе итераций, превышающем 30 n, может быть прервано с сообщением Solution will not converge (решение не сходится).

Проблему собственных значений можно решать и для матрич-

21

ного полинома (A0+λA1+λ2A2+...+λpAp)X=0 командой [R,d]=polyeig(A0,A1,.., Ap), где R- матрица размера n×(n×p) собственных векторов. При р=0 эта функция тождественна eig(A0) , при р=1 – eig(A0,-A1) и в случае матриц с n=1 (скаляров) – roots(Ap,...,A1,A0), то есть ищет корни уравнения Apλp+...+ A2λ2 +A1λ +A0=0.

Для решения ряда задач используется функция сингулярного разложения матрицы в формах s=svd(A), [U,S,V]=svd(A), [U,S,V]=svd(A,0) - матрица А размерности m×n (m≥n) представляется в виде A=U*S*V′, где U′*U=V*V′=E, S=diag(s1,s2,...,sn). Здесь U со-

стоит из n собственных векторов для n наибольших собственных значений матрицы АА′, а V – из ортонормированных собственных векторов матрицы А′А; на диагонали матрицы S - значения квадратных корней из собственных чисел матрицы А′А (сингулярные числа).

6. Операции над полиномами

Рассмотрим полином вида Pn(x)=p1xn+p2xn-1+...+pnx+an+1 . Соответственно будем обозначать Р – n+1-мерный вектор коэффициентов,

Х – массив значений аргумента.

При вычислении значений полинома для элементов массива можно использовать функцию polyval(P,X):

j =max( 1,k +1−n )

Деление полиномов можно реализовать командой [C,R]=deconv(A,B) , где С –частное и R – остаток от деления А на В.

22

Вычисление производных от полинома, произведения и от-

ношения полиномов производится соответственно командами dp=polyder(P), dc=polyder(A,B) и [f,g]=polyder(A,B):

Вычисление корней полинома реализуется функцией roots(P), а построение полинома по его корням – функцией poly(R).

В приложениях, особенно связанных с преобразованием Лапласа при решении дифференциальных уравнений, оперируют с отношениями полиномов и представлениями их в виде простых дробей:

где sk - простые корни полинома Qn(s); если некоторый корень sj имеет кратность m, то соответствующее слагаемое представляется в виде

Команда [r,s,f] =residue(P,Q) дает разложение отношения по-

линомов на простые дроби (в случае близких корней возможна значительная погрешность). В случае кратного корня пользуются функцией rj=resi2(P,Q,sj,m,j), где j –номер вычисляемого коэффициента (по умолчанию j=m); по умолчанию m=1 (простой корень).

Команда [P,Q] =residue(r,s,f) выполняет обратное действие свертки разложения в отношение полиномов.

Выполнив действия

23

видим, что

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

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

hadamard (n) – матрица Адамара (ортогональная матрица из 1 и –1); n должно быть целым и n, n/12 или n/20 должны быть степенью 2:

» hadamard (4)

hilb (n), invhilb(n) – матрица Гильберта hij=1/(i+j-1) и ей обратная (пример матрицы, плохо обусловленной :к обращению):

magic(n) – магический квадрат (квадратная матрица с элементами от 1 до n2 с равными суммами элементов по строкам и столбцам):

24

pascal(n) – матрица Паскаля - симметрическая матрица из коэффициентов разложения бинома (1+х)j (треугольника Паскаля)

1

rosser - матрица Рессера (матрица 8-го порядка, служащая тестом для алгоритмов решения симметричной проблемы собственных значений ):

имеет два кратных значения, три близких, нулевое и малое ненулевое):

1000 1000 1020.049 1020.000 1019.902 0.098 0 -1020.049 toepliz(Х) – симметрическая матрица Теплица, определяющая

перестановки элементов вектора Х ;

toepliz(Х,Y) – несимметрическая матрица Теплица, первый столбец которой совпадает с вектором X и первая строка с вектором Y (если х1≠у1, возникает конфликт на главной диагонали с предпочтением для X):