matrix.cs

矩阵类

using System;
using System.Text;

namespace GoldenChicken.MathUtils
{	
	public class Matrix
	{
		private int m_rows;
		private int m_cols;
		private double[,] m_data;

		#region 构造函数

		public Matrix(int size)
		{
			this.m_rows = size;
			this.m_rows = size;
			this.m_data = new double[size,size];			
		}

		public Matrix(int rows,int cols)
		{
			this.m_rows = rows;
			this.m_cols = cols;
			this.m_data = new double[rows,cols];			
		}
		
		public Matrix(double[,] data)
		{			
			this.m_rows = data.GetLength(0);
			this.m_cols = data.GetLength(1);
			this.m_data = data;			
		}
		
		#endregion
		
		#region overrides

		public override string ToString()
		{
			int maxSpace = 0;
			
			//find the string length of the longest number in the matrix
			for(int i = 0; i < m_rows; i++)
			{
				for(int j = 0; j < m_cols; j++)
				{
					int currentLen = this[i,j].ToString().Length;

					if(maxSpace < currentLen)
					{
						maxSpace = currentLen;
					}
				}
			}

			//max space needed is one char more than the longest number
			maxSpace++;

			//calculate an approximate value for the string builder length
			StringBuilder sb = new StringBuilder(maxSpace + (m_rows * m_cols));

			for(int i = 0; i < m_rows; i++)
			{
				for(int j = 0; j < m_cols; j++)
				{
					string currentEle = this[i,j].ToString();

					sb.Append(' ', maxSpace - currentEle.Length);
					sb.Append(currentEle);
				}

				sb.Append("\n");
			}

			return sb.ToString();
		}

		public override int GetHashCode()
		{
			double result = 0;

			for(int i = 0; i < m_rows; i++)
			{
				for(int j = 0; j < m_cols; j++)
				{
					result += this[i,j];
				}
			}

			return (int)result;
		}

		public override bool Equals(Object obj)
		{
			Matrix mtx = (Matrix)obj;

			if(this.rows != mtx.rows || this.cols != mtx.cols)
				return false;

			for(int i = 0; i < this.rows; i++)
			{
				for(int j = 0; j < this.cols; j++)
				{
					if(this[i,j] != mtx[i,j])
						return false;
				}
			}

			return true;
		}

		#endregion

		#region 操作符重载
 
		/// <summary>
		/// 两矩阵是否相等
		/// </summary>
		/// <param name="lmtx"></param>
		/// <param name="rmtx"></param>
		/// <returns></returns>
		public static bool operator == (Matrix lmtx, Matrix rmtx)
		{
			return Equals(lmtx, rmtx);
		}

		/// <summary>
		/// 两矩阵是否不等
		/// </summary>
		/// <param name="lmtx"></param>
		/// <param name="rmtx"></param>
		/// <returns></returns>
		public static bool operator != (Matrix lmtx, Matrix rmtx)
		{
			return !(lmtx == rmtx);
		}

		/// <summary>
		/// 两矩阵相加
		/// </summary>
		/// <param name="lmtx"></param>
		/// <param name="rmtx"></param>
		/// <returns></returns>
		public static Matrix operator + (Matrix lmtx, Matrix rmtx)
		{

			if(lmtx.rows != rmtx.rows || lmtx.cols != rmtx.cols)
				throw new MatrixException("矩阵的维数不匹配");
			//return null;

			Matrix result = new Matrix(lmtx.rows,lmtx.cols);

			for(int i = 0; i < lmtx.rows; i++)
			{
				for(int j = 0; j < lmtx.cols; j++)
				{
					result[i,j] = lmtx[i,j] + rmtx[i,j];
				}
			}

			return result;
		}

		/// <summary>
		/// 两矩阵相减
		/// </summary>
		/// <param name="lmtx"></param>
		/// <param name="rmtx"></param>
		/// <returns></returns>
		public static Matrix operator - (Matrix lmtx, Matrix rmtx)
		{

			if(lmtx.rows != rmtx.rows || lmtx.cols != rmtx.cols)
			{
				throw new MatrixException("矩阵的维数不匹配");
			}

			Matrix result = new Matrix(lmtx.rows,lmtx.cols);

			for(int i = 0; i < lmtx.rows; i++)
			{
				for(int j = 0; j < lmtx.cols; j++)
				{
					result[i,j] = lmtx[i,j] - rmtx[i,j];
				}
			}

			return result;
		}

		/// <summary>
		/// 两矩阵相乘
		/// </summary>
		/// <param name="lmtx"></param>
		/// <param name="rmtx"></param>
		/// <returns></returns>
		public static Matrix operator * (Matrix lmtx, Matrix rmtx)
		{

			if(lmtx.cols != rmtx.rows)
			{
				throw new MatrixException("左矩阵的列数与右矩阵的行数不等");
			}

			Matrix result = new Matrix(lmtx.rows,rmtx.cols);

			for(int i = 0; i < lmtx.rows; i++)
			{
				for(int j = 0; j < rmtx.cols; j++)
				{
					for(int k = 0; k < rmtx.rows; k++)
					{
						result[i,j] += lmtx[i,k] * rmtx[k,j];
					}
				}
			}

			return result;
		}

		/// <summary>
		/// 矩阵的数乘
		/// </summary>
		/// <param name="mtx"></param>
		/// <param name="val"></param>
		/// <returns></returns>
		public static Matrix operator * (Matrix mtx, double val)
		{

			Matrix result = new Matrix(mtx.rows,mtx.cols);

			for(int i = 0; i < mtx.rows; i++)
			{
				for(int j = 0; j < mtx.cols; j++)
				{
					result[i,j] = mtx[i,j] * val;
				}
			}

			return result;
		}

		public static Matrix operator * (double val, Matrix mtx)
		{
			return mtx * val;
		}

		/// <summary>
		/// 矩阵除法
		/// </summary>
		/// <param name="mtx"></param>
		/// <param name="val"></param>
		/// <returns></returns>
		public static Matrix operator / (Matrix mtx, double val)
		{

			if(val == 0)
			{
				throw new MatrixException("分母为零");
			}

			Matrix result = new Matrix(mtx.rows,mtx.cols);

			for(int i = 0; i < mtx.rows; i++)
			{
				for(int j = 0; j < mtx.cols; j++)
				{
					result[i,j] = mtx[i,j] / val;
				}
			}

			return result;
		}

		/// <summary>
		/// 矩阵的幂
		/// </summary>
		/// <param name="mtx"></param>
		/// <param name="val"></param>
		/// <returns></returns>
		public static Matrix operator ^ (Matrix mtx, int val)
		{
			
			if(mtx.rows != mtx.cols)
			{
				throw new MatrixException("矩阵不是方阵");
			}

			Matrix result = mtx;

			for(int i = 0; i < val; i++)
			{
				result = result * mtx;
			}

			return result;
		}

		/// <summary>
		/// 转置矩阵
		/// </summary>
		/// <param name="mtx"></param>
		/// <returns></returns>
		public static Matrix operator ~ (Matrix mtx)
		{

			Matrix result = new Matrix(mtx.cols,mtx.rows);

			for(int i = 0; i < mtx.rows; i++)
			{
				for(int j = 0; j < mtx.cols; j++)
				{
					result[j,i] = mtx[i,j];
				}
			}

			return result;
		}

		/// <summary>
		/// 逆矩阵
		/// </summary>
		/// <param name="mtx"></param>
		/// <returns></returns>
		public static Matrix operator ! (Matrix mtx)
		{
			
			if(mtx.Determinant() == 0)
			{
				throw new MatrixException("矩阵为奇异矩阵");
			}

			//inverse of a 2x2 matrix
			if(mtx.rows == 2 && mtx.cols == 2)
			{
				Matrix tempMtx = new Matrix(2,2);

				tempMtx[0,0] = mtx[1,1];
				tempMtx[0,1] = -mtx[0,1];
				tempMtx[1,0] = -mtx[1,0];
				tempMtx[1,1] = mtx[0,0];

				return tempMtx / mtx.Determinant();
			}

			return mtx.Adjoint()/mtx.Determinant();

		}

		#endregion

		#region 属性
		
		/// <summary>
		/// 矩阵的行数
		/// </summary>
		public int rows
		{
			get
			{
				return this.m_rows;
			}
			set
			{
				this.m_rows = value;
			}
		}

		/// <summary>
		/// 矩阵的列数
		/// </summary>
		public int cols
		{
			get
			{
				return this.m_cols;
			}
			set
			{
				this.m_cols = value;
			}
		}

		/// <summary>
		/// 矩阵元素
		/// </summary>
		public  double this[int i, int j]
		{
			get
			{
				return this.m_data[i, j]; 
			}
			set
			{
				this.m_data[i, j] = value;
			}
		}
		
		/// <summary>
		/// 矩阵是否为方阵
		/// </summary>
		public bool IsSquare
		{
			get
			{
				return (this.m_rows == this.m_cols);
			}
		}

		/// <summary>
		/// 矩阵是否为单位矩阵
		/// </summary>
		public bool IsIdentity
		{
			get
			{
				if(this.m_rows != this.m_cols)
					return false;

				for(int i = 0; i < this.m_rows; i++)
				{
					for(int j = 0; j < this.m_cols; j++)
					{
						if(i == j)
						{
							if(this[i,j] != 1.0) return false;
						}
						else
						{
							if(this[i,j] != 0.0) return false;
						}
					}
				}

				return true;
			}
		}

		/// <summary>
		/// 矩阵是否可逆
		/// </summary>
		public bool IsInvertible
		{
			get
			{
				return (this.Determinant() == 0);
			}
		}		
					
		#endregion

		#region 公用方法

		/// <summary>
		/// 矩阵的行列式
		/// </summary>
		/// <returns></returns>
		public double Determinant()
		{

			double determinent = 0;

			if(this.rows != this.cols)
				throw new MatrixException("Attempt to find the determinent of a non square matrix");
			//return 0;

			//get the determinent of a 2x2 matrix
			if(this.rows == 2 && this.cols == 2)
			{
				determinent = (this[0,0] * this[1,1]) - (this[0,1] * this[1,0]);   
				return determinent;
			}

			Matrix tempMtx = new Matrix(this.rows - 1, this.cols - 1);
 
			//find the determinent with respect to the first row
			for(int j = 0; j < this.cols; j++)
			{
				
				tempMtx = this.Minor(0, j);

				//recursively add the determinents
				determinent += (int)Math.Pow(-1, j) * this[0,j] * tempMtx.Determinant();
				
			}

			return determinent;
		}

		/// <summary>
		/// 矩阵的伴随矩阵
		/// </summary>
		/// <returns></returns>
		public Matrix Adjoint()
		{

			if(this.rows < 2 || this.cols < 2)
			{
				throw new MatrixException("矩阵行、列数均小于2，不存在伴随矩阵");
			}

			Matrix tempMtx = new Matrix(this.rows - 1 , this.cols - 1);
			Matrix adjMtx = new Matrix (this.cols , this.rows);

			for(int i = 0; i < this.rows; i++)
			{
				for(int j = 0; j < this.cols;j++)
				{

					tempMtx = this.Minor(i, j);

					//put the determinent of the minor in the transposed position
					adjMtx[j,i] = (int)Math.Pow(-1,i+j) * tempMtx.Determinant();
				}
			}

			return adjMtx;
		}

		//returns a minor of a matrix with respect to an element
		/// <summary>
		/// 矩阵的指定行、列的代数余子式
		/// </summary>
		/// <param name="row"></param>
		/// <param name="column"></param>
		/// <returns></returns>
		public Matrix Minor(int row, int column)
		{

			if(this.rows < 2 || this.cols < 2)
			{
				throw new MatrixException("矩阵行、列数均小于2，不存在代数余子式");
			}

			int i, j = 0;

			Matrix minorMtx = new Matrix(this.rows - 1, this.cols - 1);
 
			//find the minor with respect to the first element
			for(int k = 0; k < minorMtx.rows; k++)
			{

				if(k >= row)
					i = k + 1;
				else
					i = k;

				for(int l = 0; l < minorMtx.cols; l++)
				{
					if(l >= column)
						j = l + 1;
					else
						j = l;

					minorMtx[k,l] = this[i,j];
				}

			}

			return minorMtx;
		}

		/// <summary>
		/// 重置矩阵为单位矩阵
		/// </summary>
		public void Reset()
		{
			if(this.m_rows != this.m_cols)
			{
				throw new MatrixException("矩阵不是方阵");
			}

			for(int i = 0; i < this.m_rows; i++)
			{
				for(int j= 0; j < 5; j++)
				{
					if(i == j)
					{
						this[i,j] = 1.0;
					}
					else
					{
						this[i,j] = 0.0;
					}
				}
			}
		}

		/// <summary>
		/// 清除矩阵中的数据使之为零矩阵
		/// </summary>
		public void Clear()
		{
			for(int i = 0; i < this.m_rows; i++)
			{
				for(int j = 0; j < this.m_cols; i++)
				{
					this[i,j] = 0.0;
				}
			}
		}

		/// <summary>
		/// 解线性方程
		/// </summary>
		/// <param name="COF"></param>
		/// <param name="CON"></param>
		/// <returns></returns>
		public Matrix LinSolve(Matrix COF, Matrix CON)
		{
			return (!COF) * CON;		
		}

		#endregion

	}

	public class MatrixException : Exception
	{
		public MatrixException(string message) : base(message)
		{

		}
	}
}