pius
pius

Reputation: 119

Gauss-Jordan Elimination in C++

I was using Gauss-Jordan elimination in C++ to solve a system of linear equations. Code works fine. Was wondering why Lines 1,2,3 in void gauss() can't be replaced by Line 4 (getting incorrect output after doing so)?

#include <iostream>
using namespace std;
class Gauss
{
    float a[50][50];
    int n;
public:
    void accept()
    {
        cout<<"Enter no. of variables: ";
        cin>>n;
        for(int i=0;i<n;i++)
        {
            for(int j=0;j<n+1;j++)
            {
                if(j==n)
                    cout<<"Constant no."<<i+1<<" = ";
                else
                    cout<<"a["<<i+1<<"]["<<j+1<<"] = ";
                cin>>a[i][j];
            }
        }
    }
    void display()
    {
        for(int i=0;i<n;i++)
        {
            cout<<"\n";
            for(int j=0;j<n+1;j++)
            {
                if(j==n)
                    cout<<" ";
                cout<<a[i][j]<<"\t";
            }
        }
    }

    void gauss()//converting augmented matrix to row echelon form
    {
        float temp;//Line 1
        for(int i=0;i<n;i++)
        {
            for(int j=i+1;j<n;j++)
            {
                temp=a[j][i]/a[i][i];//Line 2
                for(int k=i;k<n+1;k++)
                {
                      a[j][k]-=temp*a[i][k];//Line 3
                    //a[j][k]-=a[j][i]*a[i][k]/a[i][i];//Line 4
                }
            }
        }
    }

    void EnterJordan()//converting to reduced row echelon form
    {
        float temp;
        for(int i=n-1;i>=0;i--)
        {

            for(int j=i-1;j>=0;j--)
            {
                temp=a[j][i]/a[i][i];
                for(int k=n;k>=i;k--)
                {
                    a[j][k]-=temp*a[i][k];
                }
            }
        }

        float x[n];
        for(int i=0;i<n;i++)//making leading coefficients zero
            x[i]=0;
        for(int i=0;i<n;i++)
        {
            for(int j=0;j<n+1;j++)
            {
                if(x[i]==0&&j!=n)
                    x[i]=a[i][j];
                if(x[i]!=0)
                    a[i][j]/=x[i];
            }
        }
    }
    void credits()
    {
        for(int i=0;i<n;i++)
        {
            cout<<"\nx"<<i+1<<" = "<<a[i][n]<<endl;
        }
    }

};

int main()
{
    Gauss obj;
    obj.accept();
    cout<<"\n\nAugmented matrix: \n\n\n";
    obj.display();
    obj.gauss();
    cout<<"\n\nRow Echelon form: \n\n\n";
    obj.display();
    obj.EnterJordan();
    cout<<"\n\nReduced row echelon form:\n\n\n";
    obj.display();
    cout<<"\n\nSolution: \n\n\n";
    obj.credits();
    return 0;
}

Note: My code doesn't take into consideration the problem of division when the pivot is zero (I'm choosing the diagonal element as the pivot every time). For the particular example I tried however, such a case was not encountered.

Augmented matrix is :

 2   1  -1    8 
-3  -1   2   -11    
-2   1   2   -3

The output matrix is :

1   0   0    2  
0   1   0    3  
0   0   1    -1 

and the solution is :

x1 = 2

x2 = 3

x3 = -1

Using Line 4, the output matrix is :

1   0   0    -0.75  
0   1   -0   8  
0   0   1    -1.5   

and the solution is :

x1 = -0.75

x2 = 8

x3 = -1.5

Upvotes: 6

Views: 12178

Answers (1)

Ben Voigt
Ben Voigt

Reputation: 283684

Your Line #4 reads from a[j][i] many times, and the first time through the inner loop, when k == i, changes a[j][i] to 0.0f, thus breaking the next n-i iterations.

Reordering reads of a variable with a write to the same location is not safe.

Upvotes: 3

Related Questions