【原創】開源Math.NET基礎數學類庫使用(15)C#計算矩陣行列式

資料之巔發表於2015-04-23
C#矩陣

               本部落格所有文章分類的總目錄：【總目錄】本部落格博文總目錄-實時更新 

開源Math.NET基礎數學類庫使用總目錄：【目錄】開源Math.NET基礎數學類庫使用總目錄

　　上個月對Math.NET的基本使用進行了介紹，主要內容有矩陣，向量的相關操作，解析資料格式，數值積分，資料統計，相關函式，求解線性方程組以及隨機數發生器的相關內容。這個月接著深入發掘Math.NET的各種功能，並對原始碼進行分析，使得大家可以儘可能的使用Math.NET在.NET平臺下輕易的開發數學計算相關的，或者可以將其中的原始碼快速移植到自己的系統中去(有時候並不需要所有的功能，只需要其中的部分功能程式碼)，今天要介紹的是Math.NET中利用C#計算矩陣行列式的功能。

本文原文地址：http://www.cnblogs.com/asxinyu/p/4304289.html

1.行列式概念與性質

　　行列式是關於方陣的元素所定義的一種運算，其運算的結果是一個數，稱為方陣的行列式值，簡稱為方陣的行列式。

     

　　行列式的概念最初是伴隨著方程組的求解而發展起來的。行列式的提出可以追溯到十七世紀，最初的雛形由日本數學家關孝和與德國數學家戈特弗裡德·萊布尼茨各自獨立得出，時間大致相同。日本數學家關孝和提出來的，他在1683年寫了一部名為解伏題之法的著作，意思是“解行列式問題的方法”，書中對行列式的概念和它的展開已經有了清楚的敘述。歐洲第一個提出行列式概念的是德國數學家，微積分學奠基人之一萊布尼茨。

　　行列式是線性代數裡面的一個基本概念，我們可以從其定義和性質中瞭解一下其作用：

1.行列式與它的轉置行列式相等；

2.互換行列式的兩行(列)，行列式變號；

3.行列式的某一行(列)的所有的元素都乘以同一數k，等於用數k乘此行列式；

4.行列式如果有兩行(列)元素成比例，則此行列式等於零；

5.若行列式的某一列(行)的元素都是兩數之和，則這個行列式是對應兩個行列式的和；

6.把行列式的某一列(行)的各元素乘以同一數然後加到另一列(行)對應的元素上去，行列式不變

2.Math.NET計算行列式的實現

　　Math.NET在對行列式的計算過程中，只是把其作為矩陣計算的一個小部分功能，作為屬性新增在各個矩陣分解演算法的抽象和實現類中。因為Math.NET中矩陣的泛型型別的相關實現主要是支援Double和Float型別，所以基本上與泛型相關的類都實現了2個版本，在實際使用時再進行呼叫。而矩陣分解演算法如：Cholesky，LU，QR，Svd等都有一個抽象泛型基類。在這些抽象類中都定義好了矩陣分解相關的計算，如行列式，方程求解等功能，然後對類進行繼承，如Cholesky分解演算法，的抽象 基類：

 1 internal abstract class Cholesky : Cholesky<double>
 2 {
 3     protected Cholesky(Matrix<double> factor)
 4         : base(factor)
 5     {
 6     }
 7 
 8     /// <summary>
 9     /// Gets the determinant of the matrix for which the Cholesky matrix was computed.
10     /// </summary>
11     public override double Determinant
12     {
13         get
14         {
15             var det = 1.0;
16             for (var j = 0; j < Factor.RowCount; j++)
17             {
18                 var d = Factor.At(j, j);
19                 det *= d*d;
20             }
21 
22             return det;
23         }
24     }
25 
26     /// <summary>
27     /// Gets the log determinant of the matrix for which the Cholesky matrix was computed.
28     /// </summary>
29     public override double DeterminantLn
30     {
31         get
32         {
33             var det = 0.0;
34             for (var j = 0; j < Factor.RowCount; j++)
35             {
36                 det += 2*Math.Log(Factor.At(j, j));
37             }
38 
39             return det;
40         }
41     }
42 }

這個基類就是繼承實現的Doule型別的版本，然後DenseCholesky繼承該類，實現更多的計算功能：

 1 internal sealed class DenseCholesky : Cholesky
 2 {       
 3     public static DenseCholesky Create(DenseMatrix matrix)
 4     {
 5         if (matrix.RowCount != matrix.ColumnCount)
 6         {
 7             throw new ArgumentException(Resources.ArgumentMatrixSquare);
 8         }
 9         var factor = (DenseMatrix) matrix.Clone();
10         Control.LinearAlgebraProvider.CholeskyFactor(factor.Values, factor.RowCount);
11         return new DenseCholesky(factor);
12     }
13 
14     DenseCholesky(Matrix<double> factor): base(factor)  {      }
15     
16     public override void Solve(Matrix<double> input, Matrix<double> result)
17     {
18         if (result.RowCount != input.RowCount)
19         {
20             throw new ArgumentException(Resources.ArgumentMatrixSameRowDimension);
21         }
22 
23         if (result.ColumnCount != input.ColumnCount)
24         {
25             throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension);
26         }
27 
28         if (input.RowCount != Factor.RowCount)
29         {
30             throw Matrix.DimensionsDontMatch<ArgumentException>(input, Factor);
31         }
32 
33         var dinput = input as DenseMatrix;
34         if (dinput == null)
35         {
36             throw new NotSupportedException("Can only do Cholesky factorization for dense matrices at the moment.");
37         }
38 
39         var dresult = result as DenseMatrix;
40         if (dresult == null)
41         {
42             throw new NotSupportedException("Can only do Cholesky factorization for dense matrices at the moment.");
43         }
44         Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length*Constants.SizeOfDouble);
45         var dfactor = (DenseMatrix) Factor;
46         Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Values, dfactor.RowCount, dresult.Values, dresult.ColumnCount);
47     }
48 
49     public override void Solve(Vector<double> input, Vector<double> result)
50     {
51         if (input.Count != result.Count)
52         {
53             throw new ArgumentException(Resources.ArgumentVectorsSameLength);
54         }
55 
56         if (input.Count != Factor.RowCount)
57         {
58             throw Matrix.DimensionsDontMatch<ArgumentException>(input, Factor);
59         }
60 
61         var dinput = input as DenseVector;
62         if (dinput == null)
63         {
64             throw new NotSupportedException("Can only do Cholesky factorization for dense vectors at the moment.");
65         }
66         var dresult = result as DenseVector;
67         if (dresult == null)
68         {
69             throw new NotSupportedException("Can only do Cholesky factorization for dense vectors at the moment.");
70         }
71 
72         Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length*Constants.SizeOfDouble);
73 
74         var dfactor = (DenseMatrix) Factor;
75         Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Values, dfactor.RowCount, dresult.Values, 1);
76     }
77 }

這裡只是介紹了具體行列式計算的實現，其實在Math.NET中這種實現只要在除錯的時候搞懂了其中一個，其他相關的都好懂了。呼叫的時候，由於矩陣的型別裡面都有相關的屬性，可以直接繼續計算，下面就演示一下如何呼叫計算行列式。

3.矩陣行列式計算的呼叫程式碼

 　　在矩陣類Martix中，已經有一個屬性Determinant可以直接獲取矩陣的行列式，所以計算也非常簡單：

 1 // 格式設定
 2 var formatProvider = (CultureInfo)CultureInfo.InvariantCulture.Clone();
 3 formatProvider.TextInfo.ListSeparator = " ";
 4 
 5 //建立一個隨機的矩陣
 6 var matrix = new DenseMatrix(5);
 7 var rnd = new Random(1); 
 8 for (var i = 0; i < matrix.RowCount; i++)
 9 {
10     for (var j = 0; j < matrix.ColumnCount; j++)
11     {
12         matrix[i, j] = rnd.NextDouble();
13     }
14 }
15 // 計算行列式
16 Console.WriteLine(@"1. 行列式 結果為");
17 Console.WriteLine(matrix.Determinant());
18 Console.WriteLine();

　　輸出結果為：

1. 行列式 結果為
0.0349730711309253

　　如果要使用特殊的分解演算法類計算行列式，也可以單獨計算，例如下面的程式碼，先利用matrix物件生成一個Cholesky分解演算法的物件，然後用它來計算行列式：

 1 var formatProvider = (CultureInfo)CultureInfo.InvariantCulture.Clone();
 2 formatProvider.TextInfo.ListSeparator = " ";
 3 
 4 //建立矩陣
 5 var matrix = DenseMatrix.OfArray(new[,] { { 2.0, 1.0 }, { 1.0, 2.0 } });
 6 Console.WriteLine(@"原始矩陣");
 7 Console.WriteLine(matrix.ToString("#0.00\t", formatProvider));
 8 Console.WriteLine();
 9 
10 //cholesky分解物件
11 var cholesky = matrix.Cholesky();
12 
13 //計算行列式
14 Console.WriteLine(@"矩陣行列式");
15 Console.WriteLine(cholesky.Determinant);
16 Console.WriteLine();

結果如下：

1 原始矩陣
2 DenseMatrix 2x2-Double
3 2.00      1.00
4 1.00      2.00
5 
6 矩陣行列式
7 3

4.資源

　　資源包括原始碼以及案例都可以去官網下載，下載地址本系列文章的目錄中第一篇文章:http://www.cnblogs.com/asxinyu/p/4264638.html,有介紹。由於原始碼很大，如果找不到相應的案例，可以進行搜尋，可以比較快的找到相應的程式碼。

相關文章