Preconditioners in SciPy - Practice Problems & Coding Challenges

Challenge - 5 Problems

❓ Predict Output

intermediate

Output of applying Jacobi preconditioner

Given the sparse matrix A and vector b, what is the output of applying the Jacobi preconditioner to b using scipy?

SciPy

import numpy as np
from scipy.sparse import diags
from scipy.sparse.linalg import LinearOperator

A = diags([4, -1, -1], [0, -1, 1], shape=(3,3))
b = np.array([1, 2, 3])

D_inv = diags(1 / A.diagonal())
M = LinearOperator(A.shape, matvec=lambda x: D_inv.dot(x))
result = M.matvec(b)
print(result)

A [0.25 0.5 0.75]

B [4 2 3]

C [0.5 1. 1.5]

D [1 2 3]

💡 Hint

Remember the Jacobi preconditioner uses the inverse of the diagonal elements of A.

The Jacobi preconditioner multiplies the vector by the inverse of the diagonal of A. Since A's diagonal is [4,4,4], the inverse diagonal is [0.25, 0.25, 0.25]. Multiplying by b = [1,2,3] gives [0.25, 0.5, 0.75].

🧠 Conceptual

intermediate

Purpose of Preconditioners in Iterative Solvers

What is the main purpose of using a preconditioner in iterative methods for solving linear systems?

A To convert the matrix into a diagonal matrix

B To increase the size of the matrix for better accuracy

C To change the solution vector to a random guess

D To reduce the condition number of the matrix and speed up convergence

💡 Hint

Think about how preconditioners affect the matrix properties.

Preconditioners improve the matrix's condition number, making iterative solvers converge faster by transforming the system into an easier one to solve.

🔧 Debug

advanced

Identify the error in incomplete LU preconditioner code

What error will this code raise when trying to create an incomplete LU preconditioner with scipy.sparse.linalg.spilu on a singular matrix?

SciPy

import numpy as np
from scipy.sparse import csr_matrix
from scipy.sparse.linalg import spilu

A = csr_matrix([[0, 1], [0, 0]])
ilu = spilu(A)

A ValueError: matrix is not square

B TypeError: unsupported operand type(s)

C RuntimeError: factor is exactly singular

D No error, returns a valid preconditioner

💡 Hint

Consider what happens if the matrix cannot be factorized.

The matrix is singular (zero row), so spilu cannot factorize it and raises a RuntimeError indicating the factor is exactly singular.

❓ data_output

advanced

Output shape of preconditioned vector

If you apply an incomplete LU preconditioner created by spilu on a vector of length 5, what will be the shape of the output vector?

SciPy

import numpy as np
from scipy.sparse import diags
from scipy.sparse.linalg import spilu

A = diags([1, 2, 3, 4, 5], 0)
ilu = spilu(A)
v = np.array([1, 2, 3, 4, 5])
result = ilu.solve(v)
print(result.shape)

A (5,)

B (1, 5)

C (5, 1)

D (25,)

💡 Hint

The preconditioner solves a system of the same size as the input vector.

The preconditioner returns a vector of the same shape as the input vector, which is (5,).

🚀 Application

expert

Choosing a preconditioner for a large sparse symmetric positive definite matrix

You have a large sparse symmetric positive definite matrix A. Which preconditioner is generally the best choice to speed up conjugate gradient solver convergence?

A Jacobi preconditioner

B Incomplete Cholesky factorization

C Incomplete LU factorization without symmetry

D No preconditioner

💡 Hint

Consider the matrix properties and the solver used.

For symmetric positive definite matrices, incomplete Cholesky is preferred as it preserves symmetry and positive definiteness, improving conjugate gradient convergence.

Practice

1.

What is the main purpose of a preconditioner in scipy when solving linear systems?

easy

A. To speed up the convergence of iterative solvers

B. To increase the size of the matrix

C. To change the solution of the system

D. To make the matrix non-square

Solution

1. Step 1: Understand the role of preconditioners

   Preconditioners are used to improve the efficiency of iterative methods by transforming the system into an easier one to solve.

2. Step 2: Identify the effect on convergence

   They help iterative solvers like Conjugate Gradient converge faster by approximating the inverse of the matrix.

3. Final Answer:

   To speed up the convergence of iterative solvers -> Option A

4. Quick Check:

   Preconditioner purpose = speed up convergence [OK]

Hint: Preconditioners help iterative solvers run faster [OK]

Common Mistakes:

• Thinking preconditioners change the solution
• Believing preconditioners increase matrix size
• Confusing preconditioners with matrix transformations that alter shape

2.

Which of the following is the correct way to create a simple Jacobi preconditioner using scipy.sparse.linalg.LinearOperator?

import numpy as np
from scipy.sparse.linalg import LinearOperator

A = np.array([[4, 1], [1, 3]])
M = LinearOperator(shape=A.shape, matvec=lambda x: ...)

easy

A. matvec=lambda x: x / np.diag(A)

B. matvec=lambda x: np.dot(A, x)

C. matvec=lambda x: x * np.diag(A)

D. matvec=lambda x: np.linalg.solve(A, x)

Solution

1. Step 1: Recall Jacobi preconditioner definition

   Jacobi preconditioner uses the inverse of the diagonal elements of matrix A.

2. Step 2: Implement matvec for Jacobi

   Applying the preconditioner means dividing each element of x by the corresponding diagonal element of A.

3. Final Answer:

   matvec=lambda x: x / np.diag(A) -> Option A

4. Quick Check:

   Jacobi preconditioner = divide by diagonal [OK]

Hint: Jacobi preconditioner divides vector by matrix diagonal [OK]

Common Mistakes:

• Using matrix multiplication instead of division
• Trying to solve full system instead of diagonal scaling
• Multiplying by diagonal instead of dividing

3.

Given the following code, what will be the output of print(M.matvec(b))?

import numpy as np
from scipy.sparse.linalg import LinearOperator

A = np.array([[2, 0], [0, 5]])
b = np.array([4, 10])
M = LinearOperator(shape=A.shape, matvec=lambda x: x / np.diag(A))
print(M.matvec(b))

medium

A. [0.5, 2.0]

B. [8.0, 50.0]

C. [2.0, 2.0]

D. [4.0, 10.0]

Solution

1. Step 1: Calculate diagonal of A

   Diagonal elements are [2, 5].

2. Step 2: Apply matvec function

   Divide each element of b by corresponding diagonal: [4/2, 10/5] = [2.0, 2.0].

3. Final Answer:

   [2.0, 2.0] -> Option C

4. Quick Check:

   Vector divided by diagonal = [2.0, 2.0] [OK]

Hint: Divide vector elements by diagonal elements to get output [OK]

Common Mistakes:

• Multiplying instead of dividing
• Confusing vector and matrix multiplication
• Using wrong diagonal values

4.

Identify the error in the following code that attempts to create a Jacobi preconditioner:

import numpy as np
from scipy.sparse.linalg import LinearOperator

A = np.array([[3, 1], [1, 4]])
M = LinearOperator(shape=A.shape, matvec=lambda x: np.diag(A) * x)
print(M.matvec(np.array([1, 2])))

medium

A. Using np.diag(A) incorrectly as a matrix

B. Multiplying by diagonal instead of dividing

C. Shape of LinearOperator is wrong

D. Input vector has wrong size

Solution

1. Step 1: Understand Jacobi preconditioner operation

   Jacobi preconditioner divides vector elements by diagonal elements of A.

2. Step 2: Check given matvec function

   Code multiplies vector by diagonal instead of dividing, which is incorrect.

3. Final Answer:

   Multiplying by diagonal instead of dividing -> Option B

4. Quick Check:

   Jacobi requires division, not multiplication [OK]

Hint: Jacobi preconditioner divides vector by diagonal, not multiply [OK]

Common Mistakes:

• Confusing multiplication with division
• Ignoring element-wise operations
• Not verifying mathematical definition

5.

You want to speed up solving a large sparse system Ax = b using Conjugate Gradient in scipy. Which approach best uses a preconditioner?

from scipy.sparse.linalg import cg, LinearOperator
import numpy as np

# A is large sparse matrix
# b is known vector

# Option 1: Use identity preconditioner
M1 = LinearOperator(A.shape, matvec=lambda x: x)

# Option 2: Use Jacobi preconditioner
diag = A.diagonal()
M2 = LinearOperator(A.shape, matvec=lambda x: x / diag)

# Option 3: Use incomplete Cholesky (not shown)

x, info = cg(A, b, M=M2)

Why is Option 2 preferred over Option 1?

hard

A. Because identity preconditioner is not a LinearOperator

B. Because identity preconditioner changes the solution

C. Because Jacobi preconditioner makes matrix larger

D. Because Jacobi preconditioner approximates inverse and speeds convergence

Solution

1. Step 1: Understand identity preconditioner effect

   Identity preconditioner does nothing; it returns the vector unchanged, so no speedup.

2. Step 2: Understand Jacobi preconditioner effect

   Jacobi approximates the inverse of the diagonal, improving convergence speed of iterative solver.

3. Final Answer:

   Because Jacobi preconditioner approximates inverse and speeds convergence -> Option D

4. Quick Check:

   Jacobi preconditioner = faster convergence [OK]

Hint: Jacobi preconditioner speeds up solver by approximating inverse [OK]

Common Mistakes:

• Thinking identity preconditioner changes solution
• Believing Jacobi increases matrix size
• Confusing LinearOperator requirements