Segmentation Fault: Implementation of Gauss Elimination Algorithm

[zak100]

using namespace std;


void GaussianElimination(double **,double *b ,double *y);
int main(int argc, char * argv[])
{
  /* values
  Row0 = 0, 2, 1; b0= -8
  Row1 = 1, -2, -3; b1 = 0
  Row3 = -1, 1, 2; b2= 3

  Ans = -4, -5, 2 */
  //int **A = new {{0.0, 2.0, 1.0, -8.0}, {1.0, -2.0, -3.0, 0.0},{-1.0, 1.0, 2.0,3.0}};
  
  
  double *A[n]={NULL, NULL, NULL} ;
  double *b = NULL;
  double *y= NULL;

  for (int i=0; i < n; i++){
    A[i] = new double [n];
    if (A[i] == NULL){
       cout<< " problem with memory Allocation:";
       return -2;
    }
  }
  b = new double[n];
  if(b== NULL){
    cout<< " Null pointer assignment:";
    return -2;
  }
 
  y = new double[n];
  if(y== NULL){
    cout<< " Null pointer assignment:";
    return -2;
  }
  
  A[0][0] =0.0;A[0][1] =2.0; A[0][2] = 1.0; //A[0][3] = -8.0;
  A[1][0] =1.0;A[1][1] =-2.0; A[1][2] =-3.0; //A[1][3] = 0.0;
  A[2][0] =-1.0;A[0][1] =1.0; A[0][2] = 2.0; //A[0][3] = 3.0;
  b[0] = -8.0;
  b[1] = 0.0;
  b[2] = 3.0;

 
  
  
  GaussianElimination(A, b ,y);
  delete [] y;
  delete [ ]b;
  for(int i = 0; i < n; i++)
    delete [] A[i];
 // delete [ ] A; No need to delete A because A is not created dynamically

  
  
}
  void GaussianElimination(double **A, double *b,double *y) {
     for(int k=0; k<n-1; ++k) {
        for(int j= k+1; j <n-1; ++j)
           A[k][j] = A[k][j]/A[k][k]; 
        y[k] = b[k]/A[k][k];
        for(int i=k+1; k<n-1; i++){
           for(int j=k+1; j<n-1; j++)
              A[i][j] = A[i][j] -A[i][k] *A[k][j];
           b[i]= b[i] -A[i][k] * y[k];
           A[i][k] = 0;
        }
     }

     for(int i=0; i<n-1; i++)
        cout<<y[i];
}        

The Algorithm is given below:

/*Aginlgorithm */
/* 1. procedure Gaussian_Elimination(A, b, y)
2. begin
3. for k:= 0 to n-1 do
4. begin
5. for j:=k + 1 to n-1 do
6. A[k,j] := A[k,j] /A[k,k];
7. y[k] : = b[k]/A[k,k];
8. A[k,k] :=1;
9. for i:= k + 1 to n-1 do
10. begin
11. for j:= k+1 to n-1 do
12. A[i,j]:= A[i,j] - A[i,k] * A[k, j];
13. b[i] := b[i] -A[i, k] * y[k];
14. A[i, k] :=0;
15. endfor;*//* Line 9 */
/*16. endfor; *//*line 3*/
/*17. end

Execution and compilation

$ g++ GES1.cpp
$ ./a.out
Segmentation fault (core dumped)
$

Some body please guide me.

Zulfi.

[lastchance]

Well, it segfaults because of your line
for(int i=k+1; k<n-1; i++){
Look at the middle item.

Even when you've fixed that, your algorithm is wrong (for example, in the upper limits of the loops). I think you would do well to work through this by hand and understand the row operations involved rather than just trying to follow some prescription from Wikipedia.

And go away and look up "partial pivoting" or you will run a very real risk of dividing by 0.

PLEASE INCLUDE COMPILEABLE CODE - you are missing several lines at the top.

[zak100]

Hi,

Thanks for your reply. I have added the code required for compilation, but i cant figure out why I am getting segmentation fault.

for(int i=k+1; k<n -1; i++)

I can't figure out what's the error in the above line.

Some body please guide me.

Zulfi.
Zulfi.

[zak100]

Hi,
You are right. i am getting a divide by 0 error because initially A[0][0] = 0, so when k =0, we would get a divide by zero error.

I am just following the algorithm. This means that the algorithm is for dense simultaneous equations.

In the GaussianElimination, I have tried following:

void GaussianElimination(double **A, double *b,double *y) {
     cout<<"Inside Gaussian 1";
     for(int k=0; k<=n-1; ++k) {
        cout<<"k =" <<k<<endl;
        for(int j= k+1; j <=n-1; ++j){
           cout<<"A["<<k<<"]["<<k<<"] =" <<A[k][k] <<endl;
           A[k][j] = A[k][j]/A[k][k]; 
         }
      
       
        y[k] = b[k]/A[k][k];
        cout<<"Reached";
        for(int i=k+1; i<=n-1; i++){
           for(int j=k+1; j<=n-1; j++)
              A[i][j] = A[i][j] -A[i][k] *A[k][j];
           b[i]= b[i] -A[i][k] * y[k];
           A[i][k] = 0;
        }
     }

     for(int i=0; i<n-1; i++)
        cout<<y[i];
}   

I am getting :

Testing0Testing1Inside Gaussian 1k =0
A[0][0] =0
A[0][0] =0
Segmentation fault (core dumped)

It can't print "Reached". I think I have to change the data.

I have corrected the middle term.

Am i right?
<I think you would do well to work through this by hand and understand the row operations involved rather than just trying to follow some prescription from Wikipedia.>

If your statement is true then computer scientist can't work in other fields even if they get the right algorithm.

Zulfi.

[lastchance]

Zulfi, look up "partial pivoting". It wil save you from the initial divide-by-zero error. There should be no need to "change the data". It will also assist in reducing round-off error.

There is also nothing in your code to deal with cases where the matrix A is singular (and hence there is no solution).

With regard to your last comment, your algorithm isn't "right", so it doesn't back up your point. A computer scientist should endeavour to understand what he/she is using. In this instance, doing such a small problem by hand will clearly show you how the algorithm works.

[zak100]

Thanks for your advice.

I would try to search the algorithm of partial pivoting.

Actually I am following the book.It starts with this algorithm and then converts into a pipelined algorithm.

Thanks for your comments.

Zulfi.