Truncate Values in an array (matrix) below a certain threshold

If you perform numerical calculations like the QR-Decomposition you can come across numerical rounding errors.

In this code snippet we calculate the first householder vector and matrix to perform a QR Decomposition of random integer 4x4 Matrix.

The matrix looks likes this:

$$ A = \begin{bmatrix} -1 & -6 & 5 & -10 \ 7 & 6 & 7 & -2 \ -1 & -10 & 0 & -2 \ -6 & 9 & 6 & -6 \ \end{bmatrix} $$

import numpy as np
import scipy
from scipy.linalg import norm

# Calculate Norm of first column vector
np.random.seed(10)
A = np.random.randint(-10, 10, (4,4))
a_1 = np.array([[-1],[-1],[1],[2]])
n_a1 = norm(a_1)

# Unit Vector
e_1 = np.zeros((4,1))
e_1[0,0] = 1

# Calculate || a_1 - ||a_1||_2 e_1 ||_2
n = norm(a_1 - n_a1*e_1)

# Calculate u_1
u_1 = 1/n * (a_1 - n_a1*e_1)

# Calculate H_1
H_1 = (np.eye(4) - 2*u_1 * u_1.transpose())

# New Submatrix A1
A1 = H_1.dot(A)

Rounding Errors

The Matrix A1 has the following values:

$$ A1 = \begin{bmatrix} -7.181325 & 3.023716 & -8.881784e-16 & -0.755929 \ 5.304513 & 8.475132 & 5.628541e+00 & 0.535574 \ 0.695487 & -12.475132 & 1.371459e+00 & -4.535574 \ -2.609025 & 4.049736 & 8.742919e+00 & -11.071147 \ \end{bmatrix} $$

You see that we have a really small value at the position (0,2). This is an example for a rounding error.

Solution

In order to remove those small values you use the np.where function this replaces the small values with 0 and keeps the other values.

# Calculate A^1
# All values below 1E-15 will be replaced.
A_1 = np.where( np.abs(_) < 1e-15, 0, _)

$$ A_1 = \begin{bmatrix} -7.181325 & 3.023716 & 0.000000 & -0.755929 \ 5.304513 & 8.475132 & 5.628541 & 0.535574 \ 0.695487 & -12.475132 & 1.371459 & -4.535574 \ -2.609025 & 4.049736 & 8.742919 & -11.071147 \ \end{bmatrix} $$

This looks better now and allows to continue with further calculations.