Бисопряженные градиенты стабилизировали метод
x = bicgstab (A, b)
bicgstab (A, b, tol)
bicgstab (A, b, tol, maxit)
bicgstab (A, b, tol, maxit, M)
bicgstab (A, b, tol, maxit, M1, M2)
bicgstab (A, b, tol, maxit, M1, M2, x0)
[x, флаг] = bicgstab (A, b...)
[x, флаг, relres] = bicgstab (A, b...)
[x, флаг, relres, проход] = bicgstab (A, b...)
[x, флаг, relres, проход, resvec] = bicgstab (A, b...)
x = bicgstab(A,b) пытается решить систему линейных уравнений A*x=b для x. n-by-n матрица коэффициентов A должен быть квадратным и должен быть большим и разреженным. Вектор - столбец b должен иметь длину n. A может быть указателем на функцию, afun, таким, что afun(x) возвращает A*x.
Параметризация Функций объясняет, как предоставить дополнительные параметры функциональному afun, а также функцию перед формирователем mfun, описанный ниже, при необходимости.
Если bicgstab сходится, сообщение к тому эффекту отображено. Если bicgstab не удается сходиться после максимального количества итераций или остановов по любой причине, предупреждающее сообщение распечатано, отобразив относительный остаточный norm(b-A*x)/norm(b) и номер итерации в который метод, остановленный или не пройдено.
bicgstab(A,b,tol) задает допуск метода. Если tol является [], то bicgstab использует значение по умолчанию, 1e-6.
bicgstab(A,b,tol,maxit) задает максимальное количество итераций. Если maxit является [], то bicgstab использует значение по умолчанию, min(n,20).
bicgstab(A,b,tol,maxit,M) и предварительный формирователь использования bicgstab(A,b,tol,maxit,M1,M2) M или M = M1*M2 и эффективно решают систему inv(M)*A*x = inv(M)*b для x. Если M является [] затем, bicgstab не применяет предварительного формирователя. M может быть указателем на функцию mfun, такой, что mfun(x) возвращает M\x.
bicgstab(A,b,tol,maxit,M1,M2,x0) задает исходное предположение. Если x0 является [], то bicgstab использует значение по умолчанию, весь нулевой вектор.
[x,flag] = bicgstab(A,b,...) также возвращает флаг сходимости.
Флаг | Сходимость |
|---|---|
|
|
|
|
| Предварительный формирователь |
|
|
| Один из скаляров, вычисленных во время |
Каждый раз, когда flag не является 0, решение, возвращенный x то, что с минимальной невязкой нормы, вычисленной по всем итерациям. Никакие сообщения не отображены, если flag вывод задан.
[x,flag,relres] = bicgstab(A,b,...) также возвращает относительный остаточный norm(b-A*x)/norm(b). Если flag является 0, relres <= tol.
[x,flag,relres,iter] = bicgstab(A,b,...) также возвращает номер итерации, в котором x был вычислен, где 0 <= iter <= maxit. iter может быть целочисленный + 0.5, указав на сходимость на полпути через итерацию.
[x,flag,relres,iter,resvec] = bicgstab(A,b,...) также возвращает вектор остаточных норм в каждой половине итерации, включая norm(b-A*x0).
Этот пример сначала решает Ax = b путем обеспечения A и предварительному формирователю M1 непосредственно в качестве аргументов.
Код:
A = gallery('wilk',21);
b = sum(A,2);
tol = 1e-12;
maxit = 15;
M1 = diag([10:-1:1 1 1:10]);
x = bicgstab(A,b,tol,maxit,M1);отображения сообщение:
bicgstab converged at iteration 12.5 to a solution with relative residual 2e-014.
Этот пример заменяет матричный A в предыдущем примере с указателем на функцию матричного векторного произведения afun и предварительный формирователь M1 с указателем на функцию backsolve mfun. Пример содержится в файле run_bicgstab это
Вызовы bicgstab с указателем на функцию @afun в качестве его первого аргумента.
Содержит afun и mfun как вложенные функции, так, чтобы все переменные в run_bicgstab были доступны afun и mfun.
Следующее показывает код для run_bicgstab:
function x1 = run_bicgstab
n = 21;
b = afun(ones(n,1));
tol = 1e-12;
maxit = 15;
x1 = bicgstab(@afun,b,tol,maxit,@mfun);
function y = afun(x)
y = [0; x(1:n-1)] + ...
[((n-1)/2:-1:0)'; (1:(n-1)/2)'].*x + ...
[x(2:n); 0];
end
function y = mfun(r)
y = r ./ [((n-1)/2:-1:1)'; 1; (1:(n-1)/2)'];
end
endКогда вы входите
x1 = run_bicgstab;
MATLAB отображает сообщение
bicgstab converged at iteration 12.5 to a solution with relative residual 2e-014.
Этот пример демонстрирует использование предварительного формирователя.
Загрузите west0479, действительное 479 479 несимметричная разреженная матрица.
load west0479;
A = west0479;Задайте b так, чтобы истинное решение было вектором из всех единиц.
b = full(sum(A,2));
Установите допуск и максимальное количество итераций.
tol = 1e-12; maxit = 20;
Используйте bicgstab, чтобы найти решение в требуемом допуске и количестве итераций.
[x0,fl0,rr0,it0,rv0] = bicgstab(A,b,tol,maxit);
fl0 равняется 1, потому что bicgstab не сходится к требуемому допуску 1e-12 в требуемых 20 итерациях. На самом деле поведение bicgstab так плохо, что исходное предположение (x0 = zeros(size(A,2),1)) является лучшим решением и возвращено, как обозначено it0 = 0. MATLAB® хранит остаточную историю в rv0.
Постройте график поведения bicgstab.
semilogy(0:0.5:maxit,rv0/norm(b),'-o'); xlabel('Iteration number'); ylabel('Relative residual');
График показывает, что решение не сходится. Можно использовать предварительный формирователь, чтобы улучшить результат.
Создайте предварительный формирователь с ilu, поскольку A несимметричен.
[L,U] = ilu(A,struct('type','ilutp','droptol',1e-5));
Error using ilu There is a pivot equal to zero. Consider decreasing the drop tolerance or consider using the 'udiag' option.
MATLAB не может создать неполный LU, когда это привело бы к сингулярному фактору, который бесполезен как предварительный формирователь.
Можно попробовать еще раз с уменьшенным допуском отбрасывания, как обозначено сообщением об ошибке.
[L,U] = ilu(A,struct('type','ilutp','droptol',1e-6)); [x1,fl1,rr1,it1,rv1] = bicgstab(A,b,tol,maxit,L,U);
fl1 0, потому что bicgstab управляет относительной невязкой к 5.9829e-014 (значение rr1). Относительная невязка является меньше, чем предписанный допуск 1e-12 в третьей итерации (значение it1), когда предобусловленный неполной LU-факторизацией с допуском отбрасывания 1e-6. Выводом rv1(1) является norm(b), и выводом rv1(7) является norm(b-A*x2), поскольку bicgstab использует половину итераций.
Можно следовать, прогресс bicgstab путем графического изображения относительных невязок в каждой итерации, начинающей с первоначальной сметы (выполните итерации номера 0).
semilogy(0:0.5:it1,rv1/norm(b),'-o'); xlabel('Iteration Number'); ylabel('Relative Residual');
[1] Барретт, Ягода R., M., Канал T.F., и др., Шаблоны для Решения Линейных систем: Стандартные блоки для Итеративных Методов, SIAM, Филадельфия, 1994.
[2] ван дер Ворст, H.A., "BI-CGSTAB: быстрый и гладко сходящийся вариант CG BI для решения несимметричных линейных систем", SIAM J. Научный Закон Comput., март 1992, Издание 13, № 2, стр 631–644.