📄 matrix.cs
字号:
* @return Matrix型，指定矩阵与value相乘之积
*/
public Matrix Multiply(double value)
{
// 构造目标矩阵
Matrix result = new Matrix(this) ; // copy ourselves
// 进行数乘
for (int i = 0 ; i < numRows ; ++i)
{
for (int j = 0 ; j < numColumns; ++j)
result.SetElement(i, j, result.GetElement(i, j) * value) ;
}
return result ;
}
/**
* 实现矩阵的乘法
*
* @param other - 与指定矩阵相乘的矩阵
* @return Matrix型，指定矩阵与other相乘之积
* @如果矩阵的行/列数不匹配，则会抛出异常
*/
public Matrix Multiply(Matrix other)
{
// 首先检查行列数是否符合要求
if (numColumns != other.GetNumRows())
throw new Exception("矩阵的行/列数不匹配。");
// ruct the object we are going to return
Matrix result = new Matrix(numRows, other.GetNumColumns());
// 矩阵乘法，即
//
// [A][B][C] [G][H] [A*G + B*I + C*K][A*H + B*J + C*L]
// [D][E][F] * [I][J] = [D*G + E*I + F*K][D*H + E*J + F*L]
// [K][L]
//
double value ;
for (int i = 0 ; i < result.GetNumRows() ; ++i)
{
for (int j = 0 ; j < other.GetNumColumns() ; ++j)
{
value = 0.0 ;
for (int k = 0 ; k < numColumns ; ++k)
{
value += GetElement(i, k) * other.GetElement(k, j) ;
}
result.SetElement(i, j, value) ;
}
}
return result ;
}
/**
* 复矩阵的乘法
*
* @param AR - 左边复矩阵的实部矩阵
* @param AI - 左边复矩阵的虚部矩阵
* @param BR - 右边复矩阵的实部矩阵
* @param BI - 右边复矩阵的虚部矩阵
* @param CR - 乘积复矩阵的实部矩阵
* @param CI - 乘积复矩阵的虚部矩阵
* @return bool型，复矩阵乘法是否成功
*/
public bool Multiply(Matrix AR, Matrix AI, Matrix BR, Matrix BI, Matrix CR, Matrix CI)
{
// 首先检查行列数是否符合要求
if (AR.GetNumColumns() != AI.GetNumColumns() ||
AR.GetNumRows() != AI.GetNumRows() ||
BR.GetNumColumns() != BI.GetNumColumns() ||
BR.GetNumRows() != BI.GetNumRows() ||
AR.GetNumColumns() != BR.GetNumRows())
return false;
// 构造乘积矩阵实部矩阵和虚部矩阵
Matrix mtxCR = new Matrix(AR.GetNumRows(), BR.GetNumColumns());
Matrix mtxCI = new Matrix(AR.GetNumRows(), BR.GetNumColumns());
// 复矩阵相乘
for (int i=0; i<AR.GetNumRows(); ++i)
{
for (int j=0; j<BR.GetNumColumns(); ++j)
{
double vr = 0;
double vi = 0;
for (int k =0; k<AR.GetNumColumns(); ++k)
{
double p = AR.GetElement(i, k) * BR.GetElement(k, j);
double q = AI.GetElement(i, k) * BI.GetElement(k, j);
double s = (AR.GetElement(i, k) + AI.GetElement(i, k)) * (BR.GetElement(k, j) + BI.GetElement(k, j));
vr += p - q;
vi += s - p - q;
}
mtxCR.SetElement(i, j, vr);
mtxCI.SetElement(i, j, vi);
}
}
CR = mtxCR;
CI = mtxCI;
return true;
}
/**
* 矩阵的转置
*
* @return Matrix型，指定矩阵转置矩阵
*/
public Matrix Transpose()
{
// 构造目标矩阵
Matrix Trans = new Matrix(numColumns, numRows);
// 转置各元素
for (int i = 0 ; i < numRows ; ++i)
{
for (int j = 0 ; j < numColumns ; ++j)
Trans.SetElement(j, i, GetElement(i, j)) ;
}
return Trans;
}
/**
* 实矩阵求逆的全选主元高斯－约当法
*
* @return bool型，求逆是否成功
*/
public bool InvertGaussJordan()
{
int i,j,k,l,u,v;
double d = 0, p = 0;
// 分配内存
int[] pnRow = new int[numColumns];
int[] pnCol = new int[numColumns];
// 消元
for (k=0; k<=numColumns-1; k++)
{
d=0.0;
for (i=k; i<=numColumns-1; i++)
{
for (j=k; j<=numColumns-1; j++)
{
l=i*numColumns+j; p=Math.Abs(elements[l]);
if (p>d)
{
d=p;
pnRow[k]=i;
pnCol[k]=j;
}
}
}
// 失败
if (d == 0.0)
{
return false;
}
if (pnRow[k] != k)
{
for (j=0; j<=numColumns-1; j++)
{
u=k*numColumns+j;
v=pnRow[k]*numColumns+j;
p=elements[u];
elements[u]=elements[v];
elements[v]=p;
}
}
if (pnCol[k] != k)
{
for (i=0; i<=numColumns-1; i++)
{
u=i*numColumns+k;
v=i*numColumns+pnCol[k];
p=elements[u];
elements[u]=elements[v];
elements[v]=p;
}
}
l=k*numColumns+k;
elements[l]=1.0/elements[l];
for (j=0; j<=numColumns-1; j++)
{
if (j != k)
{
u=k*numColumns+j;
elements[u]=elements[u]*elements[l];
}
}
for (i=0; i<=numColumns-1; i++)
{
if (i!=k)
{
for (j=0; j<=numColumns-1; j++)
{
if (j!=k)
{
u=i*numColumns+j;
elements[u]=elements[u]-elements[i*numColumns+k]*elements[k*numColumns+j];
}
}
}
}
for (i=0; i<=numColumns-1; i++)
{
if (i!=k)
{
u=i*numColumns+k;
elements[u]=-elements[u]*elements[l];
}
}
}
// 调整恢复行列次序
for (k=numColumns-1; k>=0; k--)
{
if (pnCol[k]!=k)
{
for (j=0; j<=numColumns-1; j++)
{
u=k*numColumns+j;
v=pnCol[k]*numColumns+j;
p=elements[u];
elements[u]=elements[v];
elements[v]=p;
}
}
if (pnRow[k]!=k)
{
for (i=0; i<=numColumns-1; i++)
{
u=i*numColumns+k;
v=i*numColumns+pnRow[k];
p=elements[u];
elements[u]=elements[v];
elements[v]=p;
}
}
}
// 成功返回
return true;
}
/**
* 复矩阵求逆的全选主元高斯－约当法
*
* @param mtxImag - 复矩阵的虚部矩阵，当前矩阵为复矩阵的实部
* @return bool型，求逆是否成功
*/
public bool InvertGaussJordan(Matrix mtxImag)
{
int i,j,k,l,u,v,w;
double p,q,s,t,d,b;
// 分配内存
int[] pnRow = new int[numColumns];
int[] pnCol = new int[numColumns];
// 消元
for (k=0; k<=numColumns-1; k++)
{
d=0.0;
for (i=k; i<=numColumns-1; i++)
{
for (j=k; j<=numColumns-1; j++)
{
u=i*numColumns+j;
p=elements[u]*elements[u]+mtxImag.elements[u]*mtxImag.elements[u];
if (p>d)
{
d=p;
pnRow[k]=i;
pnCol[k]=j;
}
}
}
// 失败
if (d == 0.0)
{
return false;
}
if (pnRow[k]!=k)
{
for (j=0; j<=numColumns-1; j++)
{
u=k*numColumns+j;
v=pnRow[k]*numColumns+j;
t=elements[u];
elements[u]=elements[v];
elements[v]=t;
t=mtxImag.elements[u];
mtxImag.elements[u]=mtxImag.elements[v];
mtxImag.elements[v]=t;
}
}
if (pnCol[k]!=k)
{
for (i=0; i<=numColumns-1; i++)
{
u=i*numColumns+k;
v=i*numColumns+pnCol[k];
t=elements[u];
elements[u]=elements[v];
elements[v]=t;
t=mtxImag.elements[u];
mtxImag.elements[u]=mtxImag.elements[v];
mtxImag.elements[v]=t;
}
}
l=k*numColumns+k;
elements[l]=elements[l]/d; mtxImag.elements[l]=-mtxImag.elements[l]/d;
for (j=0; j<=numColumns-1; j++)
{
if (j!=k)
{
u=k*numColumns+j;
p=elements[u]*elements[l];
q=mtxImag.elements[u]*mtxImag.elements[l];
s=(elements[u]+mtxImag.elements[u])*(elements[l]+mtxImag.elements[l]);
elements[u]=p-q;
mtxImag.elements[u]=s-p-q;
}
}
for (i=0; i<=numColumns-1; i++)
{
if (i!=k)
{
v=i*numColumns+k;
for (j=0; j<=numColumns-1; j++)
{
if (j!=k)
{
u=k*numColumns+j;
w=i*numColumns+j;
p=elements[u]*elements[v];
q=mtxImag.elements[u]*mtxImag.elements[v];
s=(elements[u]+mtxImag.elements[u])*(elements[v]+mtxImag.elements[v]);
t=p-q;
b=s-p-q;
elements[w]=elements[w]-t;
mtxImag.elements[w]=mtxImag.elements[w]-b;
}
}
}
}
for (i=0; i<=numColumns-1; i++)
{
if (i!=k)
{
u=i*numColumns+k;
p=elements[u]*elements[l];
q=mtxImag.elements[u]*mtxImag.elements[l];
s=(elements[u]+mtxImag.elements[u])*(elements[l]+mtxImag.elements[l]);
elements[u]=q-p;
mtxImag.elements[u]=p+q-s;
}
}
}
// 调整恢复行列次序
for (k=numColumns-1; k>=0; k--)
{
if (pnCol[k]!=k)
{
for (j=0; j<=numColumns-1; j++)
{
u=k*numColumns+j;
v=pnCol[k]*numColumns+j;
t=elements[u];
elements[u]=elements[v];
elements[v]=t;
t=mtxImag.elements[u];
mtxImag.elements[u]=mtxImag.elements[v];
mtxImag.elements[v]=t;
}
}
if (pnRow[k]!=k)
{
for (i=0; i<=numColumns-1; i++)
{
u=i*numColumns+k;
v=i*numColumns+pnRow[k];
t=elements[u];
elements[u]=elements[v];
elements[v]=t;
t=mtxImag.elements[u];
mtxImag.elements[u]=mtxImag.elements[v];
mtxImag.elements[v]=t;
}
}
}
// 成功返回
return true;
}
/**
* 对称正定矩阵的求逆
*
* @return bool型，求逆是否成功
*/
public bool InvertSsgj()
{
int i, j ,k, m;
double w, g;
// 临时内存
double[] pTmp = new double[numColumns];
// 逐列处理
for (k=0; k<=numColumns-1; k++)
{
w=elements[0];
if (w == 0.0)
{
return false;
}
m=numColumns-k-1;
for (i=1; i<=numColumns-1; i++)
{
g=elements[i*numColumns];
pTmp[i]=g/w;
if (i<=m)
pTmp[i]=-pTmp[i];
for (j=1; j<=i; j++)
elements[(i-1)*numColumns+j-1]=elements[i*numColumns+j]+g*pTmp[j];
}
elements[numColumns*numColumns-1]=1.0/w;
for (i=1; i<=numColumns-1; i++)
elements[(numColumns-1)*numColumns+i-1]=pTmp[i];
}
// 行列调整
for (i=0; i<=numColumns-2; i++)
for (j=i+1; j<=numColumns-1; j++)
elements[i*numColumns+j]=elements[j*numColumns+i];
return true;
}
/**
* 托伯利兹矩阵求逆的埃兰特方法
*
* @return bool型，求逆是否成功
*/
public bool InvertTrench()
{
int i,j,k;
double a,s;
// 上三角元素
double[] t = new double[numColumns];
// 下三角元素
double[] tt = new double[numColumns];
// 上、下三角元素赋值
for (i=0; i<numColumns; ++i)
{
t[i] = GetElement(0, i);
tt[i] = GetElement(i, 0);
}
// 临时缓冲区
double[] c = new double[numColumns];
double[] r = new double[numColumns];
double[] p = new double[numColumns];
// 非Toeplitz矩阵，返回
if (t[0] == 0.0)
{
return false;
}
a=t[0];
c[0]=tt[1]/t[0];
r[0]=t[1]/t[0];
for (k=0; k<=numColumns-3; k++)
{
s=0.0;
for (j=1; j<=k+1; j++)
s=s+c[k+1-j]*tt[j];
s=(s-tt[k+2])/a;
for (i=0; i<=k; i++)
p[i]=c[i]+s*r[k-i];
c[k+1]=-s;
s=0.0;
for (j=1; j<=k+1; j++)
s=s+r[k+1-j]*t[j];
s=(s-t[k+2])/a;
for (i=0; i<=k; i++)
{
r[i]=r[i]+s*c[k-i];
c[k-i]=p[k-i];
}
r[k+1]=-s;
a=0.0;
for (j=1; j<=k+2; j++)
a=a+t[j]*c[j-1];
a=t[0]-a;
// 求解失败
if (a == 0.0)
{
return false;
}
}
elements[0]=1.0/a;
for (i=0; i<=numColumns-2; i++)
{
k=i+1;
j=(i+1)*numColumns;
elements[k]=-r[i]/a;
elements[j]=-c[i]/a;
}
for (i=0; i<=numColumns-2; i++)
{
for (j=0; j<=numColumns-2; j++)
{
k=(i+1)*numColumns+j+1;
elements[k]=elements[i*numColumns+j]-c[i]*elements[j+1];
elements[k]=elements[k]+c[numColumns-j-2]*elements[numColumns-i-1];
}
}
return true;
}
/**
* 求行列式值的全选主元高斯消去法
*
* @return double型，行列式的值
*/
public double ComputeDetGauss()
{
int i,j,k,nis = 0,js = 0,l,u,v;
double f,det,q,d;
// 初值
f=1.0;
det=1.0;
// 消元
for (k=0; k<=numColumns-2; k++)
{
q=0.0;
for (i=k; i<=numColumns-1; i++)
{
for (j=k; j<=numColumns-1; j++)
{
l=i*numColumns+j;
d=Math.Abs(elements[l]);
if (d>q)
{
q=d;
nis=i;
js=j;
}
}
}
if (q == 0.0)
{
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -