ichol как замена холинку: неположительный стержень

Question

В Matlab 2012 cholinc Команда помечена как устаревшая. В предупреждающем сообщении говорится, что он должен быть заменен ichol, До сих пор я использовал cholinc(A,droptol)обычно с droptol=1E-15, В новой версии я пытался использовать ichol(A,struct('droptol',droptol,'type','ict')), который работает большую часть времени, но иногда я получаю предупреждение

Error using ichol
Encountered nonpositive pivot.

Это фундаментальная проблема (т.е. проблема cholinc также имел, но не сообщил) или есть способ сделать ichol вести себя так, как надо cholinc раньше делал?

Answer 1

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

cholinc не терпит срыва, потому что это не настоящий неполный Холецкий. Скорее это вызывает luinc без поворота выбрасывает L, а затем масштабирует результирующее U, чтобы получить некий неполный коэффициент Холецкого (см. документ для cholinc, первый абзац раздела Алгоритмы). Вы можете получить коэффициент, очень похожий на cholinc используя ilu (Обратите внимание, что luinc тоже устарела).

[L,U] = ilu(A,struct('type','ilutp','droptol',droptol,'thresh',0));
R = diag(sqrt(abs(diag(U))))\U;
% Essentially the same as cholinc.

Использование ichol настоятельно рекомендуется, если это возможно. Обратите внимание, что вы можете использовать 'diagcomp' возможность попытаться предотвратить сбои в факторизации, но найти эффективный параметр alpha может потребовать экспериментов. Смотрите документ для ichol для примера. когда ichol не разбивается, он имеет тенденцию быть намного быстрее, поскольку он использует симметрию. Кроме того, он имеет тенденцию производить более редкие факторы, чем cholinc что приводит к более быстрому применению фактора в качестве предварительного кондиционера, что может означать более быстрое время решения с pcg, Например,

M = delsq(numgrid('L',200));
tic; R1 = ichol(M,struct('type','ict','droptol',1e-2,'shape','upper')); toc
% Elapsed time is 0.013809 seconds.
nnz(R1)
% ans = 145632
tic; R = cholinc(M,1e-2); toc
% Elapsed time is 0.167155 seconds.
nnz(R)
% ans = 173851

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

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

tic; L = ichol(M,struct('type','ict','droptol',1e-2)); toc
% Elapsed time is 0.008895 seconds.

В качестве последнего замечания, допустимое отклонение 1e-15, на которое вы ссылаетесь выше, очень мало. Если вы используете именно такой допуск, вам лучше использовать полную факторизацию, такую ​​как chol, ldl, или же lu,