Вычислительная сложность разреженных операций пропорциональна nnz
, количество ненулевых элементов в матрице. Вычислительная сложность также зависит линейно от размера строки m
и размер столбца n
из матрицы, но независимо от продукта m*n
, общее количество нулевых и ненулевых элементов.
Сложность довольно сложных операций, таких как решение разреженных линейных уравнений, включает факторы как упорядоченное расположение и временная замена, которые обсуждены в предыдущем разделе. В общем случае однако машинное время, требуемое для операции разреженной матрицы, пропорционально количеству арифметических операций на ненулевых количествах.
Разреженные матрицы распространяют посредством расчетов согласно этим правилам:
Функции, которые принимают матрицу и возвращают скаляр или вектор постоянного размера всегда, производят выход в полном формате устройства хранения данных. Например, size
функционируйте всегда возвращает полный вектор, полон ли его вход или разрежен.
Функции, которые принимают скаляры или векторы и возвращают матрицы, такой как zeros
, ones
, rand
, и eye
, всегда возвращайте полные результаты. Это необходимо, чтобы не вводить разреженность неожиданно. Разреженный аналог zeros(m,n)
просто sparse(m,n)
. Разреженные аналоги rand
и eye
sprand
и speye
, соответственно. Нет никакого разреженного аналога для функции ones
.
Унарные функции, которые принимают матрицу и возвращают матрицу или вектор, сохраняют класс памяти операнда. Если S
разреженная матрица, затем chol(S)
также разреженная матрица и diag(S)
разреженный вектор. Постолбцовые функции такой как max
и sum
также возвратите разреженные векторы, даже при том, что эти векторы могут быть совершенно ненулевыми. Важные исключения к этому правилу sparse
и full
функции.
Бинарные операторы приводят к разреженным результатам, если оба операнда являются разреженными, и полными результатами, если оба полны. Для смешанных операндов результат полон, если операция не сохраняет разреженность. Если S
разреженно и F
полно, затем S+F
, S*F
, и F\S
полны, в то время как S.*F
и S&F
разреженны. В некоторых случаях результат может быть разреженным даже при том, что матрица имеет немного нулевых элементов.
Конкатенация матриц с помощью любого cat
функциональные или квадратные скобки приводят к разреженным результатам для смешанных операндов.
Сочетание строк и столбцов разреженной матрицы S
может быть представлен двумя способами:
Матрица перестановок P
действия на строках S
как P*S
или на столбцах как S*P'
.
Вектор сочетания p
, который является полным вектором, содержащим сочетание 1:n
, действия на строках S
как S(p,:)
, или на столбцах как S(:,p)
.
Например:
p = [1 3 4 2 5] I = eye(5,5); P = I(p,:) e = ones(4,1); S = diag(11:11:55) + diag(e,1) + diag(e,-1)
p = 1 3 4 2 5 P = 1 0 0 0 0 0 0 1 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 0 1 S = 11 1 0 0 0 1 22 1 0 0 0 1 33 1 0 0 0 1 44 1 0 0 0 1 55
Можно теперь попробовать некоторые сочетания с помощью вектора сочетания p
и матрица перестановок P
. Например, операторы S(p,:)
и P*S
возвратите ту же матрицу.
S(p,:)
ans = 11 1 0 0 0 0 1 33 1 0 0 0 1 44 1 1 22 1 0 0 0 0 0 1 55
P*S
ans = 11 1 0 0 0 0 1 33 1 0 0 0 1 44 1 1 22 1 0 0 0 0 0 1 55
Точно так же S(:,p)
и S*P'
оба производят
S*P'
ans = 11 0 0 1 0 1 1 0 22 0 0 33 1 1 0 0 1 44 0 1 0 0 1 0 55
Если P
разреженная матрица, затем оба представления используют устройство хранения данных, пропорциональное n
и можно применить любого к S
вовремя пропорциональный nnz(S)
. Векторное представление немного более компактно и эффективно, таким образом, различные стандартные программы сочетания разреженной матрицы, все возвращают полные векторы-строки за исключением вертящегося сочетания в LU (треугольная) факторизация, которая возвращает матрицу, совместимую с полной LU-факторизацией.
Преобразовывать между двумя представлениями сочетания:
n = 5; I = speye(n); Pr = I(p,:); Pc = I(:,p); pc = (1:n)*Pc
pc = 1 3 4 2 5
pr = (Pr*(1:n)')'
pr = 1 3 4 2 5
Инверсия P
просто R = P'
. Можно вычислить инверсию p
с r(p) = 1:n
.
r(p) = 1:5
r = 1 4 2 3 5
Переупорядочение столбцов матрицы может часто делать свой LU или факторы QR более разреженными. Переупорядочение строк и столбцов может часто делать свои Факторы Холецкого более разреженными. Самое простое такое переупорядочение должно отсортировать столбцы по ненулевому количеству. Это иногда - хорошее переупорядочение для матриц с очень неправильными структурами, особенно если существует большое изменение ненулевых количеств строк или столбцов.
colperm
вычисляет сочетание, которое заказывает столбцы матрицы количеством ненулей в каждом столбце от самого маленького до самого большого.
Упорядоченное расположение обратного алгоритма Катхилла-Макки предназначается, чтобы уменьшать профиль или полосу пропускания матрицы. Это, как гарантируют, не найдет самую маленькую полосу пропускания, но это обычно делает. symrcm
функция на самом деле управляет на ненулевой структуре симметрической матрицы A + A'
, но результат также полезен для несимметричных матриц. Это упорядоченное расположение полезно для матриц, которые прибывают из одномерных проблем или проблем, которые находятся в некотором смысле, длинном и тонком.
Степень узла в графике является количеством связей с тем узлом. Это совпадает с количеством недиагональных ненулевых элементов в соответствующей строке матрицы смежности. Аппроксимированный минимальный алгоритм степени генерирует упорядоченное расположение на основе того, как эти степени изменены во время Исключения Гаусса или факторизации Холесского. Это - сложный и мощный алгоритм, который обычно приводит к более разреженным факторам, чем большинство других упорядоченных расположений, включая количество столбцов и обратный алгоритм Катхилла-Макки. Поскольку отслеживание степени каждого узла очень длительно, аппроксимированный минимальный алгоритм степени использует приближение для степени, а не точной степени.
Они MATLAB® функции реализуют аппроксимированный минимальный алгоритм степени:
Смотрите Переупорядочение и Факторизацию Разреженных матриц для примера с помощью symamd
.
Можно изменить различные параметры, сопоставленные с деталями алгоритмов с помощью spparms
функция.
Для получения дополнительной информации на алгоритмах используется colamd
и symamd
, см. [5]. Аппроксимированная степень использование алгоритмов основана [1].
Как аппроксимированное минимальное упорядоченное расположение степени, вложенный алгоритм упорядоченного расположения рассечения, реализованный dissect
функционируйте переупорядочивает строки матрицы и столбцы, полагая, что матрица матрица смежности графика. Алгоритм уменьшает проблему вниз до намного меньшей шкалы путем сворачивания вместе пар вершин в графике. После переупорядочения маленького графика алгоритм затем применяет проекцию, и улучшение продвигается, чтобы расширить график назад до первоначального размера.
Вложенный алгоритм рассечения производит высококачественные переупорядочения и выполняет особенно хорошо с матрицами конечного элемента по сравнению с другими методами переупорядочения. Для получения дополнительной информации о вложенном алгоритме упорядоченного расположения рассечения, см. [7].
Если S
разреженная матрица, следующая команда возвращает три разреженных матрицы L
U
, и P
таким образом, что P*S = L*U
.
[L,U,P] = lu(S);
lu
получает факторы Исключением Гаусса с частичным поворотом. Матрица перестановок P
имеет только n
ненулевые элементы. Как с плотными матрицами, оператором [L,U] = lu(S)
возвращает переставленный модуль нижняя треугольная матрица и верхняя треугольная матрица, продуктом которой является S
. Отдельно, lu(S)
возвращает L
и U
в одной матрице без информации о центре.
Синтаксис с тремя выходами [L,U,P] = lu(S)
выбирает P
через числовой частичный поворот, но не вертится, чтобы улучшить разреженность в LU
факторы. С другой стороны, синтаксис с четырьмя выходами [L,U,P,Q] = lu(S)
выбирает P
через порог частичный поворот, и выбирает P
и Q
улучшить разреженность в LU
факторы.
Можно управлять поворотом в использовании разреженных матриц
lu(S,thresh)
где thresh
порог центра в [0,1]. Поворот происходит, когда диагональный элемент в столбце имеет величину меньше, чем thresh
времена величина любого поддиагонального элемента в том столбце. thresh = 0
поворот диагонали сил. thresh = 1
значение по умолчанию. (Значение по умолчанию для thresh
0.1
для синтаксиса с четырьмя выходами).
Когда вы вызываете lu
с тремя или меньше выходными параметрами MATLAB автоматически выделяет память, необходимую, чтобы содержать разреженный L
и U
факторы во время факторизации. За исключением синтаксиса с четырьмя выходами, MATLAB не использует символьной предварительной факторизации LU, чтобы определить требования к памяти и настроить структуры данных заранее.
Переупорядочение и факторизация разреженных матриц
Этот пример показывает эффекты переупорядочения и факторизации на разреженных матрицах.
Если вы получаете хорошее сочетание столбца p
это уменьшает временную замену, возможно, от symrcm
или colamd
, затем вычисление lu(S(:,p))
занимает меньше времени и устройства хранения данных, чем вычисление lu(S)
.
Создайте разреженную матрицу с помощью Маркерного примера мяча.
B = bucky;
B
имеет точно три ненулевых элемента в каждой строке и столбце.
Создайте два сочетания, r
и m
использование symrcm
и symamd
соответственно.
r = symrcm(B); m = symamd(B);
Эти два сочетания являются симметричным упорядоченным расположением обратного алгоритма Катхилла-Макки и симметричным аппроксимированным минимальным упорядоченным расположением степени.
Создайте графики шпиона показать три матрицы смежности Маркерного графика Мяча с этими тремя различными нумерациями. Локальная, основанная на пятиугольнике структура исходной нумерации не очевидна в других.
figure subplot(1,3,1) spy(B) title('Original') subplot(1,3,2) spy(B(r,r)) title('Reverse Cuthill-McKee') subplot(1,3,3) spy(B(m,m)) title('Min Degree')
Упорядоченное расположение обратного алгоритма Катхилла-Макки, r
, уменьшает полосу пропускания и концентрирует все ненулевые элементы около диагонали. Аппроксимированное минимальное упорядоченное расположение степени, m
, производит подобную фракталу структуру с большими блоками нулей.
Чтобы видеть временную замену, сгенерированную в LU-факторизации Маркерного мяча, используйте speye
, разреженная единичная матрица, чтобы вставить-3s на диагонали B
.
B = B - 3*speye(size(B));
Поскольку каждая сумма строки является теперь нулем, этим новым B
на самом деле сингулярно, но это все еще поучительно, чтобы вычислить его LU-факторизацию. Когда названо только одним выходным аргументом, lu
возвращает эти два треугольных множителя, L
и U
, в одной разреженной матрице. Количество ненулей в той матрице является мерой времени и устройства хранения данных, требуемого решить линейные системы, включающие B
.
Вот ненулевые счета для этих трех рассматриваемых сочетаний.
lu(B)
(Исходный): 1022
lu(B(r,r))
(Обратный алгоритм Катхилла-Макки): 968
lu(B(m,m))
(Аппроксимируйте минимальную степень): 636
Даже при том, что это - небольшой пример, результаты типичны. Исходная схема нумерации приводит к большей части временной замены. Временная замена для упорядоченного расположения обратного алгоритма Катхилла-Макки сконцентрирована в полосе, но это почти столь же обширно как первые два упорядоченных расположения. Для аппроксимированного минимального упорядоченного расположения степени относительно большие блоки нулей сохраняются во время устранения, и сумма временной замены значительно меньше что сгенерирована другими упорядоченными расположениями.
spy
графики ниже отражения характеристики каждого переупорядочения.
figure subplot(1,3,1) spy(lu(B)) title('Original') subplot(1,3,2) spy(lu(B(r,r))) title('Reverse Cuthill-McKee') subplot(1,3,3) spy(lu(B(m,m))) title('Min Degree')
Если S
симметричное (или Эрмитов), положительная определенная, разреженная матрица, оператор ниже возвратов разреженная, верхняя треугольная матрица R
так, чтобы R'*R = S
.
R = chol(S)
chol
автоматически не вертится для разреженности, но можно вычислить аппроксимированную минимальную степень и ограничивающие сочетания профиля для использования с chol(S(p,p))
.
Поскольку алгоритм Холесского не использует поворот для разреженности и не требует поворота для числовой устойчивости, chol
делает быстрое вычисление требуемого объема памяти и выделяет всю память в начале факторизации. Можно использовать symbfact
, который использует тот же алгоритм в качестве chol
, вычислять, сколько памяти выделяется.
MATLAB вычисляет полную QR-факторизацию разреженной матрицы S
с
[Q,R] = qr(S)
[Q,R,E] = qr(S)
но это часто непрактично. Унитарная матрица Q
часто сбои, чтобы иметь высокий процент нулевых элементов. Более практическая альтернатива, иногда известная как “QR-факторизацию Q-less”, доступна.
С одним разреженным входным параметром и одним выходным аргументом
R = qr(S)
возвращает только верхний треугольный фрагмент QR-факторизации. Матричный R
обеспечивает факторизацию Холесского для матрицы, сопоставленной с нормальными уравнениями:
R'*R = S'*S
Однако потеря числовой информации, свойственной от расчета S'*S
избегается.
С двумя входными параметрами, имеющими одинаковое число строк, и два выходных аргумента, оператор
[C,R] = qr(S,B)
применяет ортогональные преобразования к B
, создание C = Q'*B
не вычисляя Q
.
QR-факторизация Q-less позволяет решение разреженных проблем наименьших квадратов
с двумя шагами:
[c,R] = qr(A,b); x = R\c
Если A
является разреженным, но не квадратным, MATLAB использует эти шаги для линейного уравнения, решая оператор обратной косой черты:
x = A\b
Или, можно сделать факторизацию сами и исследовать R
для дефицита ранга.
Также возможно решить последовательность линейных систем наименьших квадратов с различными правыми сторонами, b
, это не обязательно известно когда R = qr(A)
вычисляется. Подход решает “полунормальные уравнения R'*R*x = A'*b
с
x = R\(R'\(A'*b))
и затем использует один шаг итеративного улучшения, чтобы уменьшать, округляют ошибку:
r = b - A*x; e = R\(R'\(A'*r)); x = x + e
ilu
и ichol
функции обеспечивают аппроксимированные, неполные факторизации, которые полезны как предварительные формирователи для разреженных итерационных методов.
ilu
функция производит три incomplete lower-upper (ILU) факторизации: zero-fill ILU (ILU(0)
), версия Crout ILU (ILUC(tau)
), и ILU с пороговым отбрасыванием и поворотом (ILUTP(tau)
). ILU(0)
никогда центры и получившиеся факторы только не имеют ненули в положениях, где входная матрица имела ненули. Оба ILUC(tau)
и ILUTP(tau)
, однако, сделайте пороговое отбрасывание с пользовательским допуском отбрасывания tau
.
Например:
A = gallery('neumann', 1600) + speye(1600);
ans = 7840
nnz(lu(A))
ans = 126478
показывает тот A
имеет 7 840 ненулей, и его полная LU-факторизация имеет 126 478 ненулей. С другой стороны, следующий код показывает различному ILU выходные параметры:
[L,U] = ilu(A); nnz(L)+nnz(U)-size(A,1)
ans = 7840
norm(A-(L*U).*spones(A),'fro')./norm(A,'fro')
ans = 4.8874e-17
opts.type = 'ilutp';
opts.droptol = 1e-4;
[L,U,P] = ilu(A, opts);
nnz(L)+nnz(U)-size(A,1)
ans = 31147
norm(P*A - L*U,'fro')./norm(A,'fro')
ans = 9.9224e-05
opts.type = 'crout';
[L,U,P] = ilu(A, opts);
nnz(L)+nnz(U)-size(A,1)
ans = 31083
norm(P*A-L*U,'fro')./norm(A,'fro')
ans = 9.7344e-05
Эти вычисления показывают, что заполнить нулями факторы имеют 7 840 ненулей, ILUTP(1e-4)
факторы имеют 31 147 ненулей и ILUC(1e-4)
факторы имеют 31 083 ненуля. Кроме того, относительная погрешность продукта заполнить нулями факторов является по существу нулем на шаблоне A
. Наконец, относительная погрешность факторизаций, произведенных с пороговым отбрасыванием, находится того же порядка допуска отбрасывания, несмотря на то, что это, как гарантируют, не произойдет. Смотрите ilu
страница с описанием для большего количества опций и деталей.
ichol
функция обеспечивает zero-fill incomplete Cholesky factorizations (IC(0)
) а также пороговые понижающиеся неполные факторизации Холесского (ICT(tau)
) из симметричных, положительных определенных разреженных матриц. Эти факторизации являются аналогами неполных LU-факторизаций выше и имеют многие из тех же характеристик. Например:
A = delsq(numgrid('S',200));
nnz(A)
ans = 195228
nnz(chol(A,'lower'))
ans = 7762589
A
имеет 195 228 ненулей, и его полная факторизация Холесского без переупорядочения имеет 7 762 589 ненулей. В отличие от этого:L = ichol(A); nnz(L)
ans = 117216
norm(A-(L*L').*spones(A),'fro')./norm(A,'fro')
ans = 3.5805e-17
opts.type = 'ict';
opts.droptol = 1e-4;
L = ichol(A,opts);
nnz(L)
ans = 1166754
norm(A-L*L','fro')./norm(A,'fro')
ans = 2.3997e-04
IC(0)
имеет ненули только в шаблоне более низкого треугольника A
, и на шаблоне A
, продукт соответствий факторов. Кроме того, ICT(1e-4)
факторы значительно более разреженны, чем полный Фактор Холецкого и относительная погрешность между A
and L*L'
находится того же порядка допуска отбрасывания. Важно отметить что различающийся факторы, обеспеченные chol
, факторы по умолчанию, обеспеченные ichol
являются нижними треугольными. Смотрите ichol
страница с описанием для получения дополнительной информации.
Две функции доступны, которые вычисляют несколько заданных собственных значений или сингулярных значений. svds
основан на eigs
.
Функции, чтобы вычислить несколько собственных значений или сингулярных значений
Эти функции наиболее часто используются с разреженными матрицами, но они могут использоваться с полными матрицами или даже с линейными операторами, заданными в коде MATLAB.
Оператор
[V,lambda] = eigs(A,k,sigma)
находит k
собственные значения и соответствующие собственные вектора матричного A
это является самым близким “сдвиг” sigma
. Если sigma
не использован, собственные значения, самые большие в величине, найдены. Если sigma
нуль, собственные значения, самые маленькие в величине, найдены. Вторая матрица, B
, может быть включен для обобщенной задачи о собственных значениях: Aυ = λBυ.
Оператор
[U,S,V] = svds(A,k)
находит k
самые большие сингулярные значения A
и
[U,S,V] = svds(A,k,'smallest')
находит k
самые маленькие сингулярные значения.
Числовые методы используются в eigs
и svds
описаны в [6].
В этом примере показано, как найти самое маленькое собственное значение и собственный вектор разреженной матрицы.
Настройте Дифференциальный оператор Лапласа с пятью точками на 65 65 сетка в L-образной, двумерной области.
L = numgrid('L',65);
A = delsq(L);
Определите порядок и количество ненулевых элементов.
size(A)
ans = 1×2
2945 2945
nnz(A)
ans = 14473
A
матрица порядка 2945 с 14 473 ненулевыми элементами.
Вычислите самое маленькое собственное значение и собственный вектор.
[v,d] = eigs(A,1,'smallestabs');
Распределите компоненты собственного вектора по соответствующим узлам решетки и произведите контурный график результата.
L(L>0) = full(v(L(L>0)));
x = -1:1/32:1;
contour(x,x,L,15)
axis square
[1] Amestoy, P. R., Т. А. Дэвис и я. S. Подновите, “Аппроксимированный Минимальный Алгоритм Упорядоченного расположения Степени”, SIAM Journal на Анализе матрицы и Приложения, Издание 17, № 4, октябрь 1996, стр 886-905.
[2] Барретт, R., М. Берри, Т. Ф. Чан, и др., Шаблоны для Решения Линейных систем: Базовые блоки для Итерационных методов, SIAM, Филадельфия, 1994.
[3] Дэвис, T.A., Гильберт, J. R. Larimore, S.I., Ын, E., Пейтон, B., “Столбец Аппроксимированный Минимальный Алгоритм Упорядоченного расположения Степени”, Proc. Конференция SIAM по Прикладной Линейной алгебре, октябрь 1997, p. 29.
[4] Гильберт, Джон Р., Клив Молер и Роберт Шрайбер, “Разреженные матрицы в MATLAB: Разработка и реализация”, SIAM J. Matrix Anal. Appl., Издание 13, № 1, январь 1992, стр 333-356.
[5] Larimore, S. I. аппроксимированный минимальный столбец степени, заказывая алгоритм, тезис MS, отдел информатики и вычислительной техники и разработки, Флоридского университета, Гейнсвилл, FL, 1998.
[6] Саад, Yousef, итерационные методы для разреженных линейных уравнений. PWS Publishing Company, 1996.
[7] Karypis, Джордж и Випин Кумар. "Быстрая и Высококачественная Многоуровневая Схема Разделения Неправильных Графиков". SIAM Journal на Научных вычислениях. Издание 20, Номер 1, 1999, стр 359–392.