methods.cpp File Reference

Implementation of computational functions for solving linear systems. More...

#include "methods.h"
#include <iostream>
#include <iomanip>
#include <cmath>
#include <algorithm>
#include <string>
#include <vector>

Include dependency graph for methods.cpp:

Go to the source code of this file.

Functions
int	Pivoting (const vector< vector< double > > &m, int current_row, int total_rows)
	Performs partial pivoting and returns the row index with the maximum pivot.
void	Exchange (vector< vector< double > > &m, int row1, int row2)
	Swaps two rows in the matrix and outputs the action.
bool	Eliminate (vector< vector< double > > &m, int current_row, int total_rows, int total_cols)
	Performs elimination on the matrix to form an upper triangular matrix.
int	GaussianElimination (vector< vector< double > > &m, int rows, int cols)
	Performs Gaussian elimination on the augmented matrix with partial pivoting.
bool	BackSubstitution (const vector< vector< double > > &m, int rows, int cols, vector< double > &solution)
	Performs back-substitution to find the solution vector.
int	DetermineRank (const vector< vector< double > > &m, int rows, int cols)
	Determines the rank of the coefficient matrix A (excluding augmented column).
void	ShowGeneralSolution (const vector< vector< double > > &m, int rows, int cols, int rank)
	Displays the general solution for systems with infinitely many solutions.
vector< int >	IdentifyPivots (const vector< vector< double > > &m, int rows, int cols)
	Identifies the pivot columns in the matrix.

Detailed Description

Implementation of computational functions for solving linear systems.

Author

Gilbert Young

Date

2024/09/25

This file implements key algorithms such as Gaussian elimination with partial pivoting, back-substitution, and rank determination. It also includes functionality to display the general solution when the system has infinitely many solutions.

Definition in file methods.cpp.

Function Documentation

BackSubstitution()

bool BackSubstitution	(	const vector< vector< double > > &	m,
		int	rows,
		int	cols,
		vector< double > &	solution )

Performs back-substitution to find the solution vector.

Parameters

m	The upper triangular matrix after Gaussian elimination.
rows	Number of rows in the matrix.
cols	Number of columns in the matrix (including augmented column).
solution	Reference to store the solution vector.

Returns

true If a unique solution exists.

false If the system is inconsistent.

Definition at line 137 of file methods.cpp.

{
    solution.assign(cols - 1, 0.0);
    cout << "Starting back-substitution process..." << endl;
    for (int i = rows - 1; i >= 0; i--)
    {
        // Find the first non-zero coefficient in the row
        int pivot_col = -1;
        for (int j = 0; j < cols - 1; j++)
        {
            if (fabs(m[i][j]) > 1e-12)
            {
                pivot_col = j;
                break;
            }
        }
 
        if (pivot_col == -1)
        {
            if (fabs(m[i][cols - 1]) > 1e-12)
            {
                // Inconsistent equation
                return false;
            }
            else
            {
                // 0 = 0, skip
                continue;
            }
        }
 
        double rhs = m[i][cols - 1];
        cout << "Calculating x" << pivot_col + 1 << ":" << endl;
        for (int j = pivot_col + 1; j < cols - 1; j++)
        {
            cout << "    " << fixed << setprecision(4) << m[i][j] << " * x" << j + 1
                 << " = " << m[i][j] * solution[j] << endl;
            rhs -= m[i][j] * solution[j];
        }
        cout << "    RHS after subtraction = " << rhs << endl;
        solution[pivot_col] = rhs / m[i][pivot_col];
        cout << "    x" << pivot_col + 1 << " = " << rhs << " / " << m[i][pivot_col]
             << " = " << fixed << setprecision(4) << solution[pivot_col] << endl
             << endl;
    }
    return true;
}

DetermineRank()

int DetermineRank	(	const vector< vector< double > > &	m,
		int	rows,
		int	cols )

Determines the rank of the coefficient matrix A (excluding augmented column).

Parameters

m	The augmented matrix [A|b].
rows	Number of rows in the matrix.
cols	Number of columns in the matrix (including augmented column).

Returns

int The rank of the matrix A.

Definition at line 188 of file methods.cpp.

{
    int rank = 0;
    for (int i = 0; i < rows; i++)
    {
        bool non_zero = false;
        for (int j = 0; j < cols - 1; j++)
        {
            if (fabs(m[i][j]) > 1e-12)
            {
                non_zero = true;
                break;
            }
        }
        if (non_zero)
            rank++;
    }
    return rank;
}

Eliminate()

bool Eliminate	(	vector< vector< double > > &	m,
		int	current_row,
		int	total_rows,
		int	total_cols )

Performs elimination on the matrix to form an upper triangular matrix.

Definition at line 47 of file methods.cpp.

{
    double pivot = m[current_row][current_row];
    if (fabs(pivot) < 1e-12)
    {
        // Pivot is too small, cannot eliminate
        return false;
    }
 
    for (int i = current_row + 1; i < total_rows; i++)
    {
        double factor = m[i][current_row] / pivot;
        cout << "Eliminating element in row " << i + 1 << ", column " << current_row + 1 << ":" << endl;
        cout << "Multiplying row " << current_row + 1 << " by " << fixed << setprecision(4) << factor
             << " and subtracting from row " << i + 1 << "." << endl;
        m[i][current_row] = 0.0;
        for (int j = current_row + 1; j < total_cols; j++)
        {
            m[i][j] -= factor * m[current_row][j];
        }
        cout << endl;
    }
    return true;
}

Exchange()

void Exchange	(	vector< vector< double > > &	m,
		int	row1,
		int	row2 )

Swaps two rows in the matrix and outputs the action.

Definition at line 40 of file methods.cpp.

{
    swap(m[row1], m[row2]);
    cout << "Swapping row " << row1 + 1 << " with row " << row2 + 1 << "." << endl;
}

GaussianElimination()

int GaussianElimination	(	vector< vector< double > > &	m,
		int	rows,
		int	cols )

Performs Gaussian elimination on the augmented matrix with partial pivoting.

Parameters

m	Reference to the augmented matrix [A|b] to be modified.
rows	Number of rows in the matrix.
cols	Number of columns in the matrix (including augmented column).

Returns

int Number of row exchanges performed during elimination.

Definition at line 75 of file methods.cpp.

{
    int exchange_count = 0;
    int n = min(rows, cols - 1); // Number of variables
 
    for (int k = 0; k < n; k++)
    {
        cout << "Processing column " << k + 1 << "..." << endl;
 
        // Find the row with the maximum pivot element
        int imax = Pivoting(m, k, rows);
 
        // Swap the current row with the pivot row if necessary
        if (imax != k)
        {
            Exchange(m, k, imax);
            exchange_count++;
        }
        else
        {
            cout << "No need to swap rows for column " << k + 1 << "." << endl;
        }
 
        // Check if pivot element is near zero (singular matrix)
        if (fabs(m[k][k]) < 1e-12)
        {
            cout << "Warning: Pivot element in row " << k + 1 << " is close to zero. The matrix may be singular." << endl;
            continue; // Skip elimination for this pivot
        }
 
        // Eliminate entries below the pivot
        if (!Eliminate(m, k, rows, cols))
        {
            cout << "Elimination failed for column " << k + 1 << "." << endl;
        }
 
        // Display current matrix state
        cout << "Current matrix state:" << endl;
        for (int r = 0; r < rows; r++)
        {
            for (int c = 0; c < cols; c++)
            {
                double coeff = round(m[r][c] * 1e12) / 1e12; // Handle floating-point precision
                if (fabs(coeff - round(coeff)) < 1e-12)
                {
                    cout << static_cast<long long>(round(coeff)) << "\t";
                }
                else
                {
                    cout << fixed << setprecision(2) << coeff << "\t";
                }
            }
            cout << endl;
        }
        cout << "-------------------------------------" << endl;
    }
    return exchange_count;
}

IdentifyPivots()

vector< int > IdentifyPivots	(	const vector< vector< double > > &	m,
		int	rows,
		int	cols )

Identifies the pivot columns in the matrix.

Parameters

m	The matrix after Gaussian elimination.
rows	Number of rows in the matrix.
cols	Number of columns in the matrix (including augmented column).

Returns

std::vector<int> A vector containing the indices of the pivot columns.

Definition at line 320 of file methods.cpp.

{
    vector<int> pivots;
    int n = min(rows, cols - 1);
    for (int i = 0; i < n; i++)
    {
        // Find the pivot column in the current row
        int pivot_col = -1;
        for (int j = 0; j < cols - 1; j++)
        {
            if (fabs(m[i][j]) > 1e-12)
            {
                pivot_col = j;
                break;
            }
        }
        if (pivot_col != -1)
            pivots.push_back(pivot_col);
    }
    return pivots;
}

Pivoting()

int Pivoting	(	const vector< vector< double > > &	m,
		int	current_row,
		int	total_rows )

Performs partial pivoting and returns the row index with the maximum pivot.

Definition at line 23 of file methods.cpp.

{
    int imax = current_row;
    double max_val = fabs(m[current_row][current_row]);
    for (int i = current_row + 1; i < total_rows; i++)
    {
        double val = fabs(m[i][current_row]);
        if (val > max_val)
        {
            imax = i;
            max_val = val;
        }
    }
    return imax;
}

ShowGeneralSolution()

void ShowGeneralSolution	(	const vector< vector< double > > &	m,
		int	rows,
		int	cols,
		int	rank )

Displays the general solution for systems with infinitely many solutions.

Parameters

m	The matrix after Gaussian elimination.
rows	Number of rows in the matrix.
cols	Number of columns in the matrix (including augmented column).
rank	The rank of the coefficient matrix A.

Definition at line 211 of file methods.cpp.

{
    cout << "The system has infinitely many solutions." << endl;
    cout << "Solution space dimension: " << (cols - 1 - rank) << endl;
 
    // Identify pivot columns
    vector<int> pivots = IdentifyPivots(m, rows, cols);
 
    // Identify free variables
    vector<int> free_vars;
    for (int j = 0; j < cols - 1; j++)
    {
        if (find(pivots.begin(), pivots.end(), j) == pivots.end())
        {
            free_vars.push_back(j);
        }
    }
 
    // Assign parameters to free variables
    int num_free = free_vars.size();
    vector<string> params;
    for (int i = 0; i < num_free; i++)
    {
        params.push_back("t" + to_string(i + 1));
    }
 
    // Initialize solution vector with parameters
    vector<double> particular_solution(cols - 1, 0.0);
    vector<vector<double>> basis_vectors;
 
    // Find a particular solution by setting all free variables to 0
    for (int i = rows - 1; i >= 0; i--)
    {
        // Find the first non-zero coefficient in the row
        int pivot_col = -1;
        for (int j = 0; j < cols - 1; j++)
        {
            if (fabs(m[i][j]) > 1e-12)
            {
                pivot_col = j;
                break;
            }
        }
 
        if (pivot_col == -1)
        {
            continue; // 0 = 0, skip
        }
 
        double rhs = m[i][cols - 1];
        for (int j = pivot_col + 1; j < cols - 1; j++)
        {
            rhs -= m[i][j] * particular_solution[j];
        }
        particular_solution[pivot_col] = rhs / m[i][pivot_col];
    }
 
    // Now, find basis vectors by setting each free variable to 1 and others to 0
    for (int i = 0; i < num_free; i++)
    {
        vector<double> basis(cols - 1, 0.0);
        basis[free_vars[i]] = 1.0; // Set the free variable to 1
 
        // Perform back-substitution for pivot variables
        for (int r = rank - 1; r >= 0; r--)
        {
            int pivot_col = pivots[r];
            double rhs = 0.0;
            for (int j = pivot_col + 1; j < cols - 1; j++)
            {
                rhs -= m[r][j] * basis[j];
            }
            basis[pivot_col] = rhs / m[r][pivot_col];
        }
 
        basis_vectors.push_back(basis);
    }
 
    // Display the general solution
    cout << "General solution:" << endl;
    cout << "x = [";
    for (int j = 0; j < cols - 1; j++)
    {
        cout << fixed << setprecision(4) << particular_solution[j];
        if (j < cols - 2)
            cout << ", ";
    }
    cout << "]";
 
    for (int i = 0; i < num_free; i++)
    {
        cout << " + " << params[i] << " * [";
        for (int j = 0; j < cols - 1; j++)
        {
            cout << fixed << setprecision(4) << basis_vectors[i][j];
            if (j < cols - 2)
                cout << ", ";
        }
        cout << "]";
        if (i < num_free - 1)
            cout << " + ";
    }
    cout << endl
         << endl;
}