## Posted By

Bangonkali on 03/24/12

# Gaussian Elimination

/ Published in: C++

Numerical Methods application for solving system of equation using Gaussian Elimination based on this Wikipedia article: http://j.mp/GV3PcN

`/*  * File:   main.cpp * Author: Bangonkali * * Created on March 23, 2012, 8:14 PM */#include <iostream>#include <cmath>#include <cstdlib>  // example taken from wikipedia // http://en.wikipedia.org/wiki/Gaussian_elimination#Example using namespace std; //void printMatrix(//        int N,//        float A[20][21]//)//{//    for (int i = 0; i < N; i++)//    {//        for (int j = 0; j <= N; j++)//        {//            cout << A[i][j] << '\t';//        }//        //        cout << endl;//    }//    //    cout << endl;//} int gauss(        int N, // number of unknowns        float A [20] [21], // coefficients and constants        float r[20]){//    printMatrix(N, A);    float multiplier, divider;    bool fin = false;     // forward substitution    for (int m=0; m<=N; m++)    {        for (int i=m+1; i<=N; i++)        {            multiplier = A[i][m]; // current row - root            divider = A[m][m]; // first row - root             if (divider == 0)            {                return 1;            }             for (int j=m; j<=N; j++)            {                if (i == N)                 {                    break;                    fin = true;                }                 A[i][j] = (A[m][j] * multiplier / divider) - A[i][j];//                cout << "A[" << i << "][" << j << "]  = (A[" << m << "][" << j << "] * " << multiplier << " / " << divider << ") - A[" << i << "][" << j << "]" << endl;//                cout << "Current cell " << i << j << endl;//                cout << "Current set " << m << m << endl;//                printMatrix(N, A);            }             for (int j=m; j<=N; j++)            {                A[i-1][j] = A[i-1][j] / divider;            }             if (fin)                break;        }    } //    printMatrix(N, A);     // back substitution    float s = 0;    r[N-1] = A[N-1][N];    int y = 0;    for (int i = N-2; i >= 0; i--)    {        s = 0;        y++;        for (int x = 0; x < N; x++)        {            s = s + (A[i][N-1-x] * r[N-(x+1)]);//            cout << "A[" << i << "][" << N-1-x << "] = " << A[i][N-1-x] << " ";//            cout << " r[" << N-(x+1) << "] = " << r[N-(x+1)] << " " << endl;//            cout << "s = " << s << endl;             if (y == x+1) break;        }        r[i] = A[i][N] - s;//        cout << "r[" << i << "] = A[" << i << "][" << N << "] - " << s << " = " << r[i] << endl;    } //    for (int i = 0; i < N; i++)//    {//        cout << r[i] << endl;//    }  }/* *  */int main(int argc, char** argv) {    float A[20][21];    float X[20];    bool err;    A[0][0] = 2;     A[0][1] = 1;     A[0][2] = -1;     A[0][3] = 8;     A[1][0] = -3;     A[1][1] = -1;     A[1][2] = 2;     A[1][3] = -11;      A[2][0] = -2;     A[2][1] = 1;     A[2][2] = 2;     A[2][3] = -3;      gauss(3, A, X);    for (int i=0; i<3; i++) cout << X[i] << "  ";     return 0;}`