C++ Program To Find Inverse Of A 3x3 Matrix

Understanding how to calculate the inverse of a matrix is a fundamental concept in linear algebra with wide applications. For a 3x3 matrix, this process involves several steps that can be efficiently implemented in C++.

In this article, you will learn how to compute the inverse of a 3x3 matrix using the determinant-adjoint method, complete with a practical C++ program.

Problem Statement

The inverse of a square matrix A, denoted as A⁻¹, is a matrix such that when multiplied by A, it yields the identity matrix I (A * A⁻¹ = A⁻¹ * A = I). Not all square matrices have an inverse; a matrix is invertible if and only if its determinant is non-zero. For a 3x3 matrix, finding the inverse involves calculating its determinant, cofactor matrix, and adjoint matrix.

Background & Knowledge Prerequisites

To understand and implement this solution, a basic understanding of the following is helpful:

• Linear Algebra Concepts:
• Matrix: A rectangular array of numbers.
• Determinant: A scalar value that can be computed from the elements of a square matrix. For a 3x3 matrix, it's a specific calculation involving sums and differences of products of elements.
• Minor: The determinant of the submatrix formed by removing a row and a column from the original matrix.
• Cofactor: A minor multiplied by (-1)^(i+j), where i and j are the row and column indices of the element being removed.
• Cofactor Matrix: A matrix where each element is replaced by its cofactor.
• Adjoint Matrix (or Adjugate Matrix): The transpose of the cofactor matrix.
• Identity Matrix: A square matrix with ones on the main diagonal and zeros elsewhere.
• C++ Programming:
• Arrays (multi-dimensional for matrices).
• Functions.
• Loops (for iterating through matrix elements).
• Conditional statements.

Use Cases

Finding the inverse of a 3x3 matrix is crucial in various fields:

• Computer Graphics: Used for transformations like rotations, scaling, and translations in 3D space, and for camera positioning.
• Physics and Engineering: Solving systems of linear equations, analyzing stress and strain in materials, and signal processing.
• Robotics: Kinematics, trajectory planning, and control systems often rely on matrix inversions.
• Cryptography: Some encryption algorithms utilize matrix operations, including inversion.
• Statistics and Data Science: Solving linear regression problems and eigenvalue decomposition.

Solution Approaches

The standard algebraic method to find the inverse of a 3x3 matrix A is given by the formula:

A⁻¹ = (1 / det(A)) * adj(A)

Where:

• det(A) is the determinant of matrix A.
• adj(A) is the adjoint of matrix A.

We will break this down into the following steps and implement them in C++:

1. Calculate the determinant of the 3x3 matrix.
2. Calculate the cofactor matrix.
3. Transpose the cofactor matrix to get the adjoint matrix.
4. Multiply the adjoint matrix by (1 / determinant) to get the inverse matrix.

Approach: Determinant-Adjoint Method

This method systematically computes the necessary components for the inverse.

One-line summary: Compute the determinant, then the cofactor matrix, transpose it to get the adjoint, and finally scale the adjoint by the inverse of the determinant.

Code example:

// 3x3 Matrix Inverse
#include <iostream>
#include <vector>
#include <iomanip> // For std::fixed and std::setprecision

// Function to calculate the determinant of a 2x2 matrix
double determinant2x2(double a, double b, double c, double d) {
    return a * d - b * c;
}

// Function to calculate the determinant of a 3x3 matrix
double determinant3x3(const std::vector<std::vector<double>>& matrix) {
    double det = 0;
    det += matrix[0][0] * determinant2x2(matrix[1][1], matrix[1][2], matrix[2][1], matrix[2][2]);
    det -= matrix[0][1] * determinant2x2(matrix[1][0], matrix[1][2], matrix[2][0], matrix[2][2]);
    det += matrix[0][2] * determinant2x2(matrix[1][0], matrix[1][1], matrix[2][0], matrix[2][1]);
    return det;
}

// Function to calculate the adjoint of a 3x3 matrix
std::vector<std::vector<double>> adjoint3x3(const std::vector<std::vector<double>>& matrix) {
    std::vector<std::vector<double>> adj(3, std::vector<double>(3));
    
    // Calculate cofactor matrix and then transpose it to get adjoint
    adj[0][0] = determinant2x2(matrix[1][1], matrix[1][2], matrix[2][1], matrix[2][2]);
    adj[0][1] = -determinant2x2(matrix[0][1], matrix[0][2], matrix[2][1], matrix[2][2]); // Transposed from (1,0) cofactor
    adj[0][2] = determinant2x2(matrix[0][1], matrix[0][2], matrix[1][1], matrix[1][2]);  // Transposed from (2,0) cofactor

    adj[1][0] = -determinant2x2(matrix[1][0], matrix[1][2], matrix[2][0], matrix[2][2]); // Transposed from (0,1) cofactor
    adj[1][1] = determinant2x2(matrix[0][0], matrix[0][2], matrix[2][0], matrix[2][2]);
    adj[1][2] = -determinant2x2(matrix[0][0], matrix[0][2], matrix[1][0], matrix[1][2]); // Transposed from (2,1) cofactor

    adj[2][0] = determinant2x2(matrix[1][0], matrix[1][1], matrix[2][0], matrix[2][1]);  // Transposed from (0,2) cofactor
    adj[2][1] = -determinant2x2(matrix[0][0], matrix[0][1], matrix[2][0], matrix[2][1]); // Transposed from (1,2) cofactor
    adj[2][2] = determinant2x2(matrix[0][0], matrix[0][1], matrix[1][0], matrix[1][1]);
    
    return adj;
}

// Function to find the inverse of a 3x3 matrix
std::vector<std::vector<double>> inverse3x3(const std::vector<std::vector<double>>& matrix) {
    double det = determinant3x3(matrix);
    if (det == 0) {
        std::cerr << "Determinant is zero, inverse does not exist." << std::endl;
        // Return an empty or zero matrix to signify failure
        return std::vector<std::vector<double>>(); 
    }

    std::vector<std::vector<double>> adj = adjoint3x3(matrix);
    std::vector<std::vector<double>> inv(3, std::vector<double>(3));
    
    for (int i = 0; i < 3; ++i) {
        for (int j = 0; j < 3; ++j) {
            inv[i][j] = adj[i][j] / det;
        }
    }
    return inv;
}

// Function to print a 3x3 matrix
void printMatrix(const std::vector<std::vector<double>>& matrix) {
    if (matrix.empty()) {
        std::cout << "Matrix is empty or invalid." << std::endl;
        return;
    }
    for (int i = 0; i < 3; ++i) {
        for (int j = 0; j < 3; ++j) {
            std::cout << std::fixed << std::setprecision(4) << matrix[i][j] << "\\t";
        }
        std::cout << std::endl;
    }
}

int main() {
    // Step 1: Define the 3x3 matrix
    std::vector<std::vector<double>> A = {
        {1, 2, 3},
        {0, 1, 4},
        {5, 6, 0}
    };

    std::cout << "Original Matrix A:" << std::endl;
    printMatrix(A);
    std::cout << std::endl;

    // Step 2: Calculate the determinant
    double det_A = determinant3x3(A);
    std::cout << "Determinant of A: " << det_A << std::endl;
    std::cout << std::endl;

    // Step 3: Find the inverse matrix
    if (det_A != 0) {
        std::vector<std::vector<double>> A_inv = inverse3x3(A);
        std::cout << "Inverse of Matrix A (A^-1):" << std::endl;
        printMatrix(A_inv);
        std::cout << std::endl;

        // Optional: Verify A * A^-1 = I
        std::cout << "Verification (A * A^-1):" << std::endl;
        std::vector<std::vector<double>> product(3, std::vector<double>(3, 0.0));
        for (int i = 0; i < 3; ++i) {
            for (int j = 0; j < 3; ++j) {
                for (int k = 0; k < 3; ++k) {
                    product[i][j] += A[i][k] * A_inv[k][j];
                }
            }
        }
        printMatrix(product);
    } else {
        std::cout << "Cannot compute inverse as determinant is zero." << std::endl;
    }

    return 0;
}

Run

Sample output:

Original Matrix A:
1.0000	2.0000	3.0000	
0.0000	1.0000	4.0000	
5.0000	6.0000	0.0000	

Determinant of A: 1.0000

Inverse of Matrix A (A^-1):
-24.0000	18.0000	5.0000	
20.0000	-15.0000	-4.0000	
-5.0000	4.0000	1.0000	

Verification (A * A^-1):
1.0000	0.0000	0.0000	
0.0000	1.0000	0.0000	
0.0000	0.0000	1.0000

Stepwise explanation:

1. determinant2x2 function: This helper function calculates the determinant of a 2x2 matrix, which is frequently needed when computing the 3x3 determinant and cofactors.
2. determinant3x3 function: This function takes a 3x3 matrix and computes its determinant using the cofactor expansion method along the first row. It calls determinant2x2 for the 2x2 minors.
3. adjoint3x3 function:

• This function is crucial for constructing the adjoint matrix. It computes the cofactor of each element of the original matrix.
• The formula for a cofactor C_ij is (-1)^(i+j) * M_ij, where M_ij is the minor.
• Instead of explicitly building a cofactor matrix and then transposing, the adjoint3x3 function directly computes the elements of the adjoint matrix. Remember that adj[i][j] is C_ji (the cofactor of element matrix[j][i]). The code directly calculates C_ji and stores it at adj[i][j], effectively transposing as it goes.
• The signs (-1)^(i+j) are applied by alternating positive and negative determinant2x2 calls.

1. inverse3x3 function:

• First, it calculates the determinant of the input matrix.
• It checks if the determinant is zero. If it is, the inverse does not exist, and an error message is printed.
• If the determinant is non-zero, it calls adjoint3x3 to get the adjoint matrix.
• Finally, it divides each element of the adjoint matrix by the determinant to obtain the inverse matrix.

1. printMatrix function: A utility function to display the contents of a 3x3 matrix in a readable format, using std::fixed and std::setprecision for consistent floating-point output.
2. main function:

• Initializes a sample 3x3 matrix A.
• Prints the original matrix.
• Calculates and prints the determinant of A.
• Calls inverse3x3 to get the inverse matrix A_inv and prints it.
• Includes an optional verification step: multiplying A by A_inv to confirm that the result is close to the identity matrix. Due to floating-point arithmetic, results might not be perfectly 1.0 or 0.0 but very close.

Conclusion

Calculating the inverse of a 3x3 matrix is a multi-step process that combines concepts from linear algebra and programming. The determinant-adjoint method provides a robust way to achieve this, involving the calculation of the determinant, cofactors, and adjoint matrix. This C++ implementation breaks down these mathematical steps into manageable functions, allowing for clear and efficient computation.

Summary

• A matrix inverse exists only if its determinant is non-zero.
• The inverse A⁻¹ = (1 / det(A)) * adj(A).
• Calculating the inverse involves:
• Finding the determinant of the 3x3 matrix.
• Constructing the adjoint matrix (transpose of the cofactor matrix).
• Scaling the adjoint matrix by the reciprocal of the determinant.
• The C++ implementation utilizes helper functions for 2x2 determinants, 3x3 determinants, and adjoint calculation to modularize the code.
• Verification by multiplying the original matrix by its inverse should yield an identity matrix.