2015-12-18 33 views
1

我需要编写一个从逗号分隔值文件中读取矩阵的程序,然后使用高斯消元来计算逆并将该反写写入新文件。使用高斯消除的文本文件的逆矩阵。 C++

读它是好的,正如它写回来一样。我想我理解高斯消元是如何工作的,并且能够通过使用数组来完成这一任务。

#include <iostream> 
#include <fstream> 
#include <iomanip> 
#include <vector> 
#include <string> 

using namespace std; 

// main program starts here 
int main() { 
// create a vector of a vector to store original matrix 
// and matrices after each calculation 
vector<vector<double>> data; 
int n_com = 0; 

// try to read input file 
ifstream readFile("test_data.txt"); 

if (readFile.is_open()){ 
    while (!readFile.eof()){ 
     int i; 

     // declare temporary line vector and line string 
     vector<double> vline; 
     string aLine; 

     // assign line of file to line string 
     getline(readFile, aLine); 

     // count number of commas in line string 
     n_com = count(aLine.begin(), aLine.end(), ','); 

     // define integer for start of each element 
     int start = 0; 

     // loop over all but final element 
     for (i=0; i < n_com; i++) { 

      // declare and find position of next comma 
      int comma_pos; 
      comma_pos = aLine.find(',', start); 

      // declare string for element and assign substring to it 
      string elems; 
      elems = aLine.substr(start, comma_pos - start); 

      // convert string to double 
      double elemd = atof(elems.c_str()); 

      // push back double to temporary vector 
      vline.push_back(elemd); 

      // redefine start for next iteration 
      start = comma_pos + 1; 
     } 
     // assign final element to string 
     string final_elems = aLine.substr(start, aLine.length() - start); 

     // convert final element to double and push back to vector 
     double final_elemd = atof(final_elems.c_str()); 
     vline.push_back(final_elemd); 

     // push back line vector to data vector 
     data.push_back(vline); 
    } 
} 
else { 
    // print error if unable to open file 
    printf("Error unable to open input file!\n"); 
    // exit program 
    exit(1); 
} 
// close input file 
readFile.close(); 

所以这就是我读过的原始矩阵。

这是我的代码有高斯消元法对于数组

// calculate width and length of original data (no. of rows and columns) 
int length = data.size(); 
int width = n_com + 1; 

// create new file to write to 
ofstream writeFile ("tranpose.txt"); 

// check outfile is open 
if (writeFile.is_open()){ 

    // declare indices 
    // width (columns) 
    int i; 
    // length (rows)0 
    int j; 
    // k 
    int k; 

    // declare a float? 
    float data[10][10] = {0},d; 

    // identity matrix 
    for (i=1; i <= length; i++){ 
     for (j=1; j <= 2 * length; j++){ 
      if (j == (i + length)){ 
       data[i][j] = 1; 
      } 
     } 
    } 

    // partial pivoting 
    for (i=length; i > 1; i--){ 
     if (data[i-1][1] < data[i][1]){ 
      for(j=1;j <= length * 2; j++){ 
       d = data[i][j]; 
       data[i][j] = data[i-1][j]; 
       data[i-1][j] = d; 
      } 
     } 
    } 
    cout<<"Augmented Matrix: "<<endl; 
    for (i=1; i <= length; i++){ 
     for (j=1;j <= length * 2; j++){ 
      cout<<data[i][j]<<" "; 
     } 
     cout<<endl; 
    } 

    // reducing to diagonal matrix 
    for (i=1; i <= length; i++){ 
     for (j=1; j <= length * 2; j++){ 
      if (j != i){ 
       d = data[j][i]/data [i][i]; 
       for (k=1; k<= length * 2; k++){ 
        data[j][k] = data[j][k] - (data[i][k] * d); 
       } 
      } 
     } 
    } 

    // reducing to unit matrix 
    for (i=1; i <= length; i++){ 
     d = data[i][i]; 
     for (j=1; j <= length * 2; j++){ 
      data[i][j] = data[i][j]/d; 
     } 
    } 
    // print inverse matrix in console 
    cout<<"Inverse Matrix "<<endl; 
    for (i=1; i <= length; i++){ 
     for (j = length + 1; j <= length * 2; j++){ 
      cout<<data[i][j]<<" "; 
     } 
     cout<<endl; 
    } 
    // loop over all rows 
    for (i=1; i <= length; i++){ 
     // loop over all columns 
     for (j = length + 1; j <= length * 2; j++){ 
      // print data in transposed positions excluding last value 
      // i.e. [j][i] instead of [i][j] 
      writeFile << setw(4) << fixed << setprecision(2) << data[i][j] << ","; 
     } 
     // print onto new line 
     int i = length - 1; 
     writeFile << setw(4) << fixed << setprecision(2) << data[i][j] << "\n"; 

    } 

    // close written file 
    writeFile.close(); 

我怎么可以这样写,这样它会使用从我的向量,或矩阵文件存储的数据,而不是一个正常的数组号码?

+0

不知道你是想问什么。矢量矢量在元素访问方面的行为与2D数组完全相同。 – Lucien

+0

感谢您的回复。这也是我的想法,但是这个代码只需要一个0的矩阵作为输入,然后反转,得到一个“nan”的矩阵。 我只想知道这段代码中的问题,以及如何正确编写高斯消除,并将其应用于我的输入文本文件。 –

+0

试一下注释掉这行'float data [10] [10] = {0},d;' – Lucien

回答

0

我改进了整个事情,因为你的算法看起来并不适合我。

using Matrix = std::vector<std::vector<double>>; 

Matrix inverse(Matrix mat) 
{ 
    // Use Gaussian elimination 
    // Using two matrixs instead of one agumented 
    // to improve peformance 

    auto height = mat.size(); 
    auto width = mat[0].size(); 

    // Create an identity matrix 
    Matrix result(height, Matrix::value_type(width)); 
    for (auto i = 0;i < width;++i) { 
     result[i][i] = 1; 
    } 
    cout << "Augmented Matrix: " << endl; 
    printTwo(mat, result); 

    // reduce to Echelon form 
    for (auto j = 0;j < width;++j) { 
     // partial pivoting 
     auto maxRow = j; 
     for (auto i = j;i < height;++i) { 
      maxRow = mat[i][j]>mat[maxRow][j] ? i : maxRow; 
     } 
     mat[j].swap(mat[maxRow]); 
     result[j].swap(result[maxRow]); 

     cout << "pivotted Matrix: " << endl; 
     printTwo(mat, result); 

     // Reduce row by row 
     auto pivot = mat[j][j]; 
     auto& row1L = mat[j]; 
     auto& row1R = result[j]; 
     for (auto i = j + 1;i < height;++i) { 
      auto& row2L = mat[i]; 
      auto& row2R = result[i]; 
      auto temp = row2L[j]; 
      for (auto k = 0;k < width;++k) { 
       row2L[k] -= temp/pivot*row1L[k]; 
       row2R[k] -= temp/pivot*row1R[k]; 
      } 
     } 
     // Make diaganal elements to 1 
     for (auto k = 0;k < width;++k) { 
      row1L[k] /= pivot; 
      row1R[k] /= pivot; 
     } 
     cout << "reduced Matrix: " << endl; 
     printTwo(mat, result); 
    } 

    //back subsitution 
    for (auto j = width - 1;;--j) { 
     auto& row1L = mat[j]; 
     auto& row1R = result[j]; 
     for (auto i = 0;i < j;++i) { 
      auto& row2L = mat[i]; 
      auto& row2R = result[i]; 
      auto temp = row2L[j]; 
      for (auto k = 0;k < width;++k) { 
       row2L[k] -= temp*row1L[k]; 
       row2R[k] -= temp*row1R[k]; 
      } 
     } 
     cout << "subsituted Matrix: " << endl; 
     printTwo(mat, result); 
     if (j == 0) break; 
    } 

    return result; 
} 

而这里的代码中使用的帮手:

void printTwo(const Matrix& lhs, const Matrix& rhs) 
{ 
    for (auto i = 0;i < lhs.size();++i) { 
     for (auto elm : lhs[i]) { 
      cout << setw(4) << elm << ' '; 
     } 

     for (auto elm : rhs[i]) { 
      cout << setw(4) << elm << ' '; 
     } 
     cout << endl; 
    } 
}