Используйте chol разложить на множители симметричную матрицу коэффициентов, и затем решить линейную систему с помощью Фактора Холецкого.

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

A = [1 0 1; 0 2 0; 1 0 3]

A = 3×3

     1     0     1
     0     2     0
     1     0     3

Вычислите Фактор Холецкого матрицы.

R = chol(A)

R = 3×3

    1.0000         0    1.0000
         0    1.4142         0
         0         0    1.4142

Создайте вектор для правой стороны уравнения $\mathrm{Ax}=\mathit{b}$.

b = sum(A,2);

С тех пор $\mathit{A}={\mathit{R}}^{\mathit{T}}\mathit{R}$ с разложением Холесского линейное уравнение становится ${\mathit{R}}^{\mathit{T}}\mathit{R}\mathit{x}=\mathit{b}$. Решите для x использование оператора обратной косой черты.

x = R\(R'\b)

x = 3×1

    1.0000
    1.0000
    1.0000

Вычислите верхние и более низкие факторизации Холецкого матрицы и проверьте результаты.

Создайте 6 6 симметричную положительную определенную тестовую матрицу использование gallery функция.

A = gallery('lehmer',6);

Вычислите Фактор Холецкого с помощью верхнего треугольника A.

R = chol(A)

R = 6×6

    1.0000    0.5000    0.3333    0.2500    0.2000    0.1667
         0    0.8660    0.5774    0.4330    0.3464    0.2887
         0         0    0.7454    0.5590    0.4472    0.3727
         0         0         0    0.6614    0.5292    0.4410
         0         0         0         0    0.6000    0.5000
         0         0         0         0         0    0.5528

Проверьте, что верхний треугольный множитель удовлетворяет R'*R - A = 0, в ошибке округления.

norm(R'*R - A)

ans = 2.5801e-16

Теперь задайте 'lower' опция, чтобы вычислить Фактор Холецкого с помощью более низкого треугольника A.

L = chol(A,'lower')

L = 6×6

    1.0000         0         0         0         0         0
    0.5000    0.8660         0         0         0         0
    0.3333    0.5774    0.7454         0         0         0
    0.2500    0.4330    0.5590    0.6614         0         0
    0.2000    0.3464    0.4472    0.5292    0.6000         0
    0.1667    0.2887    0.3727    0.4410    0.5000    0.5528

Проверьте, что нижний треугольный фактор удовлетворяет L*L' - A = 0, в ошибке округления.

norm(L*L' - A)

ans = 2.5801e-16

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

Создайте матрицу 5 на 5 биномиальных коэффициентов. Эта матрица симметрична положительный определенный, поэтому вычтите 1 из последнего элемента, чтобы гарантировать, что это более не положительно определенный.

A = pascal(5);
A(end) = A(end) - 1

A = 5×5

     1     1     1     1     1
     1     2     3     4     5
     1     3     6    10    15
     1     4    10    20    35
     1     5    15    35    69

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

[R,flag] = chol(A)

R = 4×4

     1     1     1     1
     0     1     2     3
     0     0     1     3
     0     0     0     1

flag = 5

Начиная с flag является ненулевым, это дает индекс центра, где факторизация перестала работать. chol может вычислить q = flag-1 = 4 строки и столбцы правильно прежде, чем перестать работать, когда это сталкивается с частью матрицы, которая изменилась.

Проверьте тот R'*R возвращает четыре строки и столбца, которые соглашаются с A(1:q,1:q).

q = flag-1;
R'*R

ans = 4×4

     1     1     1     1
     1     2     3     4
     1     3     6    10
     1     4    10    20

A(1:q,1:q)

ans = 4×4

     1     1     1     1
     1     2     3     4
     1     3     6    10
     1     4    10    20

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

Создайте разреженную положительную определенную матрицу на основе west0479 матрица.

load west0479
A = west0479;
S = A'*A;

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

[R,flag] = chol(S);
[RP,flagP,P] = chol(S);

Для каждого вычисления проверяйте тот flag = 0 подтвердить вычисление успешно.

if ~flag && ~flagP
    disp('Factorizations successful.')
else
    disp('Factorizations failed.')
end

Factorizations successful.

Сравните количество ненулей в chol(S) по сравнению с переупорядоченным матричным chol(P'*S*P). Лучшая практика должна использовать три выходных синтаксиса chol с разреженными матрицами, начиная с переупорядочения строк и столбцов может значительно сократить количество ненулей в Факторе Холецкого.

subplot(1,2,1)
spy(R)
title('Nonzeros in chol(S)')
subplot(1,2,2)
spy(RP)
title('Nonzeros in chol(P''*S*P)')

Используйте 'vector' опция chol возвратить информацию о сочетании как вектор, а не матрицу.

Создайте разреженную матрицу конечного элемента.

S = gallery('wathen',10,10);
spy(S)

Вычислите Фактор Холецкого для матрицы и задайте 'vector' опция, чтобы возвратить вектор сочетания p.

[R,flag,p] = chol(S,'vector');

Проверьте тот flag = 0, указание на вычисление успешно.

if ~flag
    disp('Factorization successful.')
else
    disp('Factorization failed.')
end

Factorization successful.

Проверьте тот S(p,p) = R'*R, в ошибке округления.

norm(S(p,p) - R'*R,'fro')

ans = 2.1039e-13