MATLAB Function Reference
pcg

Syntax

• ```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)
pcg(A,b,tol,maxit,M1,M2,x0,p1,p2,...)
[x,flag] = pcg(A,b,tol,maxit,M1,M2,x0,p1,p2,...)
[x,flag,relres] = pcg(A,b,tol,maxit,M1,M2,x0,p1,p2,...)
[x,flag,relres,iter] = pcg(A,b,tol,maxit,M1,M2,x0,p1,p2,...)
[x,flag,relres,iter,resvec] = pcg(A,b,tol,maxit,M1,M2,x0,p1,p2,...)
```

Description

```x = pcg(A,b) ``` attempts to solve the system of linear equations `A*x=b` for `x`. The `n`-by-`n` coefficient matrix `A` must be symmetric and positive definite, and should also be large and sparse. The column vector `b` must have length `n`. `A` can be a function `afun` such that `afun(x)` returns `A*x`.

If `pcg` converges, a message to that effect is displayed. If `pcg` fails to converge after the maximum number of iterations or halts for any reason, a warning message is printed displaying the relative residual `norm(b-A*x)/norm(b)` and the iteration number at which the method stopped or failed.

```pcg(A,b,tol) ``` specifies the tolerance of the method. If `tol` is `[]`, then `pcg` uses the default, `1e-6`.

```pcg(A,b,tol,maxit) ``` specifies the maximum number of iterations. If `maxit` is `[]`, then `pcg` uses the default, `min(n,20)`.

```pcg(A,b,tol,maxit,M) and pcg(A,b,tol,maxit,M1,M2) ``` use symmetric positive definite preconditioner `M` or `M = M1*M2` and effectively solve the system `inv(M)*A*x = inv(M)*b` for `x`. If `M` is `[]` then `pcg` applies no preconditioner. `M` can be a function that returns `M`\`x`.

```pcg(A,b,tol,maxit,M1,M2,x0) ``` specifies the initial guess. If `x0` is `[]`, then `pcg` uses the default, an all-zero vector.

```pcg(afun,b,tol,maxit,m1fun,m2fun,x0,p1,p2,...) ``` passes parameters `p1,p2,...` to functions `afun(x,p1,p2,...)`, `m1fun(x,p1,p2,...)`, and `m2fun(x,p1,p2,...)`.

```[x,flag] = pcg(A,b,tol,maxit,M1,M2,x0) ``` also returns a convergence flag.

 Flag Convergence `0` `pcg `converged to the desired tolerance `tol` within `maxit `iterations. `1` `pcg` iterated `maxit` times but did not converge. `2` Preconditioner `M` was ill-conditioned. `3` `pcg` stagnated. (Two consecutive iterates were the same.) `4` One of the scalar quantities calculated during `pcg `became too small or too large to continue computing.

Whenever `flag` is not `0`, the solution `x` returned is that with minimal norm residual computed over all the iterations. No messages are displayed if the `flag` output is specified.

```[x,flag,relres] = pcg(A,b,tol,maxit,M1,M2,x0) ``` also returns the relative residual `norm(b-A*x)/norm(b)`. If `flag` is `0`, `relres <= tol`.

```[x,flag,relres,iter] = pcg(A,b,tol,maxit,M1,M2,x0) ``` also returns the iteration number at which `x` was computed, where `0 <= iter <= maxit`.

```[x,flag,relres,iter,resvec] = pcg(A,b,tol,maxit,M1,M2,x0) ``` also returns a vector of the residual norms at each iteration including `norm(b-A*x0)`.

Examples

Example 1.

• ```A = gallery('wilk',21);
b = sum(A,2);
tol = 1e-12;
maxit = 15;
M = diag([10:-1:1 1 1:10]);
[x,flag,rr,iter,rv] = pcg(A,b,tol,maxit,M);
```

Alternatively, use this one-line matrix-vector product function

• ```function y = afun(x,n)
y = [0;
x(1:n-1)] + [((n-1)/2:-1:0)';
(1:(n-1)/2)'].*x + [x(2:n);
0];
```

and this one-line preconditioner backsolve function

• ```function y = mfun(r,n)
y = r ./ [((n-1)/2:-1:1)'; 1; (1:(n-1)/2)'];
```

as inputs to `pcg`

• ```[x1,flag1,rr1,iter1,rv1] = pcg(@afun,b,tol,maxit,@mfun,...
[],[],21);
```

Example 2.

• ```A = delsq(numgrid('C',25));
b = ones(length(A),1);
[x,flag] = pcg(A,b)
```

`flag` is `1` because `pcg` does not converge to the default tolerance of `1e-6` within the default 20 iterations.

• ```R = cholinc(A,1e-3);
[x2,flag2,relres2,iter2,resvec2] = pcg(A,b,1e-8,10,R',R)
```

`flag2` is `0` because `pcg` converges to the tolerance of `1.2e-9` (the value of `relres2`) at the sixth iteration (the value of `iter2`) when preconditioned by the incomplete Cholesky factorization with a drop tolerance of `1e-3`. `resvec2(1) = norm(b)` and `resvec2(7) = norm(b-A*x2)`. You can follow the progress of `pcg` by plotting the relative residuals at each iteration starting from the initial estimate (iterate number 0).

• ``````semilogy(0:iter2,resvec2/norm(b),'-o')
```xlabel('iteration number')
ylabel('relative residual')

```

`bicg`, `bicgstab`, `cgs`, `cholinc`, `gmres`, `lsqr`, `minres`, `qmr`, `symmlq`
`@` (function handle), `\` (backslash)