Sparse Linear System Solver

Build a GPU-accelerated solver for sparse linear systems using iterative methods.

Problem Description

Create a Python module that solves sparse linear systems of equations Ax = b where A is a sparse matrix. The solver should leverage GPU acceleration for improved performance on large systems.

The module should support:

• Solving sparse linear systems using an iterative approach
• Optional preconditioning to improve convergence
• Configurable convergence tolerance
• Ability to track solution progress

Requirements

Input Format

The solver should accept:

• A sparse matrix A in COO (coordinate) format, represented as:
  • row_indices: array of row indices
  • col_indices: array of column indices
  • values: array of non-zero values
  • shape: tuple (n_rows, n_cols)
• A right-hand side vector b (dense array)
• Optional tolerance for convergence (default: 1e-6)
• Optional maximum iterations (default: 1000)
• Optional preconditioner type: "identity" or "diagonal" (default: "identity")

Output Format

The solver should return:

• Solution vector x
• Number of iterations performed
• Final residual norm

Functionality

The implementation should:

1. Transfer sparse matrix data and vectors to GPU
2. Apply iterative solver with optional preconditioning
3. Check convergence based on residual norm
4. Return solution and convergence information

Test Cases

• Given a 3x3 sparse matrix A = [[4, 1, 0], [1, 3, 1], [0, 1, 2]] with b = [1, 2, 3], solve for x and verify the residual norm is below tolerance @test

• Given a 5x5 diagonal sparse matrix with values [2, 3, 4, 5, 6] and b = [2, 6, 12, 20, 30], solve for x and verify x equals [1, 2, 3, 4, 5] within tolerance @test

• Given a sparse matrix and applying diagonal preconditioning, verify that the solver converges in fewer iterations than without preconditioning @test

@generates

API

def solve_sparse_system(
    row_indices,
    col_indices,
    values,
    shape,
    b,
    tolerance=1e-6,
    max_iterations=1000,
    preconditioner="identity"
):
    """
    Solves the sparse linear system Ax = b using GPU-accelerated iterative methods.

    Args:
        row_indices: Array of row indices for sparse matrix (COO format)
        col_indices: Array of column indices for sparse matrix (COO format)
        values: Array of non-zero values for sparse matrix
        shape: Tuple (n_rows, n_cols) defining matrix dimensions
        b: Right-hand side vector (dense array)
        tolerance: Convergence tolerance for residual norm (default: 1e-6)
        max_iterations: Maximum number of iterations (default: 1000)
        preconditioner: Type of preconditioner - "identity" or "diagonal" (default: "identity")

    Returns:
        Tuple of (solution_vector, iterations, final_residual)
        - solution_vector: Solution x as numpy array
        - iterations: Number of iterations performed (int)
        - final_residual: Final residual norm (float)

    Raises:
        ValueError: If matrix is not square or dimensions don't match
        RuntimeError: If solver fails to converge within max_iterations
    """
    pass

Dependencies { .dependencies }

pycuda { .dependency }

Provides GPU-accelerated sparse matrix operations and iterative solvers.

@satisfied-by