我一直在使用 C++ 程序来计算 nbyn 的行列式 方阵。它适用于较低的订单（我得到了 n=6 的精确结果）。我需要计算 21by21 矩阵的行列式，但程序需要无限的时间来计算行列式。
我正在写下我使用的代码：
# include<iostream>
# include<fstream>
# include<cmath>
using namespace std;
float determinant( float matrix2[21][21], int n2)
{
float det = 0;
float submatrix2[21][21];
if (n2 == 2)
return ((matrix2[1][1] * matrix2[0][0]) - (matrix2[0][1] * matrix2[1][0]));
else {
for (int l = 0; l < n2; l++) {
int subi = 0;
for (int m = 1; m < n2; m++) {
int subj = 0;
for (int p = 0; p < n2; p++) {
if (p == l)
continue;
submatrix2[subi][subj] = matrix2[m][p];
subj++;
}
subi++;
}
det = det + (pow(-1, l) * matrix2[0][l] * determinant( submatrix2, n2 - 1 ));
}
}
return det;
}
int main()
{
ifstream aa("matrix21.dat");
float a[21][21];
int p,q;
for(p=0;p<21;p++)
{
for(q=0;q<21;q++)
{
aa>>a[p][q];
}
}
for(p=0;p<21;p++)
{
for(q=0;q<21;q++)
{
cout<<a[p][q]<<" ";
}
cout<<endl;
}
float det;
det=determinant(a,21);
cout<<"The determinant is :"<<det<<endl;
return(0);
}
即使我尝试在输入文件中给出对角矩阵，但程序无法提供结果。
如果你尝试在未成年人中使用扩展，那么 21 个未成年人将分别触发 20 个未成年人的扩展，这将分别触发 19 个未成年人的扩展，......等等。即使有一些明智的记忆，这也是站不住脚的。
相反，您可以执行高斯消去法的步骤，将矩阵简化为上三角形式。一旦主对角线下方有零，那么行列式就是对角线元素的乘积。
请注意，部分旋转（交换行以将最大可能的元素放置在枢轴位置）通常对于数值稳定性是可取的。每行交换都会将行列式的符号更改 -1。
请注意，21x21 矩阵可能具有相当大的行列式；注意溢出。
这是一个 C++ 代码及其 21x21 矩阵的输出（在底部给出）。确保您的代码知道它期望从文件中读取多大的矩阵。
我用Python检查了它（通过作弊并使用numpy.linalg库）。
c++代码
#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include <algorithm>
using namespace std;
const double SMALL = 1.0E-30;
using matrix = vector< vector<double> >;
//======================================================================
matrix readMatrix( int N, string filename ) // It is assumed you know how big your matrix is!
{
matrix M( N, vector<double>( N ) );
ifstream in( filename );
for ( int i = 0; i < N; i++ )
{
for ( int j = 0; j < N; j++ ) in >> M[i][j];
}
return M;
}
//======================================================================
double determinant( matrix A )
{
int n = A.size();
double det = 1;
// Reduce to row-echelon form by pivoting from rows i = 0, ,,,, n - 2 (n-1 not needed)
// Better numerics by swapping rows to put largest possible element as pivot
// The determinant changes sign when you swap rows, so keep a tally!
// After reduction to RE form, the determinant is just the product of the diagonal elements
for ( int i = 0; i < n - 1; i++ )
{
// Partial pivot: find row r below with largest element in column i
int r = i;
double maxA = abs( A[i][i] );
for ( int k = i + 1; k < n; k++ )
{
double val = abs( A[k][i] );
if ( val > maxA )
{
r = k;
maxA = val;
}
}
if ( r != i )
{
for ( int j = i; j < n; j++ ) swap( A[i][j], A[r][j] );
det = -det; // Swap rows - flip sign of determinant
}
// Row operations to make upper-triangular
double pivot = A[i][i];
if ( abs( pivot ) < SMALL ) return 0.0; // Singular matrix
for ( int r = i + 1; r < n; r++ ) // Remove multiples from lower rows
{
double multiple = A[r][i] / pivot; // Multiple of row i to clear element in ith column below
for ( int j = i; j < n; j++ ) A[r][j] -= multiple * A[i][j];
}
det *= pivot; // Determinant is product of diagonal elements
}
det *= A[n-1][n-1];
return det;
}
//======================================================================
int main()
{
matrix M = readMatrix( 21, "matrix.dat" );
cout << "determinant = " << determinant( M ) << '\n';
}
//======================================================================
样本矩阵的输出（如下）：
determinant = -0.0178605
Python 库代码（待检查）：
import numpy as np
M = np.loadtxt( "matrix.dat" )
print( np.linalg.det( M ) )
样本矩阵的输出（如下）：
-0.017860521764216902
示例 21x21 矩阵（使用下面的生成器随机生成）：
0.72 0.11 0.79 0.51 0.51 0.85 0.4 0.75 0.44 0.03 0.95 0.03 0.65 0.21 0.95 0.2 0.16 0.86 0.2 0.17 0.75
0.38 0.74 0.73 0.06 0.07 0.56 0.83 0.11 0.78 0.66 0.19 0.32 0.7 0.28 0.73 0.3 0.06 0.12 0.04 0.25 0.47
0.34 0.8 0.53 0.28 0.61 0.54 0.35 0.65 0.81 0.19 0.86 0.09 0.14 0.96 0.39 0.84 0.92 0.7 0.95 0.84 0.37
0.06 0.36 0.44 0.35 0.9 0.47 0.57 0.22 0.79 0.48 0.59 0.2 0.18 0.1 0.84 0.4 0.78 0.05 0.67 0.19 0.92
0.53 0.88 0.49 0.95 0.1 0.73 0.27 0.93 0.36 0.36 0.73 0.81 0.71 0.95 0.67 0.61 0.07 0.48 0.76 0.79 0.68
0.31 0.1 0.04 0.42 0.65 0.23 0.72 0.11 0.3 0.43 0.32 0.37 0.45 0.52 0.73 0.53 0.54 0.78 0.24 0.48 0.6
0.54 0.71 0.86 0.69 0.09 0.91 0.28 0.5 0.81 0.38 0.38 0.31 0.16 0.63 0.64 0.15 0.08 0.07 0.64 0.76 0.81
0.14 0.13 0.54 0.96 0.31 0.02 0.23 0.49 0.44 0.83 0.82 0.6 0.26 0.55 0.96 0.33 0.55 0.05 0.03 0.8 0.53
0.69 0.15 0.02 0.26 0.73 0.1 0.25 0.44 0.04 0.39 0.49 0.51 0.29 0.85 0.64 0.56 0.85 0.25 0.31 0.05 0.75
0.28 0.43 0.18 0.52 0.54 0.63 0.38 0.47 0.08 0.54 0.03 0.01 0.52 0.19 0.34 0.71 0.9 0.11 0.34 0.68 0.41
0.47 0.3 0.18 0.76 0.35 0.61 0.89 0.86 0.02 0.66 0.38 0.2 0.86 0.12 0.09 0.69 0.79 0.64 0.36 0.03 0.87
0.64 0.06 0.92 0.37 0.34 0.13 0.95 0.38 0.9 0.51 0.48 0.31 0.12 0.19 0.94 0.28 0.82 0.62 0.83 0.24 0.96
0.92 0.17 0.83 0.64 0.52 0.51 0.08 0.64 0.38 0.21 0.99 0.02 0.45 0.37 0.78 0.79 0.24 0.98 0.14 0.75 0.13
0.46 0.4 0 0.84 0.74 0.32 0.78 0.84 0.41 0.82 0.34 0.39 0.71 0.92 0.17 0.07 0.23 0.39 0.11 0.93 0.98
0.45 0.38 0.45 0.69 0.55 0.33 0.46 0.69 0.82 0.82 0.07 0.98 0.03 0.03 0.78 0.5 0.71 0.01 0.6 0.66 0.44
0.29 0.25 0.7 0.1 0.08 0.53 0.6 0.4 0.84 0.56 0.67 0.87 0.19 0.58 0.54 0.95 0.63 0.3 0.45 0.59 0.91
0.8 0.2 0 0.37 0.76 0.21 0.41 0.92 0.25 0.16 0.7 0.72 0.41 0.01 0.36 0.64 0.39 0.52 0.03 0.04 0.55
0.58 0.19 0.47 0.33 0.38 0.16 0.78 0.72 0.63 0.3 0.41 0.4 0.76 0.99 0.93 0.98 0.83 0.06 0.83 0.39 0.94
0.31 0.92 0.56 0.82 0.4 0.34 0.04 0.02 0.04 0.3 0.34 0.37 0.42 0.18 0.45 0.19 0.2 0.97 0.81 0.89 0.44
0.29 0.27 0.46 0.2 0.36 0.08 0.36 0.97 0.97 0.94 0.55 0.61 0.93 0.16 0.32 0.6 0.62 0.48 0.66 0.39 0.99
0.82 0.98 0.24 0.81 0.34 0.4 0.65 0.09 0.76 0.93 0.32 0.38 0.38 0.9 0.34 0.69 0.19 0.08 0.55 0.15 0.75
简单矩阵生成器：
#include <iostream>
#include <iomanip>
#include <fstream>
#include <cstdlib>
#include <ctime>
using namespace std;
int main()
{
const int N = 21;
srand( time( 0 ) );
ofstream out( "matrix.dat" );
for ( int i = 0; i < N; i++ )
{
for ( int j = 0; j < N; j++ ) out << setw( 5 ) << 0.01 * ( rand() % 100) << " ";
out << '\n';
}
}