Implement a Matrix Solver

Thomas J. Kennedy

Contents:

1 Overview

In this module we discussed Least Squares approximation. Part of the discussion included designing pseudocode for matrix multiplication. Revisit the pseudocode.

Example 1: Matrix Multiplication Pseudocode
def multiply(lhs, rhs):

    # let result be a lhs.rows by rhs.columns matrix

    for i in 0...result.rows
        for j in 0...result.cols
            for k in 0...lhs.cols
              result[i][j] += lhs[i][k] * rhs[k][j]

    return result

If you did not finish implementing this code either implement it… or make use of NumPy’s np.array.

In this ungraded practice exercise you will:

  1. Select a programming language (Python 3.11+ is strongly encouraged).
  2. Reformulate pseudocode for use with your selected language.
  3. Implement the refined pseudocode in your selected language.

2 Purpose

As a Computer Science student, you have been various software design methodologies (e.g., stepwise refinement) and mathematical problems (e.g., systems of equations as matrices). In this exercise you will take pseudocode, implement it, and analyze its complexity (temporal and spatial).

3 What to Do

  1. Start by creating a README file (either .md or .txt). Use the provided template to start:

    Example 2: Sample ReadMe.md
    # Getting Started
    
    **Language:** {selected_language_here:}
    
    # My Pseudocode
    
    {insert pseudocode here}
    
    
    # Requirements
    
    {List dependencies and libraries here}
    
    
    # Compilation & Execution Instructions
    
    {list required commands}
    

    Each {} represents a placeholder for your content.

  2. Complete, document, and test your matrix multiplication code or use the NumPy np.array.

  3. Take the Gaussian Elimination Pseudocode discussed in The Initial Pseudocode section and rewrite it to work with your selected language and data structures (e.g., np.array.

  4. Place your finalized pseudocode under # My Pseudocode in your README file.

  5. Implement a matrix solver, based on your updated pseudocode.

  6. When you are done… save this code. If you are careful, you will have a significant chunk of code you need for the semester project! A library you can use with a quick C++ #include {} or Java import {} or Python from {} import {} is your goal.

  7. Package your README file, source code, and build configuration into a single zip file or tarball.

3.1 Input

For this exercise you may hardcode the initialization of your $X$, $X^T$, and $Y$ matrices. Use the matrices from Quadratic Example as your inputs.

3.2 Output

Your output should be a single line in the form

phi_hat = c_0 + c_1 * x + c_2 * x^2

where c_0, c_1, and c_2 are the coefficients computed by your code.

4 Grading

Revisit the standard requirements for Programming Exercises (defined in the syllabus). In this assignment you will be graded on:

  1. Quality of your pseudocode
  2. Quality of your implementation
  3. Correctness of your implementation
  4. Use of NumPy (where appropriate)

5 Submission

You will submit your solution through Canvas in accordance with the requirements for Programming Exercises (defined in the syllabus).


6 The Initial Pseudocode

Example 3: Initial Pseudocode

Start

  1. Pivot
    1. Search for max $r_{i,j}$
    2. Swap row i with row largest
  2. Scale
  3. Eliminate
Example 4: Top Level Pseudocode
# Let `matrix_A` be an n by n matrix with an augmented vector in col n
def solve(matrix_A):

    # Iterate through all rows
    for every i in 0 to n-1:

        # Pivot
        idx <- find_largest_row_by_col(A, col_index, num_rows)
        swap_row(A, i, idx)

        scale(A, i, num_cols, A[i][i])
        A[i][i] = 1

        eliminate(A, i, num_rows)
Example 5: Scale
def scale_row(A, row_idx, num_cols, s):

    for every j in 0 to num_cols:
        A[row_idx][j] = A[row_idx][j] / s

Example 6: Eliminate
def eliminate(A, src_row_idx, num_cols, num_rows):

    start_col <- src_row_idx

    for every i in (src_row_idx + 1) to num_rows:

        s <- A[i][start_col]

        for every j in start_col to num_cols:
            A[i][j] = A[i][j] - s * A[src_row_idx][j]

        A[i][start_col] = 0

Example 7: Backsolve
def back_solve(matrix):

    augColIdx <- matrix.cols()  # the augmented column
    lastRow = matrix.rows() - 1

    for i in lastRow down to 1:
        for j <- (i - 1) down to 0:
            s <- matrix[j][i]

            matrix[j][i] -= (s * matrix[i][i])
            matrix[j][augCol] -= (s * matrix[i][augColIdx])