8  Python and NumPy Primer

Condensed reference for the toolset used throughout these notes.

8.1 Arrays and Basic Operations

NoteDefinition: NumPy Arrays

The ndarray is the fundamental structure.

  • Creation: np.array([1, 2]), np.zeros((m, n)), np.random.randn(m, n).

  • Indexing: 0-indexed. A[i, j] for entries; A[i, :] for rows; A[:, j] for columns.

Tip

Shapes and Rank-1 Arrays. A.shape returns (m, n). Note that a 1D array has shape (n,), which is not equivalent to a column vector (n, 1). These `rank-1 arrays'' are a common source of broadcasting bugs. Usex[:, None]orx.reshape(-1, 1)` to convert to 2D.

WarningExercise

  1. Create \(4 \times 4\) identity matrix. Verify shape and dtype.

  2. Extract row 2 and column 3 from a random \(3 \times 3\) \(\bA\).

  3. Confirm that np.zeros((3,)) fails where a (3, 1) column vector is required in a matrix-vector product.

8.2 Linear Algebra Operations

Core routines in numpy.linalg:

Tip

Performance Law: Vectorization. Explicit for loops in Python are slow. Push heavy lifting into NumPy/BLAS routines by vectorizing operations. For example, use A @ B instead of triple-nested loops.

Tip

Stability: solve vs inv. np.linalg.solve(A, b) is faster and more stable than np.linalg.inv(A) @ b. Explicit inversion is almost never desirable in scientific computation.

WarningExercise

  1. Solve a \(2 \times 2\) system and verify the residual \(\|\bA\bx - \bb\|\).

  2. Profile solve vs inv @ b for \(n=1000\).

  3. Compare np.linalg.eig and np.linalg.eigh (optimized for symmetric matrices).

8.3 Key Idioms

  • Element-wise: * and ** are entry-wise. Matrix powers use np.linalg.matrix\_power.

  • Broadcasting: Automatic dimension expansion. A * v scales columns of \(\bA\) by entries of \(v\).

  • Views vs. Copies: B = A is a view (modifying B changes A). Use B = A.copy() for independent data. Slices are also views.

WarningExercise

  1. Let \(\bA = \begin{pmatrix} 1 & 2 \\ 3 & 4 \end{pmatrix}\). Compare A * A and A @ A.

  2. Demonstrate the ``view pitfall’’ by modifying a slice of a matrix.

  3. Use broadcasting to perform diagonal scaling \(\bA\bD\) without constructing \(\bD = \text{diag}(d)\).