Введение в описание комплекса программ
для сингулярного разложения матриц

Пусть  A - вещественная матрица общего вида размера  m на n. Сингулярным разложением ( SVD ) матрицы  A  называется разложение вида

                A  =  U S VT ,

где           U и V - ортогональные матрицы порядка  m и n  соответственно, а
                      S  =  diag (s1, ..., sr )

                 r  =  min (m, n) ,   s1 ≥ ... ≥ sr ≥ 0.

Если  A - комплексная матрица, ее сингулярное разложение имеет вид

               A  =  U S VH ,

где         U и V - унитарные матрицы порядка  m и n  соответственно, а  S, также как и выше,
                           имеет вещественные диагональные элементы.

s_i называются сингулярными числами, первые  r  столбцов матрицы V - правыми сингулярными векторами, а первые  r  столбцов матрицы U - левыми сингулярными векторами матрицы A.

Сингулярные числа и сингулярные векторы удовлетворяют условиям

            A vi  =  si ui      и      AT ui  =  si vi   (или   AH ui  =  si vi ), 

где  u_i  и  v_i   являются   i - ми столбцами матриц  U  и  V  соответственно.

Подробнее, вычисление состоит из следующих этапов.

1. Матрица  A  приводится к двухдиагональному виду A = U₁ B V₁^T, если A  является вещественной (и A = U₁ B V₁^H - если A - комплексная), где U₁ и V₁ являются ортогональными матрицами (или унитарными, если  A - комплексная), а  B  является вещественной и верхней двухдиагональной при  m ≥ n (или нижней двухдиагональной при  m < n ).

2. Выполняется сингулярное разложение ( SVD ) двухдиагональной матрицы  B:

         B  =  U2 S V2T ,

где     U2 и V2 - ортогональные матрицы, а S - диагональная (как указывалось выше). 

Сингулярные числа матрицы A совпадают с сингулярными числами матрицы B. Сингулярными векторами матрицы  A  являются  U = U₁ U₂  и  V = V₁ V₂.

Приведение вещественной матрицы общего вида к двухдиагональной форме выполняется базовой подпрограммой  PDGEBRD, а сингулярное разложение матрицы  B  выполняется с использованием базовой подпрограммы  DBDSQR из пакета LAPACK [11,12].
Если требуется вычислить только сингулярные числа, то делается обращение к подпрограмме DLASQ1, которая выполняет это вычисление, используя  qd- алгоритм (см.[36]).

Базовая подпрограмма  PDGEBRD представляет  U₁ и  V₁ в виде произведений элементарных матриц отражения.
Для вещественной матрицы  A  матрицы  U₁ и  V₁ могут быть умножены на другие матрицы без явного формирования U₁  и  V₁  посредством использования базовой подпрограммы PDORMBR.

Решение проблемы сингулярного разложения для вещественных квадратных матриц выполняется с помощью целевых подпрограмм комплекса PDGESVD1, PDGESVD2, PDGESVD3.

Если  m >> n, может быть более эффективно сначала выполнить QR-разложение матрицы  A, используя базовую подпрограмму  PDGEQRF, а затем выполнить сингулярное разложение квадратной матрицы  R порядка  n, т.к. если  A = Q R  и R = U S V^T, то сингулярное разложение  A  задается формулой A = ( Q U ) S V^T. Аналогично, если  m << n, может быть более эффективно сначала выполнить  LQ факторизацию  A, используя базовую подпрограмму  DGELQF. Эти предварительные факторизации выполняются самими целевыми подпрограммами для прямоугольных матриц PDGESVD4, PDGESVD5, PDGESVD6.

Сингулярное разложение может быть использовано для нахождения решения с минимальной нормой для линейной задачи наименьших квадратов (возможно) неполного ранга. Действительный ранг  k  матрицы  A  может быть определен как количество сингулярных чисел, превосходящих некоторое заданное малое число. Пусть ^~S является ведущей подматрицей матрицы S порядка  k, а ^~V является матрицей, состоящей из первых  k  столбцов матрицы  V. Тогда решение задается формулой

                  x  =  ~V ~S - 1 ~c1 ,

где    ~c1 - состоит из первых  k  элементов

                 c  =  UT b  =  U2T U1T b . 

U₁^T b может быть вычислено с помощью подпрограммы  PDORMBR.

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