00001 00034 #include <itpp/base/fastmath.h> 00035 00036 00037 namespace itpp { 00038 00039 // m=m-v*v'*m 00040 void sub_v_vT_m(mat &m, const vec &v) 00041 { 00042 vec v2(m.cols()); 00043 double tmp, *v2p; 00044 const double *vp; 00045 int i, j; 00046 00047 it_assert(v.size() == m.rows(), "sub_v_vT_m()"); 00048 00049 v2p = v2._data(); 00050 for (j=0; j<m.cols(); j++) { 00051 tmp = 0.0; 00052 vp=v._data(); 00053 for (i=0; i<m.rows(); i++) 00054 tmp += *(vp++) * m._elem(i,j); 00055 *(v2p++) = tmp; 00056 } 00057 00058 vp=v._data(); 00059 for (i=0; i<m.rows(); i++) { 00060 v2p = v2._data(); 00061 for (j=0; j<m.cols(); j++) 00062 m._elem(i,j) -= *vp * *(v2p++); 00063 vp++; 00064 } 00065 } 00066 00067 // m=m-m*v*v' 00068 void sub_m_v_vT(mat &m, const vec &v) 00069 { 00070 vec v2(m.rows()); 00071 double tmp, *v2p; 00072 const double *vp; 00073 int i, j; 00074 00075 it_assert(v.size() == m.cols(), "sub_m_v_vT()"); 00076 00077 v2p = v2._data(); 00078 for (i=0; i<m.rows(); i++) { 00079 tmp = 0.0; 00080 vp = v._data(); 00081 for (j=0; j<m.cols(); j++) 00082 tmp += *(vp++) * m._elem(i,j); 00083 *(v2p++) = tmp; 00084 } 00085 00086 v2p = v2._data(); 00087 for (i=0; i<m.rows(); i++) { 00088 vp=v._data(); 00089 for (j=0; j<m.cols(); j++) 00090 m._elem(i,j) -= *v2p * *(vp++); 00091 v2p++; 00092 } 00093 } 00094 00095 } // namespace itpp
Generated on Thu Aug 30 02:47:18 2007 for IT++ by Doxygen 1.5.3