C#，码海拾贝（14）——实矩阵与复数矩阵“求逆”的高斯·约当算法C#源程序

1、矩阵求逆

矩阵求逆，即求矩阵的逆矩阵。矩阵是线性代数的主要内容，很多实际问题用矩阵的思想去解既简单又快捷。逆矩阵又是矩阵理论的很重要的内容，逆矩阵的求法自然也就成为线性代数研究的主要内容之一。
设A是数域上的一个n阶方阵，若在相同数域上存在另一个n阶矩B，使得： AB=BA=E。 则我们称B是A的逆矩阵，而A则被称为可逆矩阵。其中，E为单位矩阵。
典型的矩阵求逆方法有：利用定义求逆矩阵、初等变换法、伴随阵法、恒等变形法等。

2、矩阵求逆的计算过程简述

求元素为具体数字的矩阵的逆矩阵，常用初等变换法‘如果A可逆，则A’可通过初等变换，化为单位矩阵 I ，即存在初等矩阵使 P1,P2,P3,...Ps：
（1）P1 x P2 x P3 x ...x Ps x A = I；
（2）用 A-1 右乘上式两端，得：P1 x P2 x P3 x ...x Ps x I = A-1；
比较（1）、（2）两式，可以看到当A通过初等变换化为单位处阵的同时，对单位矩阵I作同样的初等变换，就化为A的逆矩阵 。
这就是求逆矩阵的初等行变换法，它是实际应用中比较简单的一种方法。需要注意的是，在作初等变换时只允许作行初等变换。同样，只用列初等变换也可以求逆矩阵。

用此方法求逆矩阵，对于小型矩阵，特别是二阶方阵求逆既方便、快捷，又有规律可循。因为二阶可逆矩阵的伴随矩阵，只需要将主对角线元素的位置互换，次对角线的元素变号即可。
若可逆矩阵是二阶或二阶以上矩阵，在求逆矩阵的过程中，需要求9个或9个以上代数余子式，还要计算一个三阶或三阶以上行列式，工作量大且中途难免出现符号及计算的差错。对于求出的逆矩阵是否正确，一般要通过 来检验。一旦发现错误，必须对每一计算逐一排查。
 

using System;namespace Zhou.CSharp.Algorithm
{/// 
/// 矩阵类/// 作者：周长发/// 改进：深度混淆/// https://blog.csdn.net/beijinghorn/// 
public partial class Matrix{/// 
/// 实矩阵求逆的全选主元高斯－约当法/// 
/// 源矩阵/// 求逆是否成功public static bool InvertGaussJordan(Matrix src){// 分配内存int[] pnRow = new int[src.Columns];int[] pnCol = new int[src.Columns];// 消元for (int k = 0; k <= src.Columns - 1; k++){double d = 0.0;for (int i = k; i <= src.Columns - 1; i++){for (int j = k; j <= src.Columns - 1; j++){int w = i * src.Columns + j;double p = Math.Abs(src[w]);if (p > d){d = p;pnRow[k] = i;pnCol[k] = j;}}}// 失败if (Math.Abs(d) < float.Epsilon){return false;}if (pnRow[k] != k){for (int j = 0; j <= src.Columns - 1; j++){int u = k * src.Columns + j;int v = pnRow[k] * src.Columns + j;double p = src[u];src[u] = src[v];src[v] = p;}}if (pnCol[k] != k){for (int i = 0; i <= src.Columns - 1; i++){int u = i * src.Columns + k;int v = i * src.Columns + pnCol[k];double p = src[u];src[u] = src[v];src[v] = p;}}{int w = k * src.Columns + k;src[w] = 1.0 / src[w];for (int j = 0; j <= src.Columns - 1; j++){if (j != k){int u = k * src.Columns + j;src[u] = src[u] * src[w];}}for (int i = 0; i <= src.Columns - 1; i++){if (i != k){for (int j = 0; j <= src.Columns - 1; j++){if (j != k){int u = i * src.Columns + j;src[u] = src[u] - src[i * src.Columns + k] * src[k * src.Columns + j];}}}}for (int i = 0; i <= src.Columns - 1; i++){if (i != k){int u = i * src.Columns + k;src[u] = -src[u] * src[w];}}}}// 调整恢复行列次序for (int k = src.Columns - 1; k >= 0; k--){if (pnCol[k] != k){for (int j = 0; j <= src.Columns - 1; j++){int u = k * src.Columns + j;int v = pnCol[k] * src.Columns + j;double p = src[u];src[u] = src[v];src[v] = p;}}if (pnRow[k] != k){for (int i = 0; i <= src.Columns - 1; i++){int u = i * src.Columns + k;int v = i * src.Columns + pnRow[k];double p = src[u];src[u] = src[v];src[v] = p;}}}// 成功返回return true;}/// 
/// 复矩阵求逆的全选主元高斯－约当法/// 
/// 源矩阵/// 复矩阵的虚部矩阵，当前矩阵为复矩阵的实部/// 求逆是否成功public static bool InvertGaussJordan(Matrix src, Matrix mtxImag){int i, j, k, z, u, v, w;double p, q, s, t, d, b;// 分配内存int[] pnRow = new int[src.Columns];int[] pnCol = new int[src.Columns];// 消元for (k = 0; k <= src.Columns - 1; k++){d = 0.0;for (i = k; i <= src.Columns - 1; i++){for (j = k; j <= src.Columns - 1; j++){u = i * src.Columns + j;p = src[u] * src[u] + mtxImag[u] * mtxImag[u];if (p > d){d = p;pnRow[k] = i;pnCol[k] = j;}}}// 失败if (Math.Abs(d) < float.Epsilon){return false;}if (pnRow[k] != k){for (j = 0; j <= src.Columns - 1; j++){u = k * src.Columns + j;v = pnRow[k] * src.Columns + j;t = src[u];src[u] = src[v];src[v] = t;t = mtxImag[u];mtxImag[u] = mtxImag[v];mtxImag[v] = t;}}if (pnCol[k] != k){for (i = 0; i <= src.Columns - 1; i++){u = i * src.Columns + k;v = i * src.Columns + pnCol[k];t = src[u];src[u] = src[v];src[v] = t;t = mtxImag[u];mtxImag[u] = mtxImag[v];mtxImag[v] = t;}}z = k * src.Columns + k;src[z] = src[z] / d;mtxImag[z] = -mtxImag[z] / d;for (j = 0; j <= src.Columns - 1; j++){if (j != k){u = k * src.Columns + j;p = src[u] * src[z];q = mtxImag[u] * mtxImag[z];s = (src[u] + mtxImag[u]) * (src[z] + mtxImag[z]);src[u] = p - q;mtxImag[u] = s - p - q;}}for (i = 0; i <= src.Columns - 1; i++){if (i != k){v = i * src.Columns + k;for (j = 0; j <= src.Columns - 1; j++){if (j != k){u = k * src.Columns + j;w = i * src.Columns + j;p = src[u] * src[v];q = mtxImag[u] * mtxImag[v];s = (src[u] + mtxImag[u]) * (src[v] + mtxImag[v]);t = p - q;b = s - p - q;src[w] = src[w] - t;mtxImag[w] = mtxImag[w] - b;}}}}for (i = 0; i <= src.Columns - 1; i++){if (i != k){u = i * src.Columns + k;p = src[u] * src[z];q = mtxImag[u] * mtxImag[z];s = (src[u] + mtxImag[u]) * (src[z] + mtxImag[z]);src[u] = q - p;mtxImag[u] = p + q - s;}}}// 调整恢复行列次序for (k = src.Columns - 1; k >= 0; k--){if (pnCol[k] != k){for (j = 0; j <= src.Columns - 1; j++){u = k * src.Columns + j;v = pnCol[k] * src.Columns + j;t = src[u];src[u] = src[v];src[v] = t;t = mtxImag[u];mtxImag[u] = mtxImag[v];mtxImag[v] = t;}}if (pnRow[k] != k){for (i = 0; i <= src.Columns - 1; i++){u = i * src.Columns + k;v = i * src.Columns + pnRow[k];t = src[u];src[u] = src[v];src[v] = t;t = mtxImag[u];mtxImag[u] = mtxImag[v];mtxImag[v] = t;}}}// 成功返回return true;}}
}

本文来自互联网用户投稿，文章观点仅代表作者本人，不代表本站立场，不承担相关法律责任。如若转载，请注明出处。 如若内容造成侵权/违法违规/事实不符，请点击【内容举报】进行投诉反馈！