Mình đang có một bài tập về lập trình một chương trình giải hệ phương trình tuyến tính Ax=b
ma trận A là vuông . vecto b .

Ma trận A sẽ được khử về dạng ma trận tam giác trên ( bước khử xuôi ) . Sau đó tìm ra các nghiệm bằng cách thế ngược lên .

Cái thế ngược mình đã làm rồi . Giờ chỉ còn cái phần chính là khử xuôi . Mình vẫn chưa giải quyết được cái trường hợp chia cho 0 và sai số trong quá trình khử xuôi

Code:
void giaiGaussXoay(vector<vector <double> > a , vector<double> &x) 
{
      unsigned i,j,k; 
      double t,m;
      unsigned n=a.size(); // a la ma tran chua he pt  ma tran [n]x[n+1]   co chua ca cot b
      
      double max;  // luu phan tu max 
      unsigned imax; // luu tru vi tri hang xoay
      
      vector<double> d; 

      for(k=0;k<n-1;k++)  // cong viec khu xuoi se hoan thanh sau n-1 buoc
      {
          cout<<"Buoc k = "<<k+1;
          cout<<endl;
       // Nhiem vu dau tien la tim phan tu a[i][k] co gia tri max
          
          max=fabs(a[k][k]);
          for(j=k+1;j<=n-1;j++)
          {
               if(fabs(a[j][k])>max)              
               {
               imax=j;// khi tim dc vi tri co a[j][k] max thi gan gia tri hang max 
               }
          }
          
          
       // doi cho tat ca cac phan tu tren hang imax cho hang k
          for(j=k;j<=n;j++) // toan bo hang j ke ca pt cua b
          {
                            t=a[imax][j];
                            a[imax][j]=a[k][j];
                            a[k][j]=t;
          }
          // den day hang thich hop da duoc de vao dung vi tri cua no
          // Ta tien hanh khu hang duoi cho hang phia tren
      
      
      // tien hanh khu ma tran ve dang tam giac tren
          for(i=k+1;i<=n-1;i++)  // bat dau tu phan tu duoi pt tru . a[k][] , a[k+1][]  ->  a[n-1][] 
          {
                       m=a[i][k]/a[k][k];
                        cout<<"m = "<<m<<"\n";
                       for(j=k;j<=n ;j++) // step 1 k=0 , i=1 ,j=k
                       {
                             
                              a[i][j] -= m*(a[k][j]); // khu pt duoi duong cheo 
                       }
                       
          }
          for(i=0;i<=n-1;i++)  // in ma tran qua moi buoc trung gian
          {
                             for(j=0;j<=n;j++)
                             {
                                              cout<<a[i][j];
                                              cout<<"\t";
                             }
                             cout<<endl;
          }

      
      }
      
      heTamGiacTren(a,x); // goi ham the nguoc giai ma tran tam giac tren
  
}