18.3 Разложение матриц ¶

: R = chol (A) ¶

: [R, p] = chol (A) ¶

: [R, p, Q] = chol (A) ¶

: [R, p, Q] = chol (A, "vector") ¶

: [L, …] = chol (…, "lower") ¶

: [R, …] = chol (…, "upper") ¶

Вычислить верхний фактор Холеского, R, вещественной симметричной или комплексной эрмитовой положительно определенной матрицы A.

Верхний фактор Холеского R вычисляется с использованием верхней треугольной части матрицы A и определяется следующим образом:

R' * R = A.

Вызов chol с необязательным флагом "upper" имеет тот же результат. В отличие от этого, используя необязательный флаг "lower", chol возвращает нижнее треугольное разложение, вычисленное с использованием нижней треугольной части матрицы A так, что

L * L' = A.

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

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

Если вызывается с тремя выходными аргументами, матрица A должна быть разреженной, и к матрице A применяется сохраняющая разреженность перестановка строк/столбцов перед факторизацией. То есть R является разложением A(Q,Q) так, что

R' * R = Q' * A * Q.

Перестановка, сохраняющая разреженность, обычно возвращается как матрица. Однако, с учётом необязательного флага "vector", Q будет возвращена как вектор, так что

R' * R = A(Q, Q).

В общем случае нижнее треугольное разложение значительно быстрее для разреженных матриц.

См. также: hess, lu, qr, qz, schur, svd, ichol, cholinv, chol2inv, cholupdate, cholinsert, choldelete, cholshift.

: Ainv = cholinv (A) ¶

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

См. также: chol, chol2inv, inv.

: Ainv = chol2inv (R) ¶

Найти обратную симметричной, положительно определенной квадратной матрицы по её разложению Холеского, R.

Обратите внимание, что R должна быть верхнетреугольной матрицей с положительными диагональными элементами. chol2inv (U) предоставляет inv (R'*R), но значительно быстрее, чем использование inv.

См. также: chol, cholinv, inv.

: [R1, info] = cholupdate (R, u, op) ¶

Обновить или пересчитать разложение Холеского.

Дана верхняя треугольная матрица R и столбец u, попытаться определить другую верхнюю треугольную матрицу R1 такую, что

• R1’*R1 = R’*R + u*u’ если op равно "+"
• R1’*R1 = R’*R - u*u’ если op равно "-"

Если op равно "-", info устанавливается в

• 0, если пересчёт был успешным,
• 1, если R’*R - u*u’ не является положительно определенной,
• 2, если R вырожденная.

Если info отсутствует, в случаях 1 и 2 выводится сообщение об ошибке.

См. также: chol, cholinsert, choldelete, cholshift.

: R1 = cholinsert (R, j, u) ¶

: [R1, info] = cholinsert (R, j, u) ¶

Обновление разложения Холеского с учётом вставки строки или столбца в исходную факторизованную матрицу.

Дано разложение Холеского вещественной симметричной или комплексной эрмитовой положительно определённой матрицы A = R’*R, R верхняя треугольная, вернуть разложение Холеского матрицы A1, где A1(p,p) = A, A1(:,j) = A1(j,:)’ = u и p = [1:j-1,j+1:n+1]. u(j) должно быть положительным.

При возврате info устанавливается в

• 0, если вставка успешна,
• 1, если A1 не является положительно определённой,
• 2, если R вырожденная.

Если info отсутствует, в случаях 1 и 2 выводится сообщение об ошибке.

См. также: chol, cholupdate, choldelete, cholshift.

: R1 = choldelete (R, j) ¶

Обновить разложение Холеского, удалив строку или столбец из исходной факторизованной матрицы.

Дано разложение Холеского вещественной симметричной или комплексной эрмитовой положительно определенной матрицы A = R’*R, R верхняя треугольная, вернуть разложение Холеского A(p,p), где p = [1:j-1,j+1:n+1].

См. также: chol, cholupdate, cholinsert, cholshift.

: R1 = cholshift (R, i, j) ¶

Обновить факторизацию Холецкого, задав диапазон столбцов для сдвига в исходной факторизованной матрице.

Учитывая факторизацию Холецкого вещественной симметричной или комплексной эрмитовой положительно определенной матрицы A = R’*R, R верхняя треугольная, вернуть факторизацию Холецкого для A(p,p), где p — перестановка
p = [1:i-1, shift(i:j, 1), j+1:n] если i < j
или
p = [1:j-1, shift(j:i,-1), i+1:n] если j < i.

См. также: chol, cholupdate, cholinsert, choldelete.

: H = hess (A) ¶

: [P, H] = hess (A) ¶

Вычислить декомпозицию Хессенберга матрицы A.

Декомпозиция Хессенберга имеет вид P * H * P' = A где P — квадратная унитарная матрица (P' * P = I, используя комплексно-сопряжённое транспонирование) и H — верхняя матрица Хессенберга (H(i, j) = 0 forall i > j+1).

Декомпозиция Хессенберга обычно используется как первый шаг в вычислении собственных значений, но также имеет и другие применения (см. Golub, Nash, и Van Loan, IEEE Transactions on Automatic Control, 1979).

См. также: eig, chol, lu, qr, qz, schur, svd.

: [L, U] = lu (A) ¶

: [L, U, P] = lu (A) ¶

: [L, U, P, Q] = lu (S) ¶

: [L, U, P, Q, R] = lu (S) ¶

: […] = lu (S, thresh) ¶

: y = lu (…) ¶

: […] = lu (…, "vector") ¶

Вычислить LU-разложение A.

Если A полная, то используются подпрограммы из LAPACK, а если A разреженная, то используется UMFPACK.

Результат возвращается в переупорядоченном виде, в соответствии с необязательным возвращаемым значением P. Например, для матрицы A = [1, 2; 3, 4],

[L, U, P] = lu (A)

возвращается

L =

  1.00000  0.00000
  0.33333  1.00000

U =

  3.00000  4.00000
  0.00000  0.66667

P =

  0  1
  1  0

Матрица не обязана быть квадратной.

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

При вызове с четвёртым выходным аргументом возвращается сохранённая по столбцам преобразование Q такое, что P * A * Q = L * U. Это предпочтительный способ вызова lu с разреженными входными матрицами.

При вызове с пятым выходным аргументом и разреженной входной матрицей, lu пытается использовать множитель масштабирования R для входной матрицы, так что P * (R \ A) * Q = L * U. Это обычно приводит к более разреженному и стабильному разложению.

Можно указать дополнительный входной аргумент thresh, определяющий порог поворота. thresh может быть скаляром, в этом случае он определяет порог поворота UMFPACK для симметричных и несимметричных случаев. Если thresh — вектор из 2 элементов, то первый элемент определяет порог поворота для несимметричной стратегии поворота UMFPACK, а второй — для симметричной. По умолчанию используются значения, определённые spparms ([0.1, 0.001]).

При передаче строкового аргумента "vector", lu возвращает значения P и Q как векторные значения, такие что для полной матрицы A(P,:) = L * U, и R(P,:) * A(:,Q) = L * U.

С двумя выходными аргументами возвращаются переупорядоченные формы верхней и нижней треугольных матриц, такие что A = L * U. С одним выходным аргументом y возвращается матрица, возвращаемая подпрограммами LAPACK. Если входная матрица разреженная, то матрица L встраивается в U, чтобы получить результат, похожий на полный случай. Для полных и разреженных матриц lu теряет информацию о перестановке.

См. также: luupdate, ilu, chol, hess, qr, qz, schur, svd.

: [L, U] = luupdate (L, U, x, y) ¶

: [L, U, P] = luupdate (L, U, P, x, y) ¶

Учитывая LU-разложение вещественной или комплексной матрицы A = L*U, L нижняя единичная трапециевидная и U верхняя трапециевидная, вернуть LU-разложение A + x*y.’, где x и y — столбцовые векторы (обновление ранга 1) или матрицы с одинаковым количеством столбцов (обновление ранга k).

Необязательно можно использовать обновление с поворотом по строкам, передав матрицу перестановки строк (поворота) P; в этом случае возвращается обновлённая матрица перестановки. Обратите внимание, если L, U, P — LU-разложение с поворотом, полученное с помощью lu:

[L, U, P] = lu (A);

то факторизацию A+x*y.' можно получить либо как

[L1, U1] = lu (L, U, P*x, y)

или

[L1, U1, P1] = lu (L, U, P, x, y)

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

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

См. также: lu, cholupdate, qrupdate.

: [Q, R] = qr (A) ¶

: [Q, R, P] = qr (A) ¶

: X = qr (A) # non-sparse A ¶

: R = qr (A) # sparse A ¶

: X = qr (A, B) # sparse A ¶

: [C, R] = qr (A, B) ¶

: […] = qr (…, 0) ¶

: […] = qr (…, "econ") ¶

: […] = qr (…, "vector") ¶

: […] = qr (…, "matrix") ¶

Вычислить QR-разложение A, используя стандартные подпрограммы LAPACK.

QR-разложение имеет вид

Q * R = A

где Q — ортогональная матрица, а R — верхняя треугольная матрица.

Например, для матрицы A = [1, 2; 3, 4],

[Q, R] = qr (A)

возвращает

Q =

  -0.31623  -0.94868
  -0.94868   0.31623

R =

  -3.16228  -4.42719
   0.00000  -0.63246

произведение которых возвращает исходную матрицу

Q * R
  ⇒
     1.0000   2.0000
     3.0000   4.0000

Если запрашивается только одно значение возврата, то это либо R, если A разреженная, или X, такое что R = triu (X) , если A полная. (Обратите внимание: в отличие от большинства команд, одно возвращаемое значение не является первым, когда запрашивается несколько значений.)

Если запрашивается третий выходной параметр P, то qr вычисляет переупорядоченное QR-разложение

Q * R = A * P

где Q — ортогональная матрица, R — верхняя треугольная матрица, а P — матрица перестановки.

Если A плотная, переупорядоченное QR-разложение обладает дополнительным свойством, что диагональные элементы R упорядочены по убыванию модуля. Другими словами, abs (diag (R)) будут упорядочены от наибольшего к наименьшему.

Если A разреженная, P — упорядочение столбцов A с целью уменьшения заполнения. В этом случае диагональные элементы R не упорядочены по убыванию модуля.

Например, для матрицы A = [1, 2; 3, 4],

[Q, R, P] = qr (A)

возвращает

Q =

  -0.44721  -0.89443
  -0.89443   0.44721

R =

  -4.47214  -3.13050
   0.00000   0.44721

P =

   0  1
   1  0

Если входная матрица A разреженная, разреженное QR-разложение вычисляется с помощью SPQR или CXSPARSE (например, если SPQR недоступен). Поскольку матрица Q в общем случае является полной матрицей, рекомендуется запрашивать только одно возвращаемое значение R. В этом случае вычисления избегают построения Q и возвращают разреженную R такую, что R = chol (A' * A).

Если A плотная, а дополнительная матрица B предоставлена и запрашиваются два возвращаемых значения, то qr возвращает C, где C = Q' * B. Это позволяет вычислить приближение наименьших квадратов A \ B как

[C, R] = qr (A, B)
X = R \ C

Если A — разреженная матрица MxN и предоставлена дополнительная матрица B, возможно одно или два возвращаемых значения. Если запрашивается одно возвращаемое значение X и M < N, то X — решение с минимальной 2-нормой A \ B. Если M ≥ N, X — приближение наименьших квадратов A \ B. Если запрашиваются два возвращаемых значения, C и R имеют тот же смысл, что и в плотном случае (C — плотная, а R — разреженная). Версия с одним параметром возврата предпочтительнее, так как она использует меньше памяти и лучше обрабатывает матрицы с вырожденным рангом.

Если последним аргументом является строка "vector", то P — вектор перестановок (столбцов A) вместо матрицы перестановок. В этом случае определяющее соотношение:

Q * R = A(:, P)

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

Если последним аргументом является скаляр 0 или строка "econ", возвращается экономичное разложение. Если исходная матрица A имеет размер MxN и M > N, то экономичное разложение вычислит только N строк в R и N столбцов в Q и опустит нули в R. Если M ≤ N, нет разницы между экономичным и стандартным разложениями. При вычислении экономичного разложения и A плотной, выходной параметр P всегда является вектором, а не матрицей. Если A разреженная, выходной параметр P — разреженная матрица перестановок.

Примечание: QR-разложение применяется при решении задач наименьших квадратов

min norm (A*x - b)

для переопределённых систем уравнений (т. е. A — высокая, узкая матрица).

Переупорядоченное QR-разложение [Q, R, P] = qr (A) позволяет построить ортогональный базис span (A).

См. также: chol, hess, lu, qz, schur, svd, qrupdate, qrinsert, qrdelete, qrshift.

: [Q1, R1] = qrupdate (Q, R, u, v) ¶

Обновление QR-разложения с помощью векторов или матриц обновлений.

Дано QR-разложение вещественной или комплексной матрицы A = Q*R, где Q унитарна, а R верхняя трапециевидная, возвращает QR-разложение A + u*v’, где u и v — столбцовые векторы (обновление ранга 1) или матрицы с одинаковым количеством столбцов (обновление ранга k). Обратите внимание, что последний случай выполняется как последовательность обновлений ранга 1; таким образом, при достаточно большом k, повторное вычисление разложения будет быстрее и точнее.

Предоставленное QR-разложение может быть полным (Q — квадратная матрица) или экономичным (R — квадратная матрица).

См. также: qr, qrinsert, qrdelete, qrshift.

: [Q1, R1] = qrinsert (Q, R, j, x, orient) ¶

Обновление QR-разложения с помощью строки или столбца для вставки в исходную факторизованную матрицу.

Дано QR-разложение вещественной или комплексной матрицы A = Q*R, где Q унитарна, а R верхняя трапециевидная, возвращает QR-разложение [A(:,1:j-1) x A(:,j:n)], где u — столбцовый вектор для вставки в A (если orient — "col"), или QR-разложение [A(1:j-1,:);x;A(:,j:n)], где x — строковый вектор для вставки в A (если orient — "row").

Значение по умолчанию для orient — "col". Если orient — "col", u может быть матрицей, а j — вектором индексов, что приводит к QR-разложению матрицы B такой, что B(:,j) даёт u и B(:,j) = [] даёт A. Обратите внимание, что последний случай выполняется как последовательность k вставок; таким образом, при достаточно большом k, повторное вычисление разложения будет быстрее и точнее.

Если orient — "col", предоставленное QR-разложение может быть полным (Q — квадратная матрица) или экономичным (R — квадратная матрица).

Если orient — "row", необходимо полное разложение.

См. также: qr, qrupdate, qrdelete, qrshift.

: [Q1, R1] = qrdelete (Q, R, j, orient) ¶

Обновить QR-разложение, удалив строку или столбец из исходной факторизованной матрицы.

Дано QR-разложение вещественной или комплексной матрицы A = Q*R, где Q унитарна, а R верхняя трапециевидная, вернуть QR-разложение [A(:,1:j-1), U, A(:,j:n)], где u — столбец, который нужно вставить в A (если orient равно "col"), или QR-разложение [A(1:j-1,:);X;A(:,j:n)], где x — строка, и orient равно "row"). Значение по умолчанию для orient равно "col".

Если orient равно "col", j может быть вектором индексов, что приводит к QR-разложению матрицы B такой, что A(:,j) = [] даёт B. Обратите внимание, что последний случай выполняется как последовательность из k удалений; таким образом, при достаточно большом k будет быстрее и точнее повторно вычислить разложение с нуля.

Если orient равно "col", предоставленное QR-разложение может быть полным (Q квадратная) или экономичным (R квадратная).

Если orient равно "row", необходимо полное разложение.

См. также: qr, qrupdate, qrinsert, qrshift.

: [Q1, R1] = qrshift (Q, R, i, j) ¶

Обновить QR-разложение, сдвинув диапазон столбцов в исходной факторизованной матрице.

Дано QR-разложение вещественной или комплексной матрицы A = Q*R, где Q унитарна, а R верхняя трапециевидная, вернуть QR-разложение A(:,p), где p — перестановка
p = [1:i-1, shift(i:j, 1), j+1:n] если i < j
или
p = [1:j-1, shift(j:i,-1), i+1:n] если j < i.

См. также: qr, qrupdate, qrinsert, qrdelete.

: [AA, BB, Q, Z, V, W] = qz (A, B) ¶

: [AA, BB, Q, Z, V, W] = qz (A, B, opt) ¶

Вычислить QZ-разложение задачи обобщённых собственных значений.

Задача обобщённых собственных значений определяется как

A x = lambda B x

Функция имеет две формы вызова:

1. [AA, BB, Q, Z, V, W, lambda] = qz (A, B)

   Вычислить комплексное QZ-разложение, обобщённые собственные векторы и обобщённые собственные значения.

   AA = Q * A * Z, BB = Q * B * Z
   A * V * diag (diag (BB)) = B * V * diag (diag (AA))
   diag (diag (BB)) * W' * A = diag (diag (AA)) * W' * B

   с AA и BB верхними треугольными, а Q и Z унитарными. Матрицы V и W соответственно содержат правые и левые обобщённые собственные векторы.

2. [AA, BB, Z {, lambda}] = qz (A, B, opt)

   Аргумент opt должен быть равен либо "real", либо "complex". Если он равен "complex", то эта форма вызова эквивалентна первой с двумя входными аргументами.

   Если opt равен "real", то вычисляется вещественное QZ-разложение. В частности, AA гарантированно является квази-верхней треугольной с 1×1 и 2×2 блоками на диагонали, а Q и Z ортогональны. Указанные выше тождества для правых и левых обобщённых собственных векторов верны только если AA является верхней треугольной (т.е. когда все обобщённые собственные значения вещественные, в этом случае вещественное и комплексное QZ совпадают).

Примечание: qz выполняет балансировку перестановки, но не масштабирование (см. balance), что может привести к менее точным результатам, чем eig. Порядок выходных аргументов выбран для совместимости с MATLAB.

См. также: eig, gsvd, balance, chol, hess, lu, qr, qzhess, schur.

: [aa, bb, q, z] = qzhess (A, B) ¶

Вычислить разложение Хессенберга-треугольной матрицы (A, B), возвращая aa = q * A * z, bb = q * B * z, с q и z ортогональными.

Например:

[aa, bb, q, z] = qzhess ([1, 2; 3, 4], [5, 6; 7, 8])
  ⇒ aa =
      -3.02244  -4.41741
       0.92998   0.69749
  ⇒ bb =
      -8.60233  -9.99730
       0.00000  -0.23250
  ⇒ q =
      -0.58124  -0.81373
      -0.81373   0.58124
  ⇒ z =
     Diagonal Matrix
       1   0
       0   1

Разложение Хессенберга-треугольной матрицы — это первый шаг в алгоритме QZ-разложения Молера и Стюарта.

Алгоритм взят из работы Голуба и Ван Лоана, Матричные вычисления, 2-е издание.

См. также: lu, chol, hess, qr, qz, schur, svd.

: S = schur (A) ¶

: S = schur (A, "real") ¶

: S = schur (A, "complex") ¶

: S = schur (A, opt) ¶

: [U, S] = schur (…) ¶

Вычислить разложение Шура матрицы A.

Разложение Шура квадратной матрицы A определяется как

S = U' * A * U

где U — унитарная матрица (U'* U — единичная матрица), а S — верхняя треугольная. Собственные значения матрицы A (и S) — диагональные элементы S. Если матрица A вещественная, вычисляется вещественное разложение Шура, в котором матрица U ортогональна, а S — блочно-верхняя треугольная с блоками размера не более 2 x 2 на диагонали.

По умолчанию для вещественных матриц используется вещественное разложение Шура. Комплексное разложение можно принудительно вычислить, передав флаг "complex".

Собственные значения можно упорядочить на диагонали в зависимости от значения opt:

opt = "a"

Переместить собственные значения с отрицательной вещественной частью в ведущий блок S. Мемоническая фраза: "a" для алгебраических уравнений Риккати, где это упорядочение полезно.

opt = "d"

Переместить собственные значения с модулем меньше единицы в ведущий блок S. Мемоническая фраза: "d" для дискретных алгебраических уравнений Риккати, где это упорядочение полезно.

opt = "u"

Без упорядочения. Нет особого упорядочения собственных значений (по умолчанию).

Ведущие k столбцов U всегда охватывают подпространство, инвариантное относительно A, соответствующее k ведущим собственным значениям S.

См. также: rsf2csf, ordschur, ordeig, lu, chol, hess, qr, qz, svd, eig.

: [U, T] = rsf2csf (UR, TR) ¶

Преобразовать вещественную верхнюю квазитреугольную форму Шура TR в комплексную верхнюю треугольную форму Шура T.

Обратите внимание, что выполняются следующие соотношения:

UR * TR * UR' = U * T * U' и U' * U — единичная матрица I.

Также обратите внимание, что U и T не уникальны.

См. также: schur.

: [UR, SR] = ordschur (U, S, select) ¶

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

Логический вектор select указывает выбранные собственные значения, как они появляются вдоль диагонали S.

Например, задана матрица A = [1, 2; 3, 4], и её разложение Шур

[U, S] = schur (A)

что возвращает

U =

  -0.82456  -0.56577
   0.56577  -0.82456

S =

  -0.37228  -1.00000
   0.00000   5.37228

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

[U, S] = ordschur (U, S, [0,1])

См. также: schur, ordeig, ordqz.

: [AR, BR, QR, ZR] = ordqz (AA, BB, Q, Z, keyword) ¶

: [AR, BR, QR, ZR] = ordqz (AA, BB, Q, Z, select) ¶

Переупорядочить QZ-разложение задачи обобщенных собственных значений.

Задача обобщенных собственных значений определяется как

A x = lambda B x

Её обобщенное разложение Шур вычисляется с помощью алгоритма qz:

[AA, BB, Q, Z] = qz (A, B)

где AA, BB, Q, и Z удовлетворяют

AA = Q * A * Z, BB = Q * B * Z

Функция ordqz вычисляет унитарное преобразование QR и ZR таким образом, чтобы порядок собственного значения на диагонали AA и BB был изменен. Результирующие переупорядоченные матрицы AR и BR удовлетворяют:

AR = QR * A * ZR, BR = QR * B * ZR

Функцию можно вызвать с аргументом keyword, который выбирает собственные значения в верхнем левом блоке AR и BR следующим образом:

"S", "udi"

small: ведущий блок содержит все |lambda| < 1

"B", "udo"

big: ведущий блок содержит все |lambda| ≥ 1

"-", "lhp"

отрицательная вещественная часть: ведущий блок содержит все собственные значения в открытой левой полуплоскости

"+", "rhp"

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

Если вместо ключевого слова задан логический вектор select, функция ordqz переупорядочивает все собственные значения k в левый блок, для которого select(k) истинно.

Примечание: Ключевые слова совместимы с теми, что в qr.

См. также: eig, ordeig, qz, schur, ordschur.

: lambda = ordeig (A) ¶

: lambda = ordeig (A, B) ¶

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

Квазитреугольная матрица A обычно является результатом факторизации Шур. Если вызвана со вторым входом B, то обобщенные собственные значения пары A, B возвращаются в порядке появления матрицы A-lambda*B. Пара A, B обычно является результатом QZ-разложения.

См. также: ordschur, ordqz, eig, schur, qz.

: angle = subspace (A, B) ¶

Определяет наибольший главный угол между двумя подпространствами, охватываемыми столбцами матриц A и B.

: s = svd (A) ¶

: [U, S, V] = svd (A) ¶

: [U, S, V] = svd (A, "econ") ¶

: [U, S, V] = svd (A, 0) ¶

Вычислить сингулярное разложение A.

Сингулярное разложение определяется соотношением

A = U*S*V'

Функция svd обычно возвращает только вектор сингулярных значений. Когда вызывается с тремя возвращаемыми значениями, она вычисляет U, S и V. Например,

svd (hilb (3))

возвращает

ans =

  1.4083189
  0.1223271
  0.0026873

и

[u, s, v] = svd (hilb (3))

возвращает

u =

  -0.82704   0.54745   0.12766
  -0.45986  -0.52829  -0.71375
  -0.32330  -0.64901   0.68867

s =

  1.40832  0.00000  0.00000
  0.00000  0.12233  0.00000
  0.00000  0.00000  0.00269

v =

  -0.82704   0.54745   0.12766
  -0.45986  -0.52829  -0.71375
  -0.32330  -0.64901   0.68867

Когда задан второй аргумент, не равный 0, svd возвращает разложение экономичного размера, исключая ненужные строки или столбцы U или V.

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

Примечания к алгоритму: При вычислении полного разложения (левые и правые сингулярные матрицы в дополнение к сингулярным значениям) существует выбор двух процедур в LAPACK. По умолчанию Octave использует процедуру gesvd. Альтернативная процедура gesdd в 5 раз быстрее, но может использовать больше памяти и быть неточной для некоторых входных матриц. Существует третья процедура gejsv, подходящая для лучшей точности при экстремальных масштабах. См. документацию по svd_driver для получения дополнительной информации о выборе драйвера.

См. также: svd_driver, svds, eig, lu, chol, hess, qr, qz.

: val = svd_driver () ¶

: old_val = svd_driver (new_val) ¶

: old_val = svd_driver (new_val, "local") ¶

Запрос или установка базового драйвера LAPACK, используемого svd.

В настоящее время распознаются значения "gesdd", "gesvd", и "gejsv". По умолчанию используется "gesvd".

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

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

Новая функция gesdd основана на алгоритме «Разделяй и властвуй» и в 5 раз быстрее альтернативной функции gesvd, основанной на QR-разложении. Однако новый алгоритм может использовать значительно больше памяти. Для матрицы MxN размер памяти составляет O(min(M,N)^2), в то время как для альтернативы — O(max(M,N)).

Функция gejsv использует алгоритм предварительно обусловленного разложения Якоби. В отличие от gesvd и gesdd, в gejsv отсутствует этап бидиагонализации, который в некоторых крайних случаях мог бы ухудшить точность. Кроме того, gejsv известна как оптимально точная в некотором смысле. Однако скорость медленнее (основной механизм однопоточный) и используется больше памяти (O(min(M,N)^2 + M + N)).

Помимо проблем со скоростью и памятью, были случаи, когда некоторые входные матрицы не были точно разложены функцией gesdd. См. активный баг https://savannah.gnu.org/bugs/?55564. До тех пор, пока проблемы точности не будут решены в новой версии библиотеки LAPACK, по умолчанию в Octave используется драйвер "gesvd".

См. также: svd.

: [housv, beta, zer] = housh (x, j, z) ¶

Вычисление вектора отражения Хаусхолдера housv для отражения x до j-ой колонки единичной матрицы, т.е.

(I - beta*housv*housv')x =  norm (x)*e(j) if x(j) < 0,
(I - beta*housv*housv')x = -norm (x)*e(j) if x(j) >= 0

Входы

x

вектор

j

индекс в векторе

z

порог для нуля (обычно должно быть число 0)

Выходы (см. Голуб и Ван Лоан):

beta

Если beta = 0, то отражение применять не нужно (устанавливается 0)

housv

вектор Хаусхолдера

: [u, h, nu] = krylov (A, V, k, eps1, pflg) ¶

Построение ортогонального базиса u блочного подпространства Кристова.

Блочное подпространство Кристова имеет следующий вид:

[v a*v a^2*v ... a^(k+1)*v]

Построение выполняется с помощью отражений Хаусхолдера для предотвращения потери ортогональности.

Если V — вектор, то h содержит матрицу Хессенберга такую, что a*u == u*h+rk*ek', в которой rk = a*u(:,k)-u*h(:,k), и ek' является вектором [0, 0, …, 1] длины k. В противном случае h бессмысленно.

Если V — вектор и k больше length (A) - 1, то h содержит матрицу Хессенберга такую, что a*u == u*h.

Значение nu — это размерность охвата подпространства Кристова (основанного на eps1).

Если b — вектор и k больше m-1, то h содержит разложение Хессенберга матрицы A.

Необязательный параметр eps1 — порог для нуля. Значение по умолчанию равно 1e-12.

Если необязательный параметр pflg отличен от нуля, используется строчный поворот для повышения числовой устойчивости. Значение по умолчанию равно 0.

Ссылка: А. Ходель, П. Мисра, Частичный поворот при вычислении подпространств Кристова больших разреженных систем, Труды 42-й конференции IEEE по принятию решений и управлению, декабрь 2003 г.