ichol

Неполная факторизация Холесского

Синтаксис

L = ichol(A)
L = ichol(A,opts)

Описание

L = ichol(A) выполняет неполную факторизацию Холесского A с заполняют нулями.

L = ichol(A,opts) выполняет неполную факторизацию Холесского A с опциями, заданными opts.

По умолчанию, ichol ссылается на более низкий треугольник A и производит нижние треугольные факторы.

Входные параметры

A

Разреженная матрица

opts

Структура максимум с пятью полями:

Имя поляСводные данныеОписание
typeТип факторизации

Указывает который разновидность неполного Холесского выполнить. Допустимыми значениями этого поля является 'nofill' и 'ict'. 'nofill' вариант выполняет incomplete Cholesky with zero-fill (IC(0)). 'ict' вариант performs incomplete Cholesky with threshold dropping (ICT). Значением по умолчанию является 'nofill'.

droptol Пропустите допуск, когда типом будет 'ict'

Неотрицательный скаляр, используемый в качестве допуска отбрасывания при выполнении ICT. Элементы, которые меньше в величине, чем локальный допуск отбрасывания, исключены из получившегося фактора за исключением диагонального элемента, который никогда не пропускается. Локальный допуск отбрасывания на шаге j из факторизации norm(A(j:end,j),1)*droptol. 'droptol' проигнорирован если 'type' 'nofill'Значение по умолчанию 0.

micholУказывает, выполнить ли модифицированный неполный Холесский

Указывает, выполняется ли modified incomplete Cholesky (MIC). Полем может быть 'on' или 'off'. При выполнении MIC диагональ компенсируется пропущенные элементы, чтобы осуществить отношение A*e = L*L'*e где   e = ones(size(A,2),1). Значением по умолчанию является 'off'.

diagcompВыполните компенсированный неполный Холесский с заданным коэффициентом

Действительный неотрицательный скаляр, используемый в качестве глобальной диагонали, переключает alpha в формировании неполного Фактора Холесского. Таким образом, вместо того, чтобы выполнить неполный Холесский на A, факторизация A + alpha*diag(diag(A)) формируется. Значение по умолчанию 0.

shapeОпределяет, на какой треугольник ссылаются и возвращают

Допустимыми значениями является 'upper' и 'lower'. Если 'upper' задан, только верхний треугольник A ссылается и R создается таким образом что A аппроксимирован R'*R. Если 'lower' задан, только более низкий треугольник A ссылается и L создается таким образом что A аппроксимирован L*L'. Значением по умолчанию является 'lower'.

Примеры

свернуть все

Этот пример генерирует неполную факторизацию Холесского.

Начните с симметричной положительной определенной матрицы, A:

N = 100;
A = delsq(numgrid('S',N));

A двумерный, дискретный отрицательный Лапласиан с пятью точками на 100 100 квадратной сетке с граничными условиями Дирихле. Размер A 98*98 = 9604 (не 10000, когда границы сетки используются, чтобы наложить условия Дирихле).

Неполная факторизация Холесского без заливок является факторизацией, которая содержит только ненули в той же позиции A содержит ненули. Эта факторизация является чрезвычайно дешевой, чтобы вычислить. Несмотря на то, что продукт L*L' обычно очень отличается от A, продукт L*L' будет совпадать с A на его шаблоне до округления.

L = ichol(A);
norm(A-L*L','fro')./norm(A,'fro')
ans = 0.0916
norm(A-(L*L').*spones(A),'fro')./norm(A,'fro')
ans = 4.9606e-17

ichol май также использоваться, чтобы сгенерировать неполные факторизации Холесского с пороговым отбрасыванием. Когда допуск отбрасывания уменьшается, фактор имеет тенденцию становиться более плотным и продукт L*L' имеет тенденцию быть лучшим приближением A. Следующие графики показывают относительную погрешность неполной факторизации, построенной против допуска отбрасывания, а также отношения плотности неполных факторов к плотности полного Фактора Холесского.

n = size(A,1);
ntols = 20;
droptol = logspace(-8,0,ntols);
nrm = zeros(1,ntols);
nz = zeros(1,ntols);
nzComplete = nnz(chol(A,'lower'));
for k = 1:ntols
    L = ichol(A,struct('type','ict','droptol',droptol(k)));
    nz(k) = nnz(L);
    nrm(k) = norm(A-L*L','fro')./norm(A,'fro');
end
figure
loglog(droptol,nrm,'LineWidth',2)
title('Drop tolerance vs norm(A-L*L'',''fro'')./norm(A,''fro'')')

figure
semilogx(droptol,nz./nzComplete,'LineWidth',2)
title('Drop tolerance vs fill ratio ichol/chol')

Относительная погрешность обычно находится на том же порядке как допуск отбрасывания, несмотря на то, что это не гарантируется.

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

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

N = 100;
A = delsq(numgrid('S',N));

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

b = ones(size(A,1),1);
tol = 1e-6;
maxit = 100;
[x0,fl0,rr0,it0,rv0] = pcg(A,b,tol,maxit);

Обратите внимание на то, что fl0 = 1 указание, что pcg не управлял относительной невязкой к требуемому допуску в максимальных позволенных итерациях. Попробуйте неполную факторизацию Холесского без заливок как предварительный формирователь.

L1 = ichol(A);
[x1,fl1,rr1,it1,rv1] = pcg(A,b,tol,maxit,L1,L1');

fl1 = 0, указание, что pcg сходившийся к требуемому допуску и сделал так в 59 итерациях (значение it1). Поскольку этой матрицей является дискретизированный Лапласиан, однако, использование изменило неполный Холесский, может создать лучший предварительный формирователь. Модифицированная неполная факторизация Холесского создает аппроксимированную факторизацию, которая сохраняет действие оператора на постоянном векторе. Таким образом, norm(A*e-L*(L'*e)) будет приблизительно нуль для e = ones(size(A,2),1) даже при том, что norm(A-L*L','fro')/norm(A,'fro') не близко к нулю. Не необходимо задать тип для этого синтаксиса начиная с nofill значение по умолчанию, но это - хорошая практика.

opts.type = 'nofill';
opts.michol = 'on';
L2 = ichol(A,opts);
e = ones(size(A,2),1);
norm(A*e-L2*(L2'*e))
ans = 3.7983e-14
[x2,fl2,rr2,it2,rv2] = pcg(A,b,tol,maxit,L2,L2');

pcg сходится (fl2 = 0) но только в 38 итерациях. Графический вывод всех трех историй сходимости показывает сходимость.

semilogy(0:maxit,rv0./norm(b),'b.');
hold on
semilogy(0:it1,rv1./norm(b),'r.');
semilogy(0:it2,rv2./norm(b),'k.');
legend('No Preconditioner','IC(0)','MIC(0)');

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

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

L3 = ichol(A, struct('type','ict','droptol',1e-1));
[x3,fl3,rr3,it3,rv3] = pcg(A,b,tol,maxit,L3,L3');
L4 = ichol(A, struct('type','ict','droptol',1e-2));
[x4,fl4,rr4,it4,rv4] = pcg(A,b,tol,maxit,L4,L4');
L5 = ichol(A, struct('type','ict','droptol',1e-3));
[x5,fl5,rr5,it5,rv5] = pcg(A,b,tol,maxit,L5,L5'); 

figure
semilogy(0:maxit,rv0./norm(b),'b-','linewidth',2);
hold on
semilogy(0:it3,rv3./norm(b),'b-.','linewidth',2);
semilogy(0:it4,rv4./norm(b),'b--','linewidth',2);
semilogy(0:it5,rv5./norm(b),'b:','linewidth',2);
legend('No Preconditioner','ICT(1e-1)','ICT(1e-2)', ...
   'ICT(1e-3)','Location','SouthEast');

Отметьте неполный предварительный формирователь Холесского, созданный с допуском отбрасывания 1e-2 обозначается как ICT(1e-2).

Как с заполнением нулями неполного Холесского, пороговая факторизация отбрасывания может извлечь выгоду из модификации (т.е. opts.michol = 'on') поскольку матрица является результатом эллиптического дифференциального уравнения с частными производными. Как с MIC(0), модифицированный пороговый понижающийся неполный Холесский сохранит действие предварительного формирователя на постоянных векторах, который является norm(A*e-L*(L'*e)) будет приблизительно нуль.

Этот пример иллюстрирует использование diagcomp опция ichol.

Неполные факторизации Холесского положительных определенных матриц не всегда существуют. Следующий код создает случайную симметричную положительную определенную матрицу и пытается решить линейную систему с помощью pcg.

S = rng('default');
A = sprandsym(1000,1e-2,1e-4,1);
rng(S);
b = full(sum(A,2));
[x0,fl0,rr0,it0,rv0] = pcg(A,b,1e-6,100);

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

L = ichol(A,struct('type','ict','droptol',1e-3));
Error using ichol
Encountered nonpositive pivot.

Если ichol ломается как выше, можно использовать diagcomp опция, чтобы создать переключенную неполную факторизацию Холесского. Таким образом, вместо того, чтобы создать L таким образом, что L*L' аппроксимирует Aichol с диагональной компенсацией создает L таким образом, что L*L' аппроксимирует M = A + alpha*diag(diag(A)) явным образом не формируя M. Когда неполные факторизации всегда существуют для по диагонали доминирующих матриц, alpha как могут находить, делает M по диагонали доминирующий.

alpha = max(sum(abs(A),2)./diag(A))-2
alpha = 62.9341
L1 = ichol(A, struct('type','ict','droptol',1e-3,'diagcomp',alpha));
[x1,fl1,rr1,it1,rv1] = pcg(A,b,1e-6,100,L1,L1');

Здесь, pcg все еще сбои, чтобы сходиться к желаемому допуску в желаемом количестве итераций, но как график ниже показов, сходимость лучше для pcg с этим предварительным формирователем, чем без предварительного формирователя. Выбор меньшего alpha может помочь. С некоторым экспериментированием мы можем обосноваться на соответствующем значении для alpha.

alpha = .1;
L2 = ichol(A, struct('type','ict','droptol',1e-3,'diagcomp',alpha));
[x2,fl2,rr2,it2,rv2] = pcg(A,b,1e-6,100,L2,L2');

Теперь pcg сходится и график может показать истории сходимости с каждым предварительным формирователем.

semilogy(0:100,rv0./norm(b),'b.');
hold on;
semilogy(0:100,rv1./norm(b),'r.');
semilogy(0:it2,rv2./norm(b),'k.');
legend('No Preconditioner','\alpha \approx 63','\alpha = .1');
xlabel('Iteration Number');
ylabel('Relative Residual');

Советы

  • Фактор, данный этой стандартной программой, может быть полезен как предварительный формирователь для системы линейных уравнений, решаемых итерационными методами, такими как pcg или minres.

  • ichol работает только на разреженные квадратные матрицы

Ссылки

[1] Саад, Yousef. “Предварительное создание условий методов”. Итерационные методы для разреженных линейных систем. PWS Publishing Company, 1996.

[2] Мантойфель, T.A. “Неполный метод факторизации для положительных определенных линейных систем”. Математика. Comput. 34, 473–497, 1980.

Смотрите также

| | |

Для просмотра документации необходимо авторизоваться на сайте