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