Matrix.h
上传用户:fxromeo
上传日期:2010-04-08
资源大小:89k
文件大小:10k
开发平台:

Visual C++

  1. // Matrix.h 矩阵模板类头文件
  2. // Ver 1.0.0.0
  3. // 版权所有(C) 何渝, 2002
  4. // 最后修改: 2002.5.31.
  5. #ifndef _MATRIX_H //避免多次编译
  6. #define _MATRIX_H
  7. #include <valarray> //模板类valarray的标准头文件
  8. #include <comm.h> //公共头文件
  9. #include <math.h> //数学头文件
  10. template <class _Ty>
  11. class matrix //矩阵类matrix
  12. {
  13. typedef matrix<_Ty> _Myt;
  14.  private:
  15. std::valarray<_Ty> m_Datas; //定义一维数组对象m_Datas
  16. size_t m_stRow; //矩阵行数变量
  17. size_t m_stCol; //矩阵列数变量
  18.  public:
  19. typedef _Ty value_type;
  20. //构造函数一(参数1, 2分别为矩阵的行与列数)
  21. /******
  22. 矩阵类matrix的构造函数一
  23. 构造函数中出现的m_Datas为valarray类的对象,申请stRow * stCol个单元,
  24. 单元内没赋值。对数组对象m_Datas使用了valarray类的构造函数:
  25. explicit valarray(size_t n)
  26.     对私有变量m_stRow和m_stCol分别赋初值stRow, stCol。
  27. ******/
  28. matrix(size_t stRow, size_t stCol)
  29. : m_Datas(stRow * stCol),
  30. m_stRow(stRow), m_stCol(stCol)
  31. {
  32. m_Datas.resize(GetRowNum() * GetColNum(), _Ty(0));
  33. /*
  34. 构造函数中出现的m_Datas为valarray类的对象,若在函数体内调用的
  35. resize(size_t n, const T& c = T())函数有两个参数(参阅valarray中的
  36. 定义),第一个参数申请有矩阵行*列那么多元素个数的存储空间,第二个
  37. 参数对这些申请的空间赋予具有模式T的初值0。如果不调用该函数,则缺省
  38. 情况下m_Datas长度为0;另外,调用该函数但不使用第二个参数,则不会对
  39. m_Datas的任何一个元素赋初值。
  40. */
  41. }
  42. //构造函数二(参数1为指向矩阵的指针,参数2, 3为矩阵行与列数)
  43. /******
  44.  矩阵类matrix的构造函数二
  45.      对私有变量m_stRow和m_stCol分别赋初值stRow, stCol。
  46.      对数组对象m_Datas使用了valarray类的构造函数:
  47. valarray(const _Ty *p, size_t n)
  48.  m_Datas初始化的第一参数为矩阵rhs指针,第二个参数为rhs的元素总个数,
  49.  即rhs行数*列数
  50. ******/
  51. matrix(const _Ty* rhs, size_t stRow, size_t stCol)
  52. : m_Datas(rhs, stRow * stCol),
  53. m_stRow(stRow), m_stCol(stCol)
  54. {
  55. }
  56. //构造函数三(拷贝构造函数,参数为对矩阵matrix的引用)
  57. /******
  58.  矩阵类matrix的构造函数三
  59.      用引用矩阵rhs的数组对象m_Datas初始化matrix所定义对象的m_Datas,
  60.  用引用矩阵rhs的行数rhs.GetRowNum()和列数rhs.GetColNum()分别初始化私
  61.  有变量m_stRow和m_stCol
  62. ******/
  63. matrix(const _Myt& rhs)
  64. : m_Datas(rhs.m_Datas),
  65. m_stRow(rhs.GetRowNum()), m_stCol(rhs.GetColNum())
  66. {
  67. }
  68. size_t GetRowNum() const //返回矩阵行数的函数
  69. {
  70. return m_stRow;
  71. }
  72. size_t GetColNum() const //返回矩阵列数的函数
  73. {
  74. return m_stCol;
  75. }
  76. /*****
  77.     重载运算符(),获取二维矩阵某元素在一维数组m_Datas的位置
  78. (matrix类定义的矩阵实际上已经被转化为一维数组m_Datas)
  79.     当数组元素(i, j)在运算符左边(被写)时,使用_Ty& operator ()
  80. 当数组元素(i, j)在运算符右边(被读)时,使用_Ty operator ()    
  81. *****/
  82. _Ty& operator () (size_t stRow, size_t stCol)
  83. {
  84. Assert(stRow < GetRowNum()); //断定stRow不超实际矩阵行值
  85. Assert(stCol < GetColNum()); //断定stCol不超实际矩阵列值
  86. return m_Datas[stRow * GetColNum() + stCol];
  87. }
  88. /******
  89.  重载()运算符,返回二维矩阵中元素(stRow, stCol)在一维数组中的位置
  90.  stRow与stCol分别为矩阵某元素行与列的位置数
  91.  当矩阵元素在运算符右边(被读)时使用
  92. ******/
  93. const _Ty operator () (size_t stRow, size_t stCol) const
  94. {
  95. Assert(stRow < GetRowNum()); //断定stRow不超实际矩阵行值
  96. Assert(stCol < GetColNum()); //断定stCol不超实际矩阵列值
  97. return m_Datas[stRow * GetColNum() + stCol];
  98. }
  99. // 赋值操作符
  100. //矩阵与矩阵的自反*, +, -运算符
  101. _Myt& operator += (const _Myt& rhs) //矩阵与矩阵的自反+
  102. {
  103. Assert(GetRowNum() == rhs.GetRowNum()); //断定两矩阵的行数相等
  104. Assert(GetColNum() == rhs.GetColNum()); //断定两矩阵的列数相等
  105. //利用类valarray定义,使其对象m_Datas定义进行两矩阵相加
  106. m_Datas += rhs.m_Datas;
  107. //结果放在左边矩阵中,返回指向左边矩阵的指针
  108. return *this;
  109. }
  110. _Myt& operator -= (const _Myt& rhs) //矩阵与矩阵的自反-
  111. {
  112. Assert(GetRowNum() == rhs.GetRowNum());
  113. Assert(GetColNum() == rhs.GetColNum());
  114. //利用类valarray定义,使其对象m_Datas定义进行两矩阵相减
  115. m_Datas -= rhs.m_Datas;
  116. //结果放在左边矩阵中,返回指向左边矩阵的指针
  117. return *this;
  118. }
  119. _Myt& operator *= (const _Myt& rhs) //矩阵与矩阵的自反*
  120. {
  121. MatrixMultiply(*this, *this, rhs); //调用矩阵乘法函数
  122. return *this; //最后结果放在左边矩阵中
  123. }
  124. //矩阵自反加、减、乘以、除以数
  125. _Myt& operator += (const _Ty& rhs) //矩阵自加数
  126. {
  127. m_Datas += rhs; //利用数组对象对每个元素加数
  128. return *this; //结果放在原矩阵(数组m_Datas)中
  129. }
  130. _Myt& operator -= (const _Ty& rhs) //矩阵自减数
  131. {
  132. m_Datas -= rhs;
  133. return *this;
  134. }
  135. _Myt& operator *= (const _Ty& rhs) //矩阵自乘数
  136. {
  137. m_Datas *= rhs;
  138. return *this;
  139. }
  140. _Myt& operator /= (const _Ty& rhs) //矩阵自除以数
  141. {
  142. m_Datas /= rhs;
  143. return *this;
  144. }
  145. // 一元操作符  对矩阵(每个元素)赋予+或-号
  146. _Myt operator + () const //赋+号
  147. {
  148. return *this; //不用作任何处理,维持原状
  149. }
  150. _Myt operator - () const //赋-号
  151. {
  152. _Myt mat(*this);
  153. mat.m_Datas = -mat.m_Datas; //每个元素赋-号
  154. return mat;
  155. }
  156. // 二元操作符
  157. //矩阵加数 mat = lhs + rhs
  158. friend _Myt operator + (const _Myt& lhs, const _Ty& rhs)
  159. {
  160. _Myt mat(lhs); //新生成一新矩阵对象mat
  161. mat.m_Datas += rhs; //对新矩阵对象每个元素加数
  162. return mat; //返回新矩阵对象
  163. }
  164. //矩阵减数 mat = lhs - rhs
  165. friend _Myt operator - (const _Myt& lhs, const _Ty& rhs)
  166. {
  167. _Myt mat(lhs);
  168. mat.m_Datas -= rhs;
  169. return mat;
  170. }
  171. //矩阵乘以数 mat = lhs * rhs
  172. friend _Myt operator * (const _Myt& lhs, const _Ty& rhs)
  173. {
  174. _Myt mat(lhs); //新生成一新矩阵对象mat
  175. mat.m_Datas *= rhs; //对新矩阵对象每个元素乘以数
  176. return mat; //返回新矩阵对象
  177. }
  178. //矩阵除以数 mat = lhs / rhs
  179. friend _Myt operator / (const _Myt& lhs, const _Ty& rhs)
  180. {
  181. _Myt mat(lhs);
  182. mat.m_Datas /= rhs;
  183. return mat;
  184. }
  185. //数加矩阵 mat = lhs + rhs
  186. friend _Myt operator + (const _Ty& lhs, const _Myt& rhs)
  187. {
  188. _Myt mat(rhs); //新生成一新矩阵对象mat
  189. mat.m_Datas += lhs; //数加上新矩阵对象的每个元素
  190. return mat;
  191. }
  192. //数减矩阵 mat = lhs - rhs
  193. friend _Myt operator - (const _Ty& lhs, const _Myt& rhs)
  194. {
  195. _Myt mat(rhs); //新生成一新矩阵对象mat
  196. mat.m_Datas -= lhs; //数减新矩阵对象的每个元素
  197. return mat;
  198. }
  199. //数乘以矩阵 mat = lhs * rhs
  200. friend _Myt operator * (const _Ty& lhs, const _Myt& rhs)
  201. {
  202. _Myt mat(rhs); //新生成一新矩阵对象mat
  203. mat.m_Datas *= lhs; //对新矩阵对象每个元素乘以数
  204. return mat;
  205. }
  206. //矩阵加矩阵 mat = lhs + rhs
  207. friend _Myt operator + (const _Myt& lhs, const _Myt& rhs)
  208. {
  209. _Myt mat(lhs); //新生成一新矩阵对象mat, 用左边阵初始化
  210. mat.m_Datas += rhs.m_Datas; //加上右边矩阵对象每个相应元素
  211. return mat; //返回新矩阵对象
  212. }
  213. //矩阵减矩阵 mat = lhs - rhs
  214. friend _Myt operator - (const _Myt& lhs, const _Myt& rhs)
  215. {
  216. _Myt mat(lhs); //新生成一新矩阵对象mat, 用左边阵初始化
  217. mat.m_Datas -= rhs.m_Datas; //减去右边矩阵对象每个相应元素
  218. return mat; //返回新矩阵对象
  219. }
  220. //矩阵乘以矩阵 mTmp = lhs * rhs
  221. friend _Myt operator * (const _Myt& lhs, const _Myt& rhs)
  222. { //生成一个矩阵新对象mTmp
  223. _Myt mTmp(lhs.GetRowNum(), rhs.GetColNum()); //没初始化
  224. MatrixMultiply(mTmp, lhs, rhs); //调用矩阵乘法函数
  225. return mTmp; //用新矩阵对象返回
  226. }
  227. //矩阵比较
  228. //比较两矩阵是否不相等,若相等返回true,否则返回false
  229. friend bool operator != (const _Myt& lhs, const _Myt& rhs)
  230. {
  231. if (lhs.GetRowNum() != rhs.GetRowNum())
  232. return true;
  233. if (lhs.GetColNum() != rhs.GetColNum())
  234. return true;
  235. for (size_t i = 0; i < lhs.m_Datas.size(); i ++)
  236. {
  237. if (lhs.m_Datas[i] != rhs.m_Datas[i])
  238. return true;
  239. }
  240. return false;
  241. }
  242. //比较两矩阵是否相等,若相等返回true,否则返回false
  243. friend bool operator == (const _Myt& lhs, const _Myt& rhs)
  244. {
  245. if (lhs.GetRowNum() != rhs.GetRowNum())
  246. return false;
  247. if (lhs.GetColNum() != rhs.GetColNum())
  248. return false;
  249. for (size_t i = 0; i < lhs.m_Datas.size(); i ++)
  250. {
  251. if (lhs.m_Datas[i] != rhs.m_Datas[i])
  252. return false;
  253. }
  254. return true;
  255. }
  256. };
  257. typedef matrix<float> matrixf;
  258. typedef matrix<double> matrixd;
  259. typedef matrix<long double> matrixld;
  260. //矩阵乘法函数
  261. template <class _Tyout, class _Tylhs, class _Tyrhs> //最后结果在mOut中
  262. matrix<_Tyout>& MatrixMultiply(matrix<_Tyout>& mOut, const matrix<_Tylhs>& lhs, const matrix<_Tyrhs>& rhs);
  263. //输出矩阵 一行一行输出全矩阵
  264. template <class _Ty>
  265. void MatrixLinePrint(const matrix<_Ty>& mOut);
  266. //输出矩阵 输出指定的某行
  267. template <class _Ty>
  268. void MatrixLinePrint(const matrix<_Ty>& mOut, size_t LineNo);
  269. //矩阵转置 原阵在mIn,转置后的矩阵在mOut
  270. template <class _Ty>
  271. void MatrixTranspose(matrix<_Ty>& mIn, matrix<_Ty>& mOut);
  272. //判断矩阵对称
  273. template <class _Ty> //对称返回true,否则返回false
  274. bool MatrixSymmetry(const matrix<_Ty>& rhs);
  275. //判断矩阵对称正定
  276. template <class _Ty>
  277. int MatrixSymmetryRegular(const matrix<_Ty>& rhs, int sym);
  278. //全选主元法求矩阵行列式
  279. template <class _Ty> //返回值为行列式值
  280. long double MatrixDeterminant(const matrix<_Ty>& rhs);
  281. //全选主元高斯(Gauss)消去法求一般矩阵的秩
  282. template <class _Ty> //返回值为秩数
  283. size_t MatrixRank(const matrix<_Ty>& rhs);
  284. //用全选主元高斯-约当(Gauss-Jordan)消去法计算实(复)矩阵的逆矩阵
  285. //矩阵类型必须是浮点型
  286. template <class _Ty>
  287. int MatrixInversionGS(matrix<_Ty>& rhs);
  288. //用“变量循环重新编号法”法求对称正定矩阵逆
  289. //矩阵类型必须是浮点型
  290. template <class _Ty>
  291. int MatrixSymmetryRegularInversion(matrix<_Ty>& rhs);
  292. //特兰持(Trench)法求托伯利兹(Toeplitz)矩阵逆
  293. //矩阵类型必须是浮点型
  294. template <class _Ty>
  295. int MatrixToeplitzInversionTrench(const valarray<_Ty>& t, 
  296.   const valarray<_Ty>& tt, matrix<_Ty>& rhs);
  297. //实矩阵LU分解
  298. //矩阵类型必须是浮点型
  299. template <class _Ty>
  300. int MatrixLU(const matrix<_Ty>& rhs, matrix<_Ty>& lhs, matrix<_Ty>& uhs);
  301. //用豪斯荷尔德(Householder)变换对一般m*n阶的实矩阵进行QR分解
  302. //矩阵类型必须是浮点型
  303. template <class _Ty>
  304. int MatrixQR(matrix<_Ty>& rhs, matrix<_Ty>& rhq);
  305. //对称正定阵的乔里斯基(Cholesky)分解及求其行列式值
  306. //矩阵类型必须是浮点型
  307. template <class _Ty>
  308. long double MatrixSymmetryRegularCholesky(const matrix<_Ty>& rhs);
  309. //一般实矩阵的奇异值分解
  310. //矩阵类型必须是浮点型
  311. template <class _Ty>
  312. int MatrixSingularValue(matrix<_Ty>& a, matrix<_Ty>& u, 
  313. matrix<_Ty>& v, _Ty eps);
  314. //广义逆的奇异值分解
  315. //矩阵类型必须是浮点型
  316. template <class _Ty>
  317. int GeneralizedInversionSingularValue(matrix<_Ty>& a, matrix<_Ty>& aa, 
  318. _Ty eps, matrix<_Ty>& u, matrix<_Ty>& v);
  319. #include "Matrix.inl" //类及相关函数的定义头文件
  320. #endif // _MATRIX_H