真正矩阵的求逆算法

1. float Inverse(CLAYMATRIX& mOut, const CLAYMATRIX& rhs)
2. {
3. CLAYMATRIX m(rhs);
4. DWORD is[4];
5. DWORD js[4];
6. float fDet = 1.0f;
7. int f = 1;
8. for (int k = 0; k < 4; k ++)
9. { float fMax = 0.0f;
10. for (DWORD i = k; i < 4; i ++)
11. {
12. for (DWORD j = k; j < 4; j ++)
13. {
14. const float f = Abs(m(i, j));
15. if (f > fMax)
16. {
17. fMax = f;
18. is[k] = i;
19. js[k] = j;
20. }
21. }
22. }
23. if (Abs(fMax) < 0.0001f)
24. return 0;
25. if (is[k] != k)
26. {
27.    f = -f;
28.    swap(m(k, 0), m(is[k], 0));
29.    swap(m(k, 1), m(is[k], 1));
30.    swap(m(k, 2), m(is[k], 2));
31.    swap(m(k, 3), m(is[k], 3));
32.   }
33. if (js[k] != k)
34.   {
35.    f = -f;
36.     swap(m(0, k), m(0, js[k]));
37.     swap(m(1, k), m(1, js[k]));
38.     swap(m(2, k), m(2, js[k]));
39.     swap(m(3, k), m(3, js[k]));
40.   }
41.    fDet *= m(k, k);
42.    m(k, k) = 1.0f / m(k, k);
43.    for (DWORD j = 0; j < 4; j ++)
44.   {
45.    if (j != k)
46.    m(k, j) *= m(k, k);
47.   }
48.    for (DWORD i = 0; i < 4; i ++)
49. {
50.   if (i != k)
51. {
52.   for (j = 0; j < 4; j ++)
53. {  
54.   if (j != k)
55.    m(i, j) = m(i, j) - m(i, k) * m(k, j);
56.    }
57.   }
58. }
59.   for (i = 0; i < 4; i ++)
60. {
61.   if (i != k)
62.   m(i, k) *= -m(k, k);
63. }
64. }
65. }
66.    for (k = 3; k >= 0; k --)
67. {
68.    if (js[k] != k)
69. {
70.   swap(m(k, 0), m(js[k], 0));
71.   swap(m(k, 1), m(js[k], 1));
72.   swap(m(k, 2), m(js[k], 2));
73.   swap(m(k, 3), m(js[k], 3));
74. }
75. if (is[k] != k)
76. {
77. swap(m(0, k), m(0, is[k]));
78. swap(m(1, k), m(1, is[k]));
79. swap(m(2, k), m(2, is[k]));
80. swap(m(3, k), m(3, is[k]));
81. }
82. }
83. mOut = m;
84. return fDet * f;
85. }

复制代码

不错.
我喜欢.

好东西居然只有一个人识货
我好伤心
                           -----------------天边浮云

呵呵
别伤心
我最近就在为数据结构而烦恼

同意

不错~
很经典~~

可是没多少人懂
我对这个只是了解那么一点
很想找个人或者几个人讨论一下
居然没有人知道
遗憾呀
                                      -------------天边浮云

慢慢来啃