C#矩阵的运算代码
#region 矩阵运算
///
/// 矩阵对应行列式的值
///
///
///
private double MatrixValue(double[,] MatrixList)
{
int Level = MatrixList.GetLength(1);
double[,] dMatrix = new double[Level, Level];
for (int i = 0; i
{
for (int j = 0; j
{
dMatrix[i, j] = MatrixList[i, j];
}
}
int sign = 1;
for (int i = 0, j = 0; i
{
//判断改行dMatrix[i, j]是否为0,若是,则寻找i 后的行(m,m>i,切dMatrix[m, j]!=0)进行交换 if (dMatrix[i, j] == 0)
{
if (i == Level - 1)
{
return 0;
}
int m = i + 1;
//获取一个dMatrix[m, j]不为为0的行
for (; dMatrix[m, j] == 0; m++)
{
if (m == Level - 1)
{
return 0;
}
}
//判断是否达到矩阵的最大行,若是,则返回0
//把i 行和m 行调换
double temp;
for (int n = j; n
{
temp = dMatrix[i, n];
dMatrix[i, n] = dMatrix[m, n];
dMatrix[m, n] = temp;
}
sign *= (-1);
}
//把当前行以后的行所对应的列变成0
double tmp;
for (int s = Level - 1; s > i; s--)
{
tmp = dMatrix[s, j];
//j行后面的所有行
for (int t = j; t
{
dMatrix[s, t] -= dMatrix[i, t] * (tmp / dMatrix[i, j]);
}
}
}
double result = 1;
for (int i = 0; i
{
if (dMatrix[i, i] != 0)
{
result *= dMatrix[i, i];
}
else
{
return 0;
}
}
return sign * result;
}
///
/// 矩阵减法
///
///
///
///
///
private double[] SubMatrix(double[] A1, double[] A2)
{
//判断矩阵的长短是否一致
int a1 = A1.GetLength(0);
int a2 = A2.GetLength(0);
if (a1 != a2)
{
return null;
}
//矩阵相减
double[] B = new double[a1];
for (int i = 0; i
{
B[i] = A1[i] - A2[i];
}
return B;
}
///
/// 矩阵乘法
///
///
///
///
private double[,] MultiplyMatrix(double[,] firstMatrix, double[,] secondMatrix)
{
double[,] resultMatrix = new double[firstMatrix.GetLength(0), secondMatrix.GetLength(1)]; //判断相乘矩阵是否合法,即第一个矩阵的列要等于第二个矩阵的行
if (firstMatrix.GetLength(1) != secondMatrix.GetLength(0))
{
return null;
}
//求结果矩阵
for (int rowIndex = 0; rowIndex
{
for (int colIndex = 0; colIndex
{
//初始化结果矩阵的元素
resultMatrix[rowIndex, colIndex] = 0;
for (int i = 0; i
{
//求结果矩阵的元素值
resultMatrix[rowIndex, colIndex] += firstMatrix[rowIndex, i] * secondMatrix[i, colIndex]; }
}
}
return resultMatrix;
}
///
/// 求逆矩阵
///
///
///
private double[,] Athwart(double[,] dMatrix)
{
//获取矩阵的行数
int Level = dMatrix.GetLength(1);
double[,] dReverseMatrix = new double[Level, 2 * Level]; //初始化矩阵Level×(2*Level)
for (int i = 0; i
{
for (int j = 0; j
{
if (j
{
dReverseMatrix[i, j] = dMatrix[i, j];
}
else
{
if (j - Level == i)
{
dReverseMatrix[i, j] = 1;
}
else
{
dReverseMatrix[i, j] = 0;
}
}
}
}
for (int i = 0, j = 0; i
if (dReverseMatrix[i, j] == 0)
{
if (i == Level - 1)
{
return null;
}
int m = i + 1;
for (; dMatrix[m, j] == 0; m++)
{
if (m == Level - 1)
{
return null;
}
if (m == Level)
{
return null;
}
else
{
//把i 行和m 行相加
for (int n = j; n
{
dReverseMatrix[i, n] += dReverseMatrix[m, n]; }
}
}
double temp = dReverseMatrix[i, j];
if (temp != 1)
{
//把i 行数据,变成以1开始的一行数据
for (int n = j; n
{
if (dReverseMatrix[i, n] != 0)
{
dReverseMatrix[i, n] /= temp;
}
}
}
//把i 行后的所有行的j 列变成0
for (int s = Level - 1; s > i; s--)
{
temp = dReverseMatrix[s, j];
for (int t = j; t
{
dReverseMatrix[s, t] -= (dReverseMatrix[i, t] * temp); }
}
}
//把矩阵Level×(2*Level)前Level×Level 转变为单位矩阵 for (int i = Level - 2; i >= 0; i--)
{
for (int j = i + 1; j
{
if (dReverseMatrix[i, j] != 0)
double tmp = dReverseMatrix[i, j];
for (int n = j; n
{
dReverseMatrix[i, n] -= (tmp * dReverseMatrix[j, n]); }
}
}
}
//返回逆矩阵
double[,] dReturn = new double[Level, Level]; for (int i = 0; i
{
for (int j = 0; j
{
dReturn[i, j] = dReverseMatrix[i, j + Level]; }
}
return dReturn;
}
#endregion
第二份: 1using System; 2using System.IO;
3using System.Diagnostics;
4
5
6namespace Adjust
7
8{ ///
9 /// Matrix 的摘要说明。
10 /// 实现矩阵的基本运算 11 ///
12 public class Matrix 1314
15 //构造方阵
16 public Matrix(int row)
17 { { 18 m_data = new double [row,row]; 19
20 }
21 public Matrix(int row, int col) 22 {
23 m_data = new double [row,col]; 24 }
25 //复制构造函数
26 public Matrix(Matrix m)
27 {
28 int row = m.Row;
29 int col = m.Col;
30 m_data = new double [row,col]; 31
33 for (int j=0;j
36 }
37
38 /*
39 //分配方阵的大小
40 //对于已含有内存的矩阵,将清空数据 41 public void SetSize(int row) 42 {
43 m_data = new double[row,row]; 44 }
45
46
47 //分配矩阵的大小
48 //对于已含有内存的矩阵,将清空数据 49 public void SetSize(int row,int col) 50 {
51 m_data = new double[row,col]; 52 }
53 */
55 //unit matrix:设为单位阵 56 public void SetUnit() 57 { 58 for (int i=0;i
61 }
62
63 //设置元素值
64 public void SetValue(double d) 65 {
66 for (int i=0;i
69 }
70
71 // Value extraction :返中行数 72 public int Row
73 {
74 get
75 {
77 return m_data.GetLength(0); 78 } 79 } 80 81 //返回列数 82 public int Col
83 {
84 get
85 {
86 return m_data.GetLength(1); 87 }
88 }
89
90 //重载索引
91 //存取数据成员
92 public double this [int row, int col] 93 {
94 get
95 {
96 return m_data[row,col]; 97 }
98 set 99 { 100 m_data[row,col] = value; 101 } 102 } 103
104 //primary change
105 // 初等变换 对调两行:rirj
106 public Matrix Exchange(int i, int j)
107 {
108 double temp;
109
110 for (int k=0;k
111 {
112 temp = m_data[i,k];
113 m_data[i,k] = m_data[j,k];
114 m_data[j,k] = temp;
115 }
116 return this ;
117 }
118
119
121 Matrix Multiple(int index, double mul) 122 { 123 for (int j=0;j
124 {
125 m_data[index,j] *= mul;
126 }
127 return this ;
128 }
129
130
131 //初等变换 第src 行乘以mul 加到第index 行
132 Matrix MultipleAdd(int index, int src, double mul)
133 {
134 for (int j=0;j
135 {
136 m_data[index,j] += m_data[src,j]*mul;
137 }
138
139 return this ;
140 }
141
143 public Matrix Transpose() 144 { 145 Matrix ret = new Matrix(Col,Row); 146
147 for (int i=0;i
148 for (int j=0;j
149 {
150 ret[j,i] = m_data[i,j];
151 }
152 return ret;
153 }
154
155 //binary addition 矩阵加
156 public static Matrix operator + (Matrix lhs,Matrix rhs)
157 {
158 if (lhs.Row != rhs.Row) //异常
159 {
160 System.Exception e = new Exception("相加的两个矩阵的行数不等"); 161 throw e;
162 }
163 if (lhs.Col != rhs.Col) //异常
164 { 165 System.Exception e = new Exception("相加的两个矩阵的列数不等"); 166 throw e; 167 } 168
169 int row = lhs.Row;
170 int col = lhs.Col;
171 Matrix ret=new Matrix(row,col);
172
173 for (int i=0;i
174 for (int j=0;j
175 {
176 double d = lhs[i,j] + rhs[i,j];
177 ret[i,j] = d;
178 }
179 return ret;
180
181 }
182
183 //binary subtraction 矩阵减
184 public static Matrix operator - (Matrix lhs,Matrix rhs)
185 {
186 if (lhs.Row != rhs.Row) //异常 187 { 188 System.Exception e = new Exception("相减的两个矩阵的行数不等"); 189 throw e;
190 }
191 if (lhs.Col != rhs.Col) //异常
192 {
193 System.Exception e = new Exception("相减的两个矩阵的列数不等"); 194 throw e;
195 }
196
197 int row = lhs.Row;
198 int col = lhs.Col;
199 Matrix ret=new Matrix(row,col);
200
201 for (int i=0;i
202 for (int j=0;j
203 {
204 double d = lhs[i,j] - rhs[i,j];
205 ret[i,j] = d;
206 }
207 return ret;
209 210 211 //binary multiple 矩阵乘 212 public static Matrix operator * (Matrix lhs,Matrix rhs)
213 {
214 if (lhs.Col != rhs.Row) //异常
215 {
216 System.Exception e = new Exception("相乘的两个矩阵的行列数不匹配"); 217 throw e;
218 }
219
220 Matrix ret = new Matrix (lhs.Row,rhs.Col);
221 double temp;
222 for (int i=0;i
223 {
224 for (int j=0;j
225 {
226 temp = 0;
227 for (int k=0;k
228 {
229 temp += lhs[i,k] * rhs[k,j];
231 ret[i,j] = temp; 232 } 233 } 234 235 return ret;
236 }
237
238
239 //binary division 矩阵除
240 public static Matrix operator / (Matrix lhs,Matrix rhs)
241 {
242 return lhs * rhs.Inverse();
243 }
244
245 //unary addition 单目加
246 public static Matrix operator + (Matrix m)
247 {
248 Matrix ret = new Matrix(m);
249 return ret;
250 }
251
252 //unary subtraction 单目减 253 public static Matrix operator - (Matrix m) 254 { 255 Matrix ret = new Matrix(m);
256 for (int i=0;i
257 for (int j= 0;j
258 {
259 ret[i,j] = -ret[i,j];
260 }
261
262 return ret;
263 }
264
265 //number multiple 数乘
266 public static Matrix operator * (double d,Matrix m)
267 {
268 Matrix ret = new Matrix(m);
269 for (int i=0;i
270 for (int j=0;j
271 ret[i,j] *= d;
272
273 return ret;
275 276 //number division 数除 277 public static Matrix operator / (double d,Matrix m)
278 {
279 return d*m.Inverse();
280 }
281
282 //功能:返回列主元素的行号
283 //参数:row 为开始查找的行号
284 //说明:在行号[row,Col)范围内查找第row 列中绝对值最大的元素,返回所在行号 285 int Pivot(int row)
286 {
287 int index=row;
288
289 for (int i=row+1;i
290 {
291 if (m_data[i,row] > m_data[index,row])
292 index=i;
293 }
294
295 return index;
297 298 //inversion 逆阵:使用矩阵的初等变换,列主元素消去法 299 public Matrix Inverse() 300 { 301 if (Row != Col) //异常, 非方阵
302 {
303 System.Exception e = new Exception("求逆的矩阵不是方阵"); 304 throw e;
305 }
306StreamWriter sw = new StreamWriter("..\\annex\\close_matrix.txt"); 307 Matrix tmp = new Matrix(this );
308 Matrix ret =new Matrix(Row); //单位阵
309 ret.SetUnit();
310
311 int maxIndex;
312 double dMul;
313
314 for (int i=0;i
315 {
316 maxIndex = tmp.Pivot(i);
317
319 { 320 System.Exception e = new Exception("求逆的矩阵的行列式的值等于0,"); 321 throw e;
322 }
323
324 if (maxIndex != i) //下三角阵中此列的最大值不在当前行,交换
325 {
326 tmp.Exchange(i,maxIndex);
327 ret.Exchange(i,maxIndex);
328
329 }
330
331 ret.Multiple(i,1/tmp[i,i]);
332
333 tmp.Multiple(i,1/tmp[i,i]);
334
335 for (int j=i+1;j
336 {
337 dMul = -tmp[j,i]/tmp[i,i];
338 tmp.MultipleAdd(j,i,dMul);
339 ret.MultipleAdd(j,i,dMul);
341 } 342sw.WriteLine("tmp=\r\n"+tmp); 343sw.WriteLine("ret=\r\n"+ret); 344 }//end for 345 346
347sw.WriteLine("**=\r\n"+ this *ret);
348
349 for (int i=Row-1;i>0;i--)
350 {
351 for (int j=i-1;j>=0;j--)
352 {
353 dMul = -tmp[j,i]/tmp[i,i];
354 tmp.MultipleAdd(j,i,dMul);
355 ret.MultipleAdd(j,i,dMul);
356 }
357 }//end for
358
359
360sw.WriteLine("tmp=\r\n"+tmp);
361sw.WriteLine("ret=\r\n"+ret);
362sw.WriteLine("***=\r\n"+ this *ret); 363sw.Close(); 364 365 return ret; 366 367 }//end Inverse
368
369
370 #region /*
371 //inversion 逆阵:使用矩阵的初等变换,列主元素消去法
372 public Matrix Inverse()
373 {
374 if(Row != Col) //异常, 非方阵
375 {
376 System.Exception e = new Exception("求逆的矩阵不是方阵");
377 throw e;
378 }
379 ///////////////
380 StreamWriter sw = new StreamWriter("..\\annex\\matrix_mul.txt"); 381 ////////////////////
382 ///
383 Matrix tmp = new Matrix(this);
384 Matrix ret =new Matrix(Row); //单位阵 385 ret.SetUnit(); 386 387 int maxIndex; 388 double dMul;
389
390 for(int i=0;i
391 {
392
393 maxIndex = tmp.Pivot(i);
394
395 if(tmp.m_data[maxIndex,i]==0)
396 {
397 System.Exception e = new Exception("求逆的矩阵的行列式的值等于0,"); 398 throw e;
399 }
400
401 if(maxIndex != i) //下三角阵中此列的最大值不在当前行,交换
402 {
403 tmp.Exchange(i,maxIndex);
404 ret.Exchange(i,maxIndex);
405
407 408 ret.Multiple(i,1/tmp[i,i]); 409 410 ///////////////////////// 411 //sw.WriteLine("nul \t"+tmp[i,i]+"\t"+ret[i,i]);
412 ////////////////
413 tmp.Multiple(i,1/tmp[i,i]);
414 //sw.WriteLine("mmm \t"+tmp[i,i]+"\t"+ret[i,i]);
415 sw.WriteLine("111111111 tmp=\r\n"+tmp);
416 for(int j=i+1;j
417 {
418 dMul = -tmp[j,i];
419 tmp.MultipleAdd(j,i,dMul);
420 ret.MultipleAdd(j,i,dMul);
421
422 }
423 sw.WriteLine("[1**********]2=\r\n"+tmp);
424
425 }//end for
426
427
429 for(int i=Row-1;i>0;i--) 430 { 431 for(int j=i-1;j>=0;j--) 432 { 433 dMul = -tmp[j,i];
434 tmp.MultipleAdd(j,i,dMul);
435 ret.MultipleAdd(j,i,dMul);
436 }
437 }//end for
438
439 //////////////////////////
440
441
442 sw.WriteLine("tmp = \r\n" + tmp.ToString());
443
444 sw.Close();
445 ///////////////////////////////////////
446 ///
447 return ret;
448
449 }//end Inverse
451*/ 452 453 #endregion 454 455 //determine if the matrix is square:方阵
456 public bool IsSquare()
457 {
458 return Row==Col;
459 }
460
461 //determine if the matrix is symmetric 对称阵
462 public bool IsSymmetric()
463
464
465 if (Row != Col)
466 return false ;
467
468 for (int i=0;i
469 for (int j=i+1;j
470 if ( m_data[i,j] != m_data[j,i])
471 return false ; {
473 return true ; 474 } 475 476 //一阶矩阵->实数 477 public double ToDouble()
478 {
479 Trace.Assert(Row==1 && Col==1);
480
481 return m_data[0,0];
482 }
483
484 //conert to string
485 public override string ToString()
486
487
488 string s="";
489 for (int i=0;i
490 { {
491 for (int j=0;j
492 s += string .Format("{0} ",m_data[i,j]);
493
495 } 496 return s;
497
498 }
499
500
501 //私有数据成员
502 private double [,] m_data;
503
504 }//end class Matrix
505}
/// 矩阵的乘
public bool MatrixMultiply(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(1) != b.GetLength(0))
return false;
if (a.GetLength(0) != c.GetLength(0) || b.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = 0;
for (int k = 0; k
{
c[i, j] += a[i, k] * b[k, j];
}
}
}
return true;
}
/// 矩阵的加
public bool MatrixAdd(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(0) != b.GetLength(0) || a.GetLength(1) != b.GetLength(1)
|| a.GetLength(0) != c.GetLength(0) || a.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = a[i, j] + b[i, j];
}
}
return true;
}
/// 矩阵的减
public bool MatrixSubtration(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(0) != b.GetLength(0) || a.GetLength(1) != b.GetLength(1)
|| a.GetLength(0) != c.GetLength(0) || a.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = a[i, j] - b[i, j];
}
}
return true;
}
/// 矩阵的行列式的值
public double MatrixSurplus(double[,] a)
{
int i, j, k, p, r, m, n;
m = a.GetLength(0);
n = a.GetLength(1);
double X, temp = 1, temp1 = 1, s = 0,s1 = 0;
if(n == 2)
{
for (i = 0; i
for (j = 0; j
if ((i + j) % 2 > 0) temp1 *= a[i, j];
else temp *= a[i, j];
X=temp-temp1;
}
else
{
for (k = 0; k
{
for (i = 0, j = k; i
temp *= a[i, j];
if (m - i > 0)
{
for (p = m - i, r = m - 1; p > 0; p--, r--)
temp *= a[r, p - 1];
}
s += temp;
temp = 1;
}
for (k = n - 1; k >= 0; k--)
{
for (i = 0, j = k; i = 0; i++, j--)
temp1 *= a[i, j];
if (m - i > 0)
{
for (p = m - 1, r = i; r
temp1 *= a[r, p];
}
s1 += temp1;
temp1 = 1;
}
X = s - s1;
}
return X;
}
/// 矩阵的转置
public bool MatrixInver(double[,] a, ref double[,] b)
{
if (a.GetLength(0) != b.GetLength(1) || a.GetLength(1) != b.GetLength(0)) return false;
for (int i = 0; i
for (int j = 0; j
b[i, j] = a[j, i];
return true;
}
/// 矩阵的逆
public bool MatrixOpp(double[,] a, ref double[,] b)
{
double X = MatrixSurplus(a);
if (X == 0) return false;
X = 1 / X;
double[,] B = new double[a.GetLength(0), a.GetLength(1)];
double[,] SP = new double[a.GetLength(0), a.GetLength(1)];
double[,] AB = new double[a.GetLength(0), a.GetLength(1)];
for (int i = 0; i
for (int j = 0; j
{
for (int m = 0; m
for (int n = 0; n
B[m, n] = a[m, n];
{
for (int x = 0; x
B[i, x] = 0;
for (int y = 0; y
B[y, j] = 0;
B[i, j] = 1;
SP[i, j] = MatrixSurplus(B);
AB[i, j] = X * SP[i, j];
}
}
MatrixInver(AB, ref b);
return true;
}
C#矩阵的运算代码
#region 矩阵运算
///
/// 矩阵对应行列式的值
///
///
///
private double MatrixValue(double[,] MatrixList)
{
int Level = MatrixList.GetLength(1);
double[,] dMatrix = new double[Level, Level];
for (int i = 0; i
{
for (int j = 0; j
{
dMatrix[i, j] = MatrixList[i, j];
}
}
int sign = 1;
for (int i = 0, j = 0; i
{
//判断改行dMatrix[i, j]是否为0,若是,则寻找i 后的行(m,m>i,切dMatrix[m, j]!=0)进行交换 if (dMatrix[i, j] == 0)
{
if (i == Level - 1)
{
return 0;
}
int m = i + 1;
//获取一个dMatrix[m, j]不为为0的行
for (; dMatrix[m, j] == 0; m++)
{
if (m == Level - 1)
{
return 0;
}
}
//判断是否达到矩阵的最大行,若是,则返回0
//把i 行和m 行调换
double temp;
for (int n = j; n
{
temp = dMatrix[i, n];
dMatrix[i, n] = dMatrix[m, n];
dMatrix[m, n] = temp;
}
sign *= (-1);
}
//把当前行以后的行所对应的列变成0
double tmp;
for (int s = Level - 1; s > i; s--)
{
tmp = dMatrix[s, j];
//j行后面的所有行
for (int t = j; t
{
dMatrix[s, t] -= dMatrix[i, t] * (tmp / dMatrix[i, j]);
}
}
}
double result = 1;
for (int i = 0; i
{
if (dMatrix[i, i] != 0)
{
result *= dMatrix[i, i];
}
else
{
return 0;
}
}
return sign * result;
}
///
/// 矩阵减法
///
///
///
///
///
private double[] SubMatrix(double[] A1, double[] A2)
{
//判断矩阵的长短是否一致
int a1 = A1.GetLength(0);
int a2 = A2.GetLength(0);
if (a1 != a2)
{
return null;
}
//矩阵相减
double[] B = new double[a1];
for (int i = 0; i
{
B[i] = A1[i] - A2[i];
}
return B;
}
///
/// 矩阵乘法
///
///
///
///
private double[,] MultiplyMatrix(double[,] firstMatrix, double[,] secondMatrix)
{
double[,] resultMatrix = new double[firstMatrix.GetLength(0), secondMatrix.GetLength(1)]; //判断相乘矩阵是否合法,即第一个矩阵的列要等于第二个矩阵的行
if (firstMatrix.GetLength(1) != secondMatrix.GetLength(0))
{
return null;
}
//求结果矩阵
for (int rowIndex = 0; rowIndex
{
for (int colIndex = 0; colIndex
{
//初始化结果矩阵的元素
resultMatrix[rowIndex, colIndex] = 0;
for (int i = 0; i
{
//求结果矩阵的元素值
resultMatrix[rowIndex, colIndex] += firstMatrix[rowIndex, i] * secondMatrix[i, colIndex]; }
}
}
return resultMatrix;
}
///
/// 求逆矩阵
///
///
///
private double[,] Athwart(double[,] dMatrix)
{
//获取矩阵的行数
int Level = dMatrix.GetLength(1);
double[,] dReverseMatrix = new double[Level, 2 * Level]; //初始化矩阵Level×(2*Level)
for (int i = 0; i
{
for (int j = 0; j
{
if (j
{
dReverseMatrix[i, j] = dMatrix[i, j];
}
else
{
if (j - Level == i)
{
dReverseMatrix[i, j] = 1;
}
else
{
dReverseMatrix[i, j] = 0;
}
}
}
}
for (int i = 0, j = 0; i
if (dReverseMatrix[i, j] == 0)
{
if (i == Level - 1)
{
return null;
}
int m = i + 1;
for (; dMatrix[m, j] == 0; m++)
{
if (m == Level - 1)
{
return null;
}
if (m == Level)
{
return null;
}
else
{
//把i 行和m 行相加
for (int n = j; n
{
dReverseMatrix[i, n] += dReverseMatrix[m, n]; }
}
}
double temp = dReverseMatrix[i, j];
if (temp != 1)
{
//把i 行数据,变成以1开始的一行数据
for (int n = j; n
{
if (dReverseMatrix[i, n] != 0)
{
dReverseMatrix[i, n] /= temp;
}
}
}
//把i 行后的所有行的j 列变成0
for (int s = Level - 1; s > i; s--)
{
temp = dReverseMatrix[s, j];
for (int t = j; t
{
dReverseMatrix[s, t] -= (dReverseMatrix[i, t] * temp); }
}
}
//把矩阵Level×(2*Level)前Level×Level 转变为单位矩阵 for (int i = Level - 2; i >= 0; i--)
{
for (int j = i + 1; j
{
if (dReverseMatrix[i, j] != 0)
double tmp = dReverseMatrix[i, j];
for (int n = j; n
{
dReverseMatrix[i, n] -= (tmp * dReverseMatrix[j, n]); }
}
}
}
//返回逆矩阵
double[,] dReturn = new double[Level, Level]; for (int i = 0; i
{
for (int j = 0; j
{
dReturn[i, j] = dReverseMatrix[i, j + Level]; }
}
return dReturn;
}
#endregion
第二份: 1using System; 2using System.IO;
3using System.Diagnostics;
4
5
6namespace Adjust
7
8{ ///
9 /// Matrix 的摘要说明。
10 /// 实现矩阵的基本运算 11 ///
12 public class Matrix 1314
15 //构造方阵
16 public Matrix(int row)
17 { { 18 m_data = new double [row,row]; 19
20 }
21 public Matrix(int row, int col) 22 {
23 m_data = new double [row,col]; 24 }
25 //复制构造函数
26 public Matrix(Matrix m)
27 {
28 int row = m.Row;
29 int col = m.Col;
30 m_data = new double [row,col]; 31
33 for (int j=0;j
36 }
37
38 /*
39 //分配方阵的大小
40 //对于已含有内存的矩阵,将清空数据 41 public void SetSize(int row) 42 {
43 m_data = new double[row,row]; 44 }
45
46
47 //分配矩阵的大小
48 //对于已含有内存的矩阵,将清空数据 49 public void SetSize(int row,int col) 50 {
51 m_data = new double[row,col]; 52 }
53 */
55 //unit matrix:设为单位阵 56 public void SetUnit() 57 { 58 for (int i=0;i
61 }
62
63 //设置元素值
64 public void SetValue(double d) 65 {
66 for (int i=0;i
69 }
70
71 // Value extraction :返中行数 72 public int Row
73 {
74 get
75 {
77 return m_data.GetLength(0); 78 } 79 } 80 81 //返回列数 82 public int Col
83 {
84 get
85 {
86 return m_data.GetLength(1); 87 }
88 }
89
90 //重载索引
91 //存取数据成员
92 public double this [int row, int col] 93 {
94 get
95 {
96 return m_data[row,col]; 97 }
98 set 99 { 100 m_data[row,col] = value; 101 } 102 } 103
104 //primary change
105 // 初等变换 对调两行:rirj
106 public Matrix Exchange(int i, int j)
107 {
108 double temp;
109
110 for (int k=0;k
111 {
112 temp = m_data[i,k];
113 m_data[i,k] = m_data[j,k];
114 m_data[j,k] = temp;
115 }
116 return this ;
117 }
118
119
121 Matrix Multiple(int index, double mul) 122 { 123 for (int j=0;j
124 {
125 m_data[index,j] *= mul;
126 }
127 return this ;
128 }
129
130
131 //初等变换 第src 行乘以mul 加到第index 行
132 Matrix MultipleAdd(int index, int src, double mul)
133 {
134 for (int j=0;j
135 {
136 m_data[index,j] += m_data[src,j]*mul;
137 }
138
139 return this ;
140 }
141
143 public Matrix Transpose() 144 { 145 Matrix ret = new Matrix(Col,Row); 146
147 for (int i=0;i
148 for (int j=0;j
149 {
150 ret[j,i] = m_data[i,j];
151 }
152 return ret;
153 }
154
155 //binary addition 矩阵加
156 public static Matrix operator + (Matrix lhs,Matrix rhs)
157 {
158 if (lhs.Row != rhs.Row) //异常
159 {
160 System.Exception e = new Exception("相加的两个矩阵的行数不等"); 161 throw e;
162 }
163 if (lhs.Col != rhs.Col) //异常
164 { 165 System.Exception e = new Exception("相加的两个矩阵的列数不等"); 166 throw e; 167 } 168
169 int row = lhs.Row;
170 int col = lhs.Col;
171 Matrix ret=new Matrix(row,col);
172
173 for (int i=0;i
174 for (int j=0;j
175 {
176 double d = lhs[i,j] + rhs[i,j];
177 ret[i,j] = d;
178 }
179 return ret;
180
181 }
182
183 //binary subtraction 矩阵减
184 public static Matrix operator - (Matrix lhs,Matrix rhs)
185 {
186 if (lhs.Row != rhs.Row) //异常 187 { 188 System.Exception e = new Exception("相减的两个矩阵的行数不等"); 189 throw e;
190 }
191 if (lhs.Col != rhs.Col) //异常
192 {
193 System.Exception e = new Exception("相减的两个矩阵的列数不等"); 194 throw e;
195 }
196
197 int row = lhs.Row;
198 int col = lhs.Col;
199 Matrix ret=new Matrix(row,col);
200
201 for (int i=0;i
202 for (int j=0;j
203 {
204 double d = lhs[i,j] - rhs[i,j];
205 ret[i,j] = d;
206 }
207 return ret;
209 210 211 //binary multiple 矩阵乘 212 public static Matrix operator * (Matrix lhs,Matrix rhs)
213 {
214 if (lhs.Col != rhs.Row) //异常
215 {
216 System.Exception e = new Exception("相乘的两个矩阵的行列数不匹配"); 217 throw e;
218 }
219
220 Matrix ret = new Matrix (lhs.Row,rhs.Col);
221 double temp;
222 for (int i=0;i
223 {
224 for (int j=0;j
225 {
226 temp = 0;
227 for (int k=0;k
228 {
229 temp += lhs[i,k] * rhs[k,j];
231 ret[i,j] = temp; 232 } 233 } 234 235 return ret;
236 }
237
238
239 //binary division 矩阵除
240 public static Matrix operator / (Matrix lhs,Matrix rhs)
241 {
242 return lhs * rhs.Inverse();
243 }
244
245 //unary addition 单目加
246 public static Matrix operator + (Matrix m)
247 {
248 Matrix ret = new Matrix(m);
249 return ret;
250 }
251
252 //unary subtraction 单目减 253 public static Matrix operator - (Matrix m) 254 { 255 Matrix ret = new Matrix(m);
256 for (int i=0;i
257 for (int j= 0;j
258 {
259 ret[i,j] = -ret[i,j];
260 }
261
262 return ret;
263 }
264
265 //number multiple 数乘
266 public static Matrix operator * (double d,Matrix m)
267 {
268 Matrix ret = new Matrix(m);
269 for (int i=0;i
270 for (int j=0;j
271 ret[i,j] *= d;
272
273 return ret;
275 276 //number division 数除 277 public static Matrix operator / (double d,Matrix m)
278 {
279 return d*m.Inverse();
280 }
281
282 //功能:返回列主元素的行号
283 //参数:row 为开始查找的行号
284 //说明:在行号[row,Col)范围内查找第row 列中绝对值最大的元素,返回所在行号 285 int Pivot(int row)
286 {
287 int index=row;
288
289 for (int i=row+1;i
290 {
291 if (m_data[i,row] > m_data[index,row])
292 index=i;
293 }
294
295 return index;
297 298 //inversion 逆阵:使用矩阵的初等变换,列主元素消去法 299 public Matrix Inverse() 300 { 301 if (Row != Col) //异常, 非方阵
302 {
303 System.Exception e = new Exception("求逆的矩阵不是方阵"); 304 throw e;
305 }
306StreamWriter sw = new StreamWriter("..\\annex\\close_matrix.txt"); 307 Matrix tmp = new Matrix(this );
308 Matrix ret =new Matrix(Row); //单位阵
309 ret.SetUnit();
310
311 int maxIndex;
312 double dMul;
313
314 for (int i=0;i
315 {
316 maxIndex = tmp.Pivot(i);
317
319 { 320 System.Exception e = new Exception("求逆的矩阵的行列式的值等于0,"); 321 throw e;
322 }
323
324 if (maxIndex != i) //下三角阵中此列的最大值不在当前行,交换
325 {
326 tmp.Exchange(i,maxIndex);
327 ret.Exchange(i,maxIndex);
328
329 }
330
331 ret.Multiple(i,1/tmp[i,i]);
332
333 tmp.Multiple(i,1/tmp[i,i]);
334
335 for (int j=i+1;j
336 {
337 dMul = -tmp[j,i]/tmp[i,i];
338 tmp.MultipleAdd(j,i,dMul);
339 ret.MultipleAdd(j,i,dMul);
341 } 342sw.WriteLine("tmp=\r\n"+tmp); 343sw.WriteLine("ret=\r\n"+ret); 344 }//end for 345 346
347sw.WriteLine("**=\r\n"+ this *ret);
348
349 for (int i=Row-1;i>0;i--)
350 {
351 for (int j=i-1;j>=0;j--)
352 {
353 dMul = -tmp[j,i]/tmp[i,i];
354 tmp.MultipleAdd(j,i,dMul);
355 ret.MultipleAdd(j,i,dMul);
356 }
357 }//end for
358
359
360sw.WriteLine("tmp=\r\n"+tmp);
361sw.WriteLine("ret=\r\n"+ret);
362sw.WriteLine("***=\r\n"+ this *ret); 363sw.Close(); 364 365 return ret; 366 367 }//end Inverse
368
369
370 #region /*
371 //inversion 逆阵:使用矩阵的初等变换,列主元素消去法
372 public Matrix Inverse()
373 {
374 if(Row != Col) //异常, 非方阵
375 {
376 System.Exception e = new Exception("求逆的矩阵不是方阵");
377 throw e;
378 }
379 ///////////////
380 StreamWriter sw = new StreamWriter("..\\annex\\matrix_mul.txt"); 381 ////////////////////
382 ///
383 Matrix tmp = new Matrix(this);
384 Matrix ret =new Matrix(Row); //单位阵 385 ret.SetUnit(); 386 387 int maxIndex; 388 double dMul;
389
390 for(int i=0;i
391 {
392
393 maxIndex = tmp.Pivot(i);
394
395 if(tmp.m_data[maxIndex,i]==0)
396 {
397 System.Exception e = new Exception("求逆的矩阵的行列式的值等于0,"); 398 throw e;
399 }
400
401 if(maxIndex != i) //下三角阵中此列的最大值不在当前行,交换
402 {
403 tmp.Exchange(i,maxIndex);
404 ret.Exchange(i,maxIndex);
405
407 408 ret.Multiple(i,1/tmp[i,i]); 409 410 ///////////////////////// 411 //sw.WriteLine("nul \t"+tmp[i,i]+"\t"+ret[i,i]);
412 ////////////////
413 tmp.Multiple(i,1/tmp[i,i]);
414 //sw.WriteLine("mmm \t"+tmp[i,i]+"\t"+ret[i,i]);
415 sw.WriteLine("111111111 tmp=\r\n"+tmp);
416 for(int j=i+1;j
417 {
418 dMul = -tmp[j,i];
419 tmp.MultipleAdd(j,i,dMul);
420 ret.MultipleAdd(j,i,dMul);
421
422 }
423 sw.WriteLine("[1**********]2=\r\n"+tmp);
424
425 }//end for
426
427
429 for(int i=Row-1;i>0;i--) 430 { 431 for(int j=i-1;j>=0;j--) 432 { 433 dMul = -tmp[j,i];
434 tmp.MultipleAdd(j,i,dMul);
435 ret.MultipleAdd(j,i,dMul);
436 }
437 }//end for
438
439 //////////////////////////
440
441
442 sw.WriteLine("tmp = \r\n" + tmp.ToString());
443
444 sw.Close();
445 ///////////////////////////////////////
446 ///
447 return ret;
448
449 }//end Inverse
451*/ 452 453 #endregion 454 455 //determine if the matrix is square:方阵
456 public bool IsSquare()
457 {
458 return Row==Col;
459 }
460
461 //determine if the matrix is symmetric 对称阵
462 public bool IsSymmetric()
463
464
465 if (Row != Col)
466 return false ;
467
468 for (int i=0;i
469 for (int j=i+1;j
470 if ( m_data[i,j] != m_data[j,i])
471 return false ; {
473 return true ; 474 } 475 476 //一阶矩阵->实数 477 public double ToDouble()
478 {
479 Trace.Assert(Row==1 && Col==1);
480
481 return m_data[0,0];
482 }
483
484 //conert to string
485 public override string ToString()
486
487
488 string s="";
489 for (int i=0;i
490 { {
491 for (int j=0;j
492 s += string .Format("{0} ",m_data[i,j]);
493
495 } 496 return s;
497
498 }
499
500
501 //私有数据成员
502 private double [,] m_data;
503
504 }//end class Matrix
505}
/// 矩阵的乘
public bool MatrixMultiply(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(1) != b.GetLength(0))
return false;
if (a.GetLength(0) != c.GetLength(0) || b.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = 0;
for (int k = 0; k
{
c[i, j] += a[i, k] * b[k, j];
}
}
}
return true;
}
/// 矩阵的加
public bool MatrixAdd(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(0) != b.GetLength(0) || a.GetLength(1) != b.GetLength(1)
|| a.GetLength(0) != c.GetLength(0) || a.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = a[i, j] + b[i, j];
}
}
return true;
}
/// 矩阵的减
public bool MatrixSubtration(double[,] a, double[,] b, ref double[,] c)
{
if (a.GetLength(0) != b.GetLength(0) || a.GetLength(1) != b.GetLength(1)
|| a.GetLength(0) != c.GetLength(0) || a.GetLength(1) != c.GetLength(1)) return false;
for (int i = 0; i
{
for (int j = 0; j
{
c[i, j] = a[i, j] - b[i, j];
}
}
return true;
}
/// 矩阵的行列式的值
public double MatrixSurplus(double[,] a)
{
int i, j, k, p, r, m, n;
m = a.GetLength(0);
n = a.GetLength(1);
double X, temp = 1, temp1 = 1, s = 0,s1 = 0;
if(n == 2)
{
for (i = 0; i
for (j = 0; j
if ((i + j) % 2 > 0) temp1 *= a[i, j];
else temp *= a[i, j];
X=temp-temp1;
}
else
{
for (k = 0; k
{
for (i = 0, j = k; i
temp *= a[i, j];
if (m - i > 0)
{
for (p = m - i, r = m - 1; p > 0; p--, r--)
temp *= a[r, p - 1];
}
s += temp;
temp = 1;
}
for (k = n - 1; k >= 0; k--)
{
for (i = 0, j = k; i = 0; i++, j--)
temp1 *= a[i, j];
if (m - i > 0)
{
for (p = m - 1, r = i; r
temp1 *= a[r, p];
}
s1 += temp1;
temp1 = 1;
}
X = s - s1;
}
return X;
}
/// 矩阵的转置
public bool MatrixInver(double[,] a, ref double[,] b)
{
if (a.GetLength(0) != b.GetLength(1) || a.GetLength(1) != b.GetLength(0)) return false;
for (int i = 0; i
for (int j = 0; j
b[i, j] = a[j, i];
return true;
}
/// 矩阵的逆
public bool MatrixOpp(double[,] a, ref double[,] b)
{
double X = MatrixSurplus(a);
if (X == 0) return false;
X = 1 / X;
double[,] B = new double[a.GetLength(0), a.GetLength(1)];
double[,] SP = new double[a.GetLength(0), a.GetLength(1)];
double[,] AB = new double[a.GetLength(0), a.GetLength(1)];
for (int i = 0; i
for (int j = 0; j
{
for (int m = 0; m
for (int n = 0; n
B[m, n] = a[m, n];
{
for (int x = 0; x
B[i, x] = 0;
for (int y = 0; y
B[y, j] = 0;
B[i, j] = 1;
SP[i, j] = MatrixSurplus(B);
AB[i, j] = X * SP[i, j];
}
}
MatrixInver(AB, ref b);
return true;
}