#include<iostream>
#include<>
#include<>
#include<>
#include<fstream>
#include <iomanip>
using namespace std;
#define ROW 25
#define COL 25
typedef struct
{
float Real;
float Image;
}Complex;
Complex add(Complex a,Complex b)
{
Complex c;
=+;
=+;
return c;
}
Complex sub(Complex a,Complex b)
{
Complex c;
=;
=;
return c;
}
Complex Mul(Complex a,Complex b)
{
Complex c;
=**;
=*+*;
return c;
}
Complex Div(Complex a,Complex b)
{
Complex c;
=(*+*)/(*+*);
=(**)/(*+*);
return c;
}
Complex matrix_det(const Complex *mat, int n)
{
int i,j,k,is,js,l,v;
Complex det;
float flag,pivot,tmp;
Complex *cpmat;
Complex tempValue;
if(mat == NULL) /* 检查输入的指针是否为空*/
{
printf("matrix pointer is NULL.\n");
}
cpmat = (Complex*)malloc(n*n*sizeof(Complex)); /* 将输入矩阵的内容拷贝一份,以免破坏*/
for(i=0; i<n*n; i++)
cpmat[i] = mat[i];
= 1.0; /* 设置行列式值初置*/
=0.0;
flag= 1.0; /* 这个变量原来记录行列式值的符号*/
for(k=0; k<n-1; k++){ /* 最多进行n-1次消去*/
pivot = 0.0; /* 选择主元*/
for(i=k; i<n; i++){
for(j=k; j<n; j++){
tmp = cpmat[i*n+j].Real*cpmat[i*n+j].Real+cpmat[i*n+j].Image*cpmat[i*n+j].Image;
if(tmp > pivot){
pivot = tmp;
is = i;
js = j;
}
}
}
if(pivot < 1e-5){ /* 如果找到的主元小于eps,则认为是0。*/
= 0.0; /*此时行列式值也是0。*/
=0.0;
return(det);
}
if(is != k){ /* 判断是否需要行交换*/
flag = -flag; /* 行交换一次,行列式值变号*/
for(j=k; j<n; j++){ /* 进行行交换*/
l = k*n + j;
v = is*n + j;
tempValue = cpmat[l];
cpmat[l] = cpmat[v];
cpmat[v] = tempValue;
}
}
if(js != k){ /* 判断是否需要列交换*/
flag = -flag; /* 列交换一次,行列式值变号*/
for(i=k; i<n; i++){ /* 进行列交换*/
l = i*n + k;
v = i*n + js;
tempValue = cpmat[v];
cpmat[v] = cpmat[l];
cpmat[l] = tempValue;
}
}
for(i=k+1; i<n; i++){ /* 进行消去*/
tempValue=Div(cpmat[i*n+k],cpmat[k*n+k]); /* 记录下此值,减少除法的次数*/
for(j=k+1; j<n; j++) { /* 消去*/
Complex aa;
aa=Mul(tempValue,cpmat[k*n+j]);
cpmat[i*n+j]=sub(cpmat[i*n+j],aa);
}
}
det = Mul(det,cpmat[k*n+k]); /*更新det的值*/
}
det = Mul(det,cpmat[k*n+k]); /* 最终更新det的值*/
=flag*;
=flag*;
free(cpmat);
return(det);
}
int main()
{
Complex *inputdata;
Complex result;
inputdata=(Complex*)malloc(sizeof(Complex)*ROW*COL);
ofstream myfile("");
srand(time(0));
for(int i=0;i<ROW*COL;i++){
inputdata[i].Real=(float)rand()/10000;
inputdata[i].Image=(float)rand()/10000;
}
//cout<<"[";
myfile<<"b=[";
for(int i=0;i<ROW;i++){
for(int j=0;j<COL;j++){
//cout<<inputdata[i*COL+j].Real<<"+"<<inputdata[i*COL+j].Image<<"i";
myfile<<inputdata[i*COL+j].Real<<"+"<<inputdata[i*COL+j].Image<<"i";
if(j<COL-1){
// cout<<",";
myfile<<",";
}
}
if(i<COL-1){
//cout<<";";
myfile<<";";
}
}
//cout<<"]"<<endl;
myfile<<"]"<<endl;;
result=matrix_det(inputdata,ROW);
//cout<<"Result:"<<endl;
//cout<<<<"+"<<<<"i"<<endl;
myfile<<setiosflags(ios::showpoint)<<<<"+"<<<<"i"<<endl;
myfile<<flush;
();
return 0;
}
计算结果与matlab计算对比,没有问题复数矩阵计算行列式
项目上需要对复矩阵的行列式计算,根据计算一般矩阵行列式的代码改成了复矩阵行列式计算。