pcg

Предобусловленный метод сопряженных градиентов

Синтаксис

x = pcg (A, b)
pcg (A, b, tol)
pcg (A, b, tol, maxit)
pcg (A, b, tol, maxit, M)
pcg (A, b, tol, maxit, M1, M2)
pcg (A, b, tol, maxit, M1, M2, x0)
[x, флаг] = pcg (A, b...)
[x, флаг, relres] = pcg (A, b...)
[x, флаг, relres, проход] = pcg (A, b...)
[x, флаг, relres, проход, resvec] = pcg (A, b...)

Описание

x = pcg(A,b) пытается решить систему линейных уравнений A*x=b для x. n-by-n матрица коэффициентов A должен быть симметричным и положительный определенный, и должен также быть большим и разреженным. Вектор - столбец b должен иметь длину n. Также можно задать A, чтобы быть указателем на функцию, afun, таким, что afun(x) возвращает A*x.

Параметризация Функций объясняет, как предоставить дополнительные параметры функциональному afun, а также функцию перед формирователем mfun, описанный ниже, при необходимости.

Если pcg сходится, сообщение к тому эффекту отображено. Если pcg не удается сходиться после максимального количества итераций или остановов по любой причине, предупреждающее сообщение распечатано, отобразив относительный остаточный norm(b-A*x)/norm(b) и номер итерации в который метод, остановленный или не пройдено.

pcg(A,b,tol) задает допуск метода. Если tol является [], то pcg использует значение по умолчанию, 1e-6.

pcg(A,b,tol,maxit) задает максимальное количество итераций. Если maxit является [], то pcg использует значение по умолчанию, min(n,20).

pcg(A,b,tol,maxit,M) и pcg(A,b,tol,maxit,M1,M2) используют симметричный положительный определенный предварительный формирователь M или M = M1*M2 и эффективно решают систему inv(M)*A*x = inv(M)*b для x. Если M является [] затем, pcg не применяет предварительного формирователя. M может быть указателем на функцию mfun, таким образом, что mfun(x) возвращает M\x.

pcg(A,b,tol,maxit,M1,M2,x0) задает исходное предположение. Если x0 является [], то pcg использует значение по умолчанию, все-нулевой вектор.

[x,flag] = pcg(A,b,...) также возвращает флаг сходимости.

Флаг

Сходимость

0

pcg сходился к желаемому допуску tol в итерациях maxit .

1

pcg выполнил итерации времен maxit, но не сходился.

2

Предварительный формирователь M был плохо обусловлен.

3

pcg застоялся. (Два последовательных выполняют итерации, было то же самое.)

4

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

Каждый раз, когда flag не является 0, решение, возвращенный x то, что с минимальной невязкой нормы, вычисленной по всем итерациям. Никакие сообщения не отображены, если flag вывод задан.

[x,flag,relres] = pcg(A,b,...) также возвращает относительный остаточный norm(b-A*x)/norm(b). Если flag является 0, relres <= tol.

[x,flag,relres,iter] = pcg(A,b,...) также возвращает номер итерации, в котором x был вычислен, где 0 <= iter <= maxit.

[x,flag,relres,iter,resvec] = pcg(A,b,...) также возвращает вектор остаточных норм в каждой итерации включая norm(b-A*x0).

Примеры

Используя pcg с Большими Матрицами

Этот пример показывает, как использовать pcg с матричным входным параметром и с указателем на функцию.

n1 = 21; 
A = gallery('moler',n1);  
b1 = sum(A,2);
tol = 1e-6;  
maxit = 15;  
M1 = spdiags((1:n1)',0,n1,n1);
[x1,flag1,rr1,iter1,rv1] = pcg(A,b1,tol,maxit,M1);

Также можно использовать следующую функцию вместо матричного A:

function y = applyMoler(x)
y = x;
y(end-1:-1:1) = y(end-1:-1:1) - cumsum(y(end:-1:2));
y(2:end) = y(2:end) - cumsum(y(1:end-1));

При помощи этой функции можно решить большие системы более эффективно, когда нет никакой потребности сохранить целый матричный A:

n2 = 21;
b2 = applyMoler(ones(n2,1));
tol = 1e-6;
maxit = 15;
M2 = spdiags((1:n2)',0,n2,n2);
[x2,flag2,rr2,iter2,rv2] = pcg(@applyMoler,b2,tol,maxit,M2);

Используя pcg с Предварительным формирователем

Этот пример демонстрирует, как использовать матрицу перед формирователем с pcg.

Создайте входную матрицу и попытайтесь решить систему с pcg.

A = delsq(numgrid('S',100));
b = ones(size(A,1),1);
[x0,fl0,rr0,it0,rv0] = pcg(A,b,1e-8,100);

fl0 является 1, потому что pcg не сходится к требуемому допуску 1e-8 в требуемых максимальных 100 итерациях. Предварительный формирователь может заставить систему сходиться более быстро.

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

L = ichol(A);
[x1,fl1,rr1,it1,rv1] = pcg(A,b,1e-8,100,L,L');

fl1 является 0, потому что pcg управляет относительной невязкой к 9.8e-09 (значение rr1), который является меньше, чем требуемый допуск 1e-8 в семьдесят седьмой итерации (значение it1), когда предобусловленный заполнением нулями неполной факторизации Холесского. rv1(1) = norm(b) и rv1(78) = norm(b-A*x1).

Предыдущая матрица представляет дискретизацию Лапласиана на 100x100 сетка с граничными условиями Дирихле. Это означает, что измененный неполный предварительный формирователь Холесского может выполнить еще лучше.

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

L = ichol(A,struct('michol','on'));
[x2,fl2,rr2,it2,rv2] = pcg(A,b,1e-8,100,L,L');

В этом случае вы достигаете сходимости только в сорока семи итерациях.

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

figure;
semilogy(0:it0,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)');
xlabel('iteration number');
ylabel('relative residual');
hold off;

Ссылки

[1] Барретт, Ягода R., M., Т. F. Чан, и др., Шаблоны для Решения Линейных систем: Стандартные блоки для Итеративных Методов, SIAM, Филадельфия, 1994.

Расширенные возможности

Представлено до R2006a

Была ли эта тема полезной?