【基础代码】C# 矩阵类
转载请注明出处：http://blog.csdn.net/CyberZHG
using System;using System.Collections.Generic;using System.Linq;using System.Text;using System.Diagnostics;namespace ZMatrix{ /// <summary> /// 类名：ZMatrix /// 日期：2012年8月19日 /// 作者：ZHG /// 说明：矩阵类，提供基本运算和统计运算。 /// 已有矩阵操作： /// 矩阵加减、乘系数、矩阵乘法、转置、增广、余子阵 /// 高斯消元、秩、行列式、逆矩阵 /// 已有统计操作： /// 求和、平均值、方差、协方差、相关系数 /// </summary> class ZMatrix { private const double eps = 1e-8; private double[,] m; private int h, w; /// <summary> /// 创建矩阵，默认为1*1的矩阵。 /// </summary> public ZMatrix() { this.setSize(1, 1); } /// <summary> /// 创建一个指定大小的矩阵。 /// </summary> /// <param name="h">矩阵的行数</param> /// <param name="w">矩阵的列数</param> public ZMatrix(int h, int w) { this.setSize(h, w); } /// <summary> /// 从已有矩阵创建矩阵，完全复制矩阵的内容。 /// </summary> /// <param name="mat">被复制的矩阵</param> public ZMatrix(ZMatrix mat) { this.setSize(mat.getHeight(), mat.getWidth()); for (int i = 0; i < this.h; ++i) { for (int j = 0; j < this.w; ++j) { this.m[i, j] = mat.m[i, j]; } } } /// <summary> /// 由数组中的数据创建矩阵 /// </summary> /// <param name="value">要设置的数据</param> public ZMatrix(double[,] value) { this.setValue(value); } /// <summary> /// 设置矩阵大小，将清空原有的数据。 /// </summary> /// <param name="h">矩阵的行数</param> /// <param name="w">矩阵的列数</param> public void setSize(int h, int w) { Trace.Assert(h > 0, "矩阵行数必须为正整数"); Trace.Assert(w > 0, "矩阵列数必须为正整数"); this.h = h; this.w = w; this.m = new double[h, w]; } /// <summary> /// 获取行的数目，和getHeight()等效。 /// </summary> /// <returns>矩阵的行数</returns> public int getRowNumber() { return this.h; } /// <summary> /// 获取列的数目，和getColumn()等效。 /// </summary> /// <returns>矩阵的列数</returns> public int getColumnNumber() { return this.w; } /// <summary> /// 获取矩阵高度，和getRowNumber()等效。 /// </summary> /// <returns>矩阵的行数</returns> public int getHeight() { return this.h; } /// <summary> /// 获取矩阵宽度，和getWidth()等效。 /// </summary> /// <returns>矩阵的列数</returns> public int getWidth() { return this.w; } /// <summary> /// 设置矩阵中元素的值。 /// </summary> /// <param name="row">要设置的行</param> /// <param name="col">要设置的列</param> /// <param name="value">要设置的数值</param> public void setValue(int row, int col, double value) { Trace.Assert(row >= 0 && row < h, "矩阵行数不在范围之内"); Trace.Assert(col >= 0 && col < w, "矩阵列数不在范围之内"); this.m[row, col] = value; } /// <summary> /// 设置矩阵中元素的值，矩阵大小自动更改为与数组的大小相同。 /// </summary> /// <param name="value">要设置的数值</param> public void setValue(double[,] value) { this.setSize(value.GetLength(0), value.GetLength(1)); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { this.m[i, j] = value[i, j]; } } } /// <summary> /// 获取元素的数值，默认为第一个元素的值。 /// </summary> /// <param name="row">要获取的行</param> /// <param name="col">要获取的列</param> /// <returns>矩阵中的数值</returns> public double getValue(int row = 0, int col = 0) { Trace.Assert(row >= 0 && row < h, "矩阵行数不在范围之内"); Trace.Assert(col >= 0 && col < w, "矩阵列数不在范围之内"); return this.m[row, col]; } /// <summary> /// 获取元素的数值。 /// </summary> /// <returns>矩阵中数值元素的数组</returns> public double[,] getValue() { return this.m; } /// <summary> /// 在控制台输出矩阵。 /// </summary> public void print() { Console.WriteLine("ZMatrix - Size(" + this.h + ", " + this.w + "): "); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { Console.Write(this.m[i, j].ToString("0.000") + " "); } Console.WriteLine(); } Console.WriteLine(); } /// <summary> /// 返回一个值全部为0的矩阵。 /// </summary> /// <param name="w">矩阵的行数</param> /// <param name="h">矩阵的列数</param> /// <returns>零矩阵</returns> public static ZMatrix ZERO(int h, int w) { ZMatrix res = new ZMatrix(h, w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] = 0.0; } } return res; } /// <summary> /// 返回一个值全部为0的方阵。 /// </summary> /// <param name="size">方阵的大小</param> /// <returns>零矩阵</returns> public static ZMatrix ZERO(int size) { return ZMatrix.ZERO(size, size); } /// <summary> /// 返回一个值全部为1的矩阵。 /// </summary> /// <param name="w">矩阵的行数</param> /// <param name="h">矩阵的列数</param> /// <returns>一矩阵</returns> public static ZMatrix ONE(int h, int w) { ZMatrix res = new ZMatrix(h, w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] = 1.0; } } return res; } /// <summary> /// 返回一个值全部为1的方阵。 /// </summary> /// <param name="size">矩阵的大小</param> /// <returns>一矩阵</returns> public static ZMatrix ONE(int size) { return ZMatrix.ONE(size, size); } /// <summary> /// 返回一个单位阵，和IDENTITY(int)等效。 /// </summary> /// <param name="size">方阵的大小</param> /// <returns>单位阵</returns> public static ZMatrix UNIT(int size) { ZMatrix res = ZMatrix.ZERO(size, size); for (int i = 0; i < size; ++i) { res.m[i, i] = 1.0; } return res; } /// <summary> /// 返回一个单位阵，和UNIT(int)等效。 /// </summary> /// <param name="size">方阵的大小</param> /// <returns>单位阵</returns> public static ZMatrix IDENTITY(int size) { return ZMatrix.UNIT(size); } /// <summary> /// 矩阵求和。 /// </summary> /// <param name="mat">要加上的矩阵</param> /// <returns>求和后的矩阵</returns> public ZMatrix add(ZMatrix mat) { Trace.Assert(h == mat.h && w == mat.w, "矩阵的大小不相同"); ZMatrix res = new ZMatrix(this); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] += mat.m[i, j]; } } return res; } /// <summary> /// 矩阵求差。 /// </summary> /// <param name="mat">要减去的矩阵</param> /// <returns>求差后的矩阵</returns> public ZMatrix subtract(ZMatrix mat) { Trace.Assert(h == mat.h && w == mat.w, "矩阵的大小不相同"); ZMatrix res = new ZMatrix(this); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] -= mat.m[i, j]; } } return res; } /// <summary> /// 将矩阵乘上一个系数。 /// </summary> /// <param name="num">要乘上的系数</param> /// <returns>乘系数后的矩阵</returns> public ZMatrix multiply(double num) { ZMatrix res = new ZMatrix(this); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] *= num; } } return res; } /// <summary> /// 矩阵相乘。 /// </summary> /// <param name="mat">要乘上的矩阵</param> /// <returns>矩阵相乘的结果</returns> public ZMatrix multiply(ZMatrix mat) { Trace.Assert(this.w == mat.h, "不满足乘法条件"); ZMatrix res = ZMatrix.ZERO(this.h, mat.w); for (int i = 0; i < this.h; ++i) { for (int j = 0; j < mat.w; ++j) { for (int k = 0; k < this.w; ++k) { res.m[i, j] += this.m[i, k] * mat.m[k, j]; } } } return res; } /// <summary> /// 矩阵转置。 /// </summary> /// <returns>转置后的矩阵</returns> public ZMatrix transpose() { ZMatrix res = new ZMatrix(this.w, this.h); for (int i = 0; i < res.h; ++i) { for (int j = 0; j < res.w; ++j) { res.m[i, j] = this.m[j, i]; } } return res; } /// <summary> /// 增广矩阵，默认在右侧添加一个单位阵。 /// </summary> /// <returns>增广后的矩阵</returns> public ZMatrix augment() { return this.augment(ZMatrix.UNIT(this.h)); } /// <summary> /// 增广矩阵的行，默认在下方添加一个单位阵。 /// </summary> /// <returns>增广后的矩阵</returns> public ZMatrix rowAugment() { return this.rowAugement(ZMatrix.UNIT(this.w)); } /// <summary> /// 增广矩阵的列，和augment()等效，默认在右侧添加一个单位阵。 /// </summary> /// <returns>增广后的矩阵</returns> public ZMatrix columnAugment() { return this.augment(); } /// <summary> /// 增广矩阵。 /// </summary> /// <param name="mat">将要放置在原矩阵右侧的矩阵</param> /// <returns>增广后的矩阵</returns> public ZMatrix augment(ZMatrix mat) { Trace.Assert(this.h == mat.h, "矩阵的行数不相同"); ZMatrix res = new ZMatrix(this.h, this.w + mat.w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] = this.m[i, j]; } for (int j = 0; j < mat.w; ++j) { res.m[i, j + w] = mat.m[i, j]; } } return res; } /// <summary> /// 增广矩阵的行。 /// </summary> /// <param name="mat">将要放置在原矩阵下方的矩阵</param> /// <returns>增广后的矩阵</returns> public ZMatrix rowAugement(ZMatrix mat) { Trace.Assert(this.w == mat.w, "矩阵的列数不相同"); ZMatrix res = new ZMatrix(this.h + mat.h, this.w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] = this.m[i, j]; } } for (int i = 0; i < mat.h; ++i) { for (int j = 0; j < w; ++j) { res.m[i + h, j] = mat.m[i, j]; } } return res; } /// <summary> /// 增广矩阵的列，与augment(ZMatrix)等效。 /// </summary> /// <param name="mat">将要放置在原矩阵右侧的矩阵</param> /// <returns>增广后的矩阵</returns> public ZMatrix columnAugement(ZMatrix mat) { return this.augment(mat); } /// <summary> /// 交换矩阵的两行，不做返回。 /// </summary> /// <param name="row1">要交换的行</param> /// <param name="row2">要交换的行</param> public void rowSwap(int row1, int row2) { Trace.Assert(row1 >= 0 && row1 < h, "矩阵行数不在范围之内"); Trace.Assert(row2 >= 0 && row2 < h, "矩阵行数不在范围之内"); if(row1 == row2) { double temp = 0.0; for (int i = 0; i < w; ++i) { temp = m[row1, i]; m[row1, i] = m[row2, i]; m[row2, i] = temp; } } } /// <summary> /// 交换矩阵的两列，不做返回。 /// </summary> /// <param name="col1">要交换的列</param> /// <param name="col2">要交换的列</param> public void columnSwap(int col1, int col2) { Trace.Assert(col1 >= 0 && col2 < w, "矩阵列数不在范围之内"); Trace.Assert(col1 >= 0 && col2 < w, "矩阵列数不在范围之内"); if(col1 == col2) { double temp = 0.0; for (int i = 0; i < h; ++i) { temp = m[i, col1]; m[i, col1] = m[i, col2]; m[i, col2] = temp; } } } /// <summary> /// 将某一行乘上一个系数，不做返回。 /// </summary> /// <param name="row">要乘系数的行</param> /// <param name="num">要乘的系数</param> public void rowMultiply(int row, double num) { Trace.Assert(row >= 0 && row < h, "矩阵行数不在范围之内"); for (int i = 0; i < w; ++i) { m[row, i] *= num; } } /// <summary> /// 将某一列乘上一个系数，不做返回。 /// </summary> /// <param name="col">要乘系数的列</param> /// <param name="num">要乘的系数</param> public void columnMultiply(int col, double num) { Trace.Assert(col >= 0 && col < w, "矩阵列数不在范围之内"); for (int i = 0; i < h; ++i) { m[i, col] *= num; } } /// <summary> /// 求余子阵。 /// </summary> /// <param name="row">行</param> /// <param name="col">列</param> /// <returns>余子阵</returns> public ZMatrix complementarySubmatrix(int row, int col) { Trace.Assert(row >= 0 && row < h, "矩阵行数不在范围之内"); Trace.Assert(col >= 0 && col < w, "矩阵列数不在范围之内"); ZMatrix res = new ZMatrix(this.h - 1, this.w - 1); for (int i = 0, ii = 0; i < h; ++i) { if (i == row) { continue; } for (int j = 0, jj = 0; j < w; ++j) { if (j == col) { continue; } res.m[ii, jj] = m[i, j]; ++jj; } ++ii; } return res; } /// <summary> /// 高斯消元法，最终得到的是一个上行梯阵式。 /// </summary> /// <param name="det">是否用于计算行列式，当值为true时将把最后一行做适当处理</param> /// <returns>高斯消元后的矩阵</returns> public ZMatrix gaussian(bool det = false) { ZMatrix res = new ZMatrix(this); bool swap = false; int i = 0; for (int j = 0; i < h && j < w; ++j) { int index = i; for (int k = i + 1; k < h; ++k) { if (Math.Abs(res.m[k, j]) > Math.Abs(res.m[index, j])) { index = k; } } if (index != i) { res.rowSwap(i, index); swap = !swap; } if (Math.Abs(res.m[i, j]) > eps) { for (int k = i + 1; k < h; ++k) { double temp = res.m[k, j] / res.m[i, j]; for (int l = j; l < w; ++l) { res.m[k, l] -= res.m[i, l] * temp; } } ++i; } } if (det) { if (swap) { res.rowMultiply(--i, -1.0); } } return res; } /// <summary> /// 求矩阵的秩。 /// </summary> /// <returns>矩阵的秩</returns> public int rank() { ZMatrix mat = this.gaussian(); int r = 0; for (int i = 0; i < w && i < h; ++i) { if (Math.Abs(mat.m[i, i]) > eps) { ++r; } else { break; } } return r; } /// <summary> /// 求矩阵的行列式，是determinate()的简称。 /// </summary> /// <returns>矩阵的行列式</returns> public double det() { return this.determinate(); } /// <summary> /// 求矩阵的行列式。 /// </summary> /// <returns>矩阵的行列式</returns> public double determinate() { Trace.Assert(w == h, "矩阵不是方阵"); ZMatrix mat = this.gaussian(true); double det = 1.0; for (int i = 0; i < w && i < h; ++i) { det *= mat.m[i, i]; } return det; } /// <summary> /// 求矩阵的逆，是inverse()的简称 /// </summary> /// <returns>矩阵的逆</returns> public ZMatrix inv() { return this.inverse(); } /// <summary> /// 求矩阵的逆 /// </summary> /// <returns>逆矩阵</returns> public ZMatrix inverse() { Trace.Assert(w == h, "矩阵不是方阵"); ZMatrix mat = this.augment().gaussian(); for (int i = h - 1, j = i; i >= 0; --i, --j) { for (int k = i - 1; k >= 0; --k) { double temp = mat.m[k, j] / mat.m[i, j]; for (int l = 0; l < mat.w; ++l) { mat.m[k, l] -= mat.m[i, l] * temp; } } } for (int i = 0; i < h; ++i) { mat.rowMultiply(i, 1 / mat.m[i, i]); } ZMatrix res = new ZMatrix(this.h, this.w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, j] = mat.m[i, j + w]; } } return res; } /// <summary> /// 将矩阵的每一行求和 /// </summary> /// <returns>求和后的矩阵</returns> public ZMatrix rowSum() { ZMatrix res = ZMatrix.ZERO(h, 1); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, 0] += m[i, j]; } } return res; } /// <summary> /// 将矩阵的每一列求和 /// </summary> /// <returns>求和后的矩阵</returns> public ZMatrix columnSum() { ZMatrix res = ZMatrix.ZERO(1, w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[0, j] += m[i, j]; } } return res; } /// <summary> /// 求矩阵中所有元素的和。 /// </summary> /// <returns>总和</returns> public double sum() { double sum = 0.0; for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { sum += m[i, j]; } } return sum; } /// <summary> /// 将矩阵每一行平均值 /// </summary> /// <returns>求平均值后的矩阵</returns> public ZMatrix rowAverage() { ZMatrix res = this.rowSum(); for (int i = 0; i < h; ++i) { res.m[i, 0] /= w; } return res; } /// <summary> /// 将矩阵每一列平均值 /// </summary> /// <returns>求平均值后的矩阵</returns> public ZMatrix columnAverage() { ZMatrix res = this.columnSum(); for (int i = 0; i < w; ++i) { res.m[0, i] /= h; } return res; } /// <summary> /// 求矩阵中所有元素的平均值。 /// </summary> /// <returns>平均值</returns> public double average() { return this.sum() / this.w / this.h; } /// <summary> /// 计算矩阵每一行的方差，是rowVariance()的简称。 /// </summary> /// <returns>求方差后的矩阵</returns> public ZMatrix rowVar() { return this.rowVariance(); } /// <summary> /// 计算矩阵每一行的方差 /// </summary> /// <returns>求方差后的矩阵</returns> public ZMatrix rowVariance() { Trace.Assert(w > 1, "矩阵列数为1，不能计算方差"); ZMatrix ave = this.rowAverage(); ZMatrix res = ZMatrix.ZERO(h, 1); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[i, 0] += (m[i, j] - ave.m[i, 0]) * (m[i, j] - ave.m[i, 0]); } } for (int i = 0; i < h; ++i) { res.m[i, 0] /= w - 1; } return res; } /// <summary> /// 计算矩阵每一行的方差，是colVariance()的简称。 /// </summary> /// <returns>求方差后的矩阵</returns> public ZMatrix colVar() { return this.columnVariance(); } /// <summary> /// 计算矩阵每一列的方差 /// </summary> /// <returns>求方差后的矩阵</returns> public ZMatrix columnVariance() { Trace.Assert(h > 1, "矩阵行数为1，不能计算方差"); ZMatrix ave = this.columnAverage(); ZMatrix res = ZMatrix.ZERO(1, w); for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { res.m[0, j] += (m[i, j] - ave.m[0, j]) * (m[i, j] - ave.m[0, j]); } } for (int i = 0; i < w; ++i) { res.m[0, i] /= h - 1; } return res; } /// <summary> /// 计算矩阵的方差，是variance()的简称。 /// </summary> /// <returns>方差</returns> public double var() { return this.variance(); } /// <summary> /// 计算矩阵的方差。 /// </summary> /// <returns>方差</returns> public double variance() { Trace.Assert(this.w * this.h > 1, "只有一个元素，不能计算方差"); double ave = this.average(); double var = 0.0; for (int i = 0; i < h; ++i) { for (int j = 0; j < w; ++j) { var += (m[i, j] - ave) * (m[i, j] - ave); } } return var / (this.w * this.h - 1); } /// <summary> /// 计算矩阵每一行的标准差，是rowStandardDeviation()的简称。 /// </summary> /// <returns>求标准差后的矩阵</returns> public ZMatrix rowDev() { return this.rowStandardDeviation(); } /// <summary> /// 计算矩阵每一行的标准差。 /// </summary> /// <returns>求标准差后的矩阵</returns> public ZMatrix rowStandardDeviation() { ZMatrix res = this.rowVariance(); for (int i = 0; i < h; ++i) { res.m[i, 0] = Math.Sqrt(res.m[i, 0]); } return res; } /// <summary> /// 计算矩阵每一列的标准差，是columnStandardDeviation()的简称。 /// </summary> /// <returns>求标准差后的矩阵</returns> public ZMatrix colDev() { return this.columnStandardDeviation(); } /// <summary> /// 计算矩阵每一列的标准差。 /// </summary> /// <returns>求标准差后的矩阵</returns> public ZMatrix columnStandardDeviation() { ZMatrix res = this.columnVariance(); for (int i = 0; i < w; ++i) { res.m[0, i] = Math.Sqrt(res.m[0, i]); } return res; } /// <summary> /// 计算矩阵的标准差，是standardDeviation()的简称。 /// </summary> /// <returns>标准差</returns> public double dev() { return this.standardDeviation(); } /// <summary> /// 计算矩阵的标准差。 /// </summary> /// <returns>标准差</returns> public double standardDeviation() { return Math.Sqrt(this.variance()); } /// <summary> /// 计算列之间的协方差矩阵。 /// </summary> /// <returns>协方差矩阵</returns> public ZMatrix cov() { return this.columnCovariance(); } /// <summary> /// 计算行之间的协方差矩阵。 /// </summary> /// <returns>协方差矩阵</returns> public ZMatrix rowCovariance() { Trace.Assert(w > 1, "矩阵列数等于1，无法计算协方差"); ZMatrix res = ZMatrix.ZERO(h, h); ZMatrix ave = this.rowAverage(); for (int i = 0; i < h; ++i) { for (int j = 0; j < h; ++j) { for (int k = 0; k < w; ++k) { res.m[i, j] += (m[i, k] - ave.m[i, 0]) * (m[j, k] - ave.m[j, 0]); } res.m[i, j] /= w - 1; } } return res; } /// <summary> /// 计算列之间的协方差矩阵。 /// </summary> /// <returns>协方差矩阵</returns> public ZMatrix columnCovariance() { Trace.Assert(h > 1, "矩阵行数等于1，无法计算协方差"); ZMatrix res = ZMatrix.ZERO(w, w); ZMatrix ave = this.columnAverage(); for (int i = 0; i < w; ++i) { for (int j = 0; j < w; ++j) { for (int k = 0; k < h; ++k) { res.m[i, j] += (m[k, i] - ave.m[0, i]) * (m[k, j] - ave.m[0, j]); } res.m[i, j] /= h - 1; } } return res; } /// <summary> /// 求列的相关系数矩阵。 /// </summary> /// <returns>相关系数矩阵</returns> public ZMatrix corrcoe() { return this.columnCorrelationCoefficient(); } /// <summary> /// 求行的相关系数矩阵 /// </summary> /// <returns>相关系数矩阵</returns> public ZMatrix rowCorrelationCoefficient() { ZMatrix cov = this.rowCovariance(); ZMatrix dev = this.rowStandardDeviation(); for (int i = 0; i < h; ++i) { for (int j = 0; j < h; ++j) { cov.m[i, j] /= dev.m[i, 0] * dev.m[j, 0]; } } return cov; } /// <summary> /// 求列的相关系数矩阵。 /// </summary> /// <returns>相关系数矩阵</returns> public ZMatrix columnCorrelationCoefficient() { ZMatrix cov = this.columnCovariance(); ZMatrix dev = this.columnStandardDeviation(); for (int i = 0; i < w; ++i) { for (int j = 0; j < w; ++j) { cov.m[i, j] /= dev.m[0, i] * dev.m[0, j]; } } return cov; } }}转载请注明出处：http://blog.csdn.net/CyberZHG