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 не завершается ошибкой. Значение 0 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.

cholinv (A)

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

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

chol2inv (U)

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

Обратите внимание, что U должна быть верхней треугольной матрицей с положительными диагональными элементами. chol2inv (U) предоставляет inv (U'*U), но это намного быстрее, чем использовать 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, and 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, thres)

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. Это обычно приводит к более разреженному и более устойчивому разложению.

Дополнительный входной аргумент thres, определяющий порог для перестановки, может быть задан. thres может быть скаляром, в этом случае он определяет порог перестановки UMFPACK для симметричных и несимметричных случаев. Если thres является вектором из 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) # non-sparse A

X = qr (A) # non-sparse A

R = qr (A) # sparse A

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

[…] = qr (…, 0)

[…] = 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 — полная матрица. (Примечание: в отличие от большинства команд, единственный возвращаемый результат не является первым возвращаемым значением, когда запрашивается несколько значений.)

Если матрица A полная, и требуется третий выход P, то qr вычисляет переставленную QR-факторизацию

Q * R = A * P

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

Переставленная QR-факторизация обладает дополнительным свойством, что диагональные элементы R упорядочены по убыванию абсолютных значений. Другими словами, abs (diag (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-факторизация вычисляется с помощью CSPARSE. Поскольку матрица Q вообще является полной матрицей, рекомендуется запросить только одно возвращаемое значение R. В этом случае вычисление избегает построения Q и возвращает R такое, что R = chol (A' * A).

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

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

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

Q * R = A(:, P)

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

Если последним аргументом является скаляр 0, то возвращается экономная QR-факторизация. Когда исходная матрица A имеет размер MxN и M > N, то экономная QR-факторизация вычислит только N строк в R и N столбцов в Q и опустит нули в R. Если M ≤ N, то нет разницы между экономной и стандартной факторизациями. При вычислении экономной факторизации выход 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 — строка, которую нужно вставить в A (если 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.

lambda = qz (A, B)

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

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

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

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

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

A x = lambda B x

Существует три варианта вызова функции:

1. lambda = qz (A, B)

   Вычислить обобщённые собственные значения lambda.

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

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

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

   с Q и Z ортогональными (унитарными для комплексного случая).

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

   Как и в форме 2 выше, но позволяет упорядочить обобщённые пары собственных значений, например, для решения дискретных алгебраических уравнений Риккати. Форма 3 недоступна для комплексных матриц и не вычисляет обобщённые собственные векторы V, W, а также ортогональную матрицу Q.

   opt

   для упорядочивания собственных значений обобщённой задачи на собственные значения. Главный блок пересмотренного определителя содержит все собственные значения, которые удовлетворяют:

   "N"

   неупорядоченные (по умолчанию)

   "S"

   малые: главный блок содержит все |lambda| < 1

   "B"

   большие: главный блок содержит все |lambda| ≥ 1

   "-"

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

   "+"

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

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

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

[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 Молера и Стюарта.

Алгоритм взят из Golub и Van Loan, Вычисления матриц, 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.

Разложение Шура определяется как

S = U' * A * U

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

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

Собственные значения, по желанию, упорядочиваются вдоль диагонали в соответствии со значением opt. opt = "a" указывает, что все собственные значения с отрицательной вещественной частью должны быть перемещены в ведущий блок S (используется в are), opt = "d" указывает, что все собственные значения с модулем меньше единицы должны быть перемещены в ведущий блок S (используется в dare), а opt = "u", значение по умолчанию, указывает, что никакого упорядочения собственных значений не должно происходить. Ведущие k столбцов U всегда охватывают A-инвариантное подпространство, соответствующее k ведущим собственным значениям S.

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

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

[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.

lambda = ordeig (A)

lambda = ordeig (A, B)

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

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

См. также: ordschur, 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 раз быстрее, но может использовать больше памяти и может быть неточной для некоторых входных матриц. См. документацию для svd_driver для получения дополнительной информации о выборе драйвера.

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

val = svd_driver ()

old_val = svd_driver (new_val)

svd_driver (new_val, "local")

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

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

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

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

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