Article

Solving Systems Using Iterative Methods

Use iterative methods to solve systems of equations where the coefficient matrix is sparse.

Overview

The code discussed here solves the following equation by using the iterative method Least Squares Minimum Residual (LSMR) to find the solution to this equation:

Equation to solve using the iterative method Least Squares Minimum Residual (LSMR).

In the equation above, A refers to the 4 x 3 matrix, and b to the right-hand-side vector. The code in this article solves the equation Ax = b by finding x.

Because A is overdetermined (it has more rows than columns), for most right-hand sides there will not be an exact solution. In this case, it is normally required to find the closest solution, in the sense that it minimizes the 2-norm of the error. That is, it solves the optimization min ‖ Ax - b ‖₂. This is known as the least-squares problem.

You could solve this problem through a direct method such as Sparse QR; however, for some problems, a faster method that provides a "good enough" solution is the iterative method LSMR. Unlike direct methods that factorize the matrix A, iterative methods only require the ability to multiply by the matrix (and its transpose, A). They move through a sequence of approximate solutions, converging to the correct answer. However, these methods run into numerical difficulties more often than direct methods. Resolving these issues requires expert knowledge, and is sometimes impossible. The most common method to improve and accelerate convergence is to use a preconditioner—an operator that approximates A¹, and is hence problem specific. For least-squares problems, just using a diagonal matrix with entries equal to the 2-norm of each column is often sufficient, and is the method discussed here.

Although you do not need the matrix explicitly to use iterative methods, the Accelerate framework offers the ability to supply a matrix explicitly if it is available, reducing the code required. The example takes advantage of this ability.

Create the Matrix

Use the code below—which is discussed in detail in Creating Sparse Matrices—to define the unsymmetric matrix A:

int rowIndices[]    = {   0,    1,   3,   0,   1,    2,   3,   1,   2 };
double values[]     = { 2.0, -0.2, 2.5, 1.0, 3.2, -0.1, 1.1, 1.4, 0.5 };
 
long columnStarts[] = {   0,              3,                   7,        9};
 
  SparseMatrix_Double A = {
    .structure = {
      .rowCount     = 4,
      .columnCount  = 3,
      .columnStarts = columnStarts,
      .rowIndices   = rowIndices,
      // Matrix meta-data
      .attributes = {
        .kind = SparseOrdinary,          
      },
      .blockSize = 1
    },
    .data = values
  };

Solve the Equation

You define the right-hand-side and solution vectors as arrays wrapped in flexible data types. The initial values of x are used as an initial guess of the solution. If no good estimate is known, initialize all the values to zero:

double bValues[] = { 1.200, 1.013, 0.205, -0.172 };
DenseVector_Double b = {
    .count = 4,    
    .data = bValues 
  };
 
double xValues[] = {   0.0,   0.0,   0.0 };
DenseVector_Double x = {
    .count = 3,     
    .data = xValues 
  };

With the matrix and vectors created, you are ready to perform the full LSMR iteration and iterate over the results in x. SparseSolve returns a status which, if equal to SparseIterativeConverged, indicates that the vector, x, contains the solution.

__auto_type status = SparseSolve( SparseLSMR(), A, b, x, SparsePreconditionerDiagScaling);
if(status!=SparseIterativeConverged) {    
    printf("Failed to converge. Returned with error %d\n", status);}
else {
    printf("x = "); 
    for(int i=0; i<x.count; i++) 
        printf(" %.2f", x.data[i]); printf("\n");
}

You should get this output: x = 0.10 0.20 0.30.

See Also

Iterative Sparse Solving Methods

Implementing Iterative Methods

Use iterative methods to solve large problems faster and with a lower memory overhead than with direct methods.

Preconditioners

Create preconditioners for iterative solves.

Sparse Iterative Methods

Select a suitable iterative method to solve a system.

Iterative Sparse Solve Functions

Solve a system using an iterative method.

Sparse Iterate Functions

Perform a single iteration of the specified iterative method.