Changeset 75
Legend:
- Unmodified
- Added
- Removed
-
bdm/math/libDC.cpp
r37 r75 162 162 double x = 0.0, sum; 163 163 int i,j; 164 vec s(v.length()); 165 vec S=L*v; 164 166 165 167 for ( i=0; i<D.length(); i++ ) { //rows of L 166 168 sum = 0.0; 167 for ( j= i; j<=i; j++ ){sum+=L( i,j )*v( i);}169 for ( j=0; j<=i; j++ ){sum+=L( i,j )*v( j );} 168 170 x +=D( i )*sum*sum; 171 s(i)=sum; 172 }; 173 return x; 174 } 175 176 double ldmat::invqform( const vec &v ) const { 177 double x = 0.0; 178 int i; 179 vec pom(v.length()); 180 181 backward_substitution(L.T(),v,pom); 182 183 for ( i=0; i<D.length(); i++ ) { //rows of L 184 x +=pom(i)*pom(i)/D(i); 169 185 }; 170 186 return x; … … 298 314 j = i + k; //change in .m 1+1=2, here 0+0+1=1 299 315 s = L( j, i ); 300 for ( m = i + 1; m < ( j - 1); m++ ) {316 for ( m = i + 1; m < ( j ); m++ ) { 301 317 s += L( m, i ) * Il( j, m ); 302 318 } -
bdm/math/libDC.h
r51 r75 33 33 public: 34 34 /*! 35 * Perfroms a rank-1 update by outer product of vectors: $V = V + w v v'$.35 * Perfroms a rank-1 update by outer product of vectors: \f$V = V + w v v'\f$. 36 36 * @param v Vector forming the outer product to be added 37 37 * @param w weight of updating; can be negative … … 46 46 virtual mat to_mat() =0; 47 47 48 /*! \brief Inplace symmetric multiplication by a SQUARE matrix $C$, i.e. $V = C*V*C'$48 /*! \brief Inplace symmetric multiplication by a SQUARE matrix \f$C\f$, i.e. \f$V = C*V*C'\f$ 49 49 @param C multiplying matrix, 50 50 */ 51 51 virtual void mult_sym ( const mat &C ) =0; 52 52 53 /*! \brief Inplace symmetric multiplication by a SQUARE transpose of matrix $C$, i.e. $V = C'*V*C$53 /*! \brief Inplace symmetric multiplication by a SQUARE transpose of matrix \f$C\f$, i.e. \f$V = C'*V*C\f$ 54 54 @param C multiplying matrix, 55 55 */ … … 64 64 65 65 /*! 66 \brief Multiplies square root of $V$ by vector $x$.66 \brief Multiplies square root of \f$V\f$ by vector \f$x\f$. 67 67 68 68 Used e.g. in generating normal samples. … … 71 71 72 72 /*! 73 \brief Evaluates quadratic form $x= v'*V*v$;73 \brief Evaluates quadratic form \f$x= v'*V*v\f$; 74 74 75 75 */ 76 76 virtual double qform (const vec &v ) const =0; 77 78 /*! 79 \brief Evaluates quadratic form \f$x= v'*inv(V)*v\f$; 80 81 */ 82 virtual double invqform (const vec &v ) const =0; 77 83 78 84 // //! easy version of the … … 112 118 void mult_sym ( const mat &C); 113 119 void mult_sym_t ( const mat &C); 114 //! store result of \c mult_sym in external matrix $U$120 //! store result of \c mult_sym in external matrix \f$U\f$ 115 121 void mult_sym ( const mat &C, fsqmat &U) const; 116 //! store result of \c mult_sym_t in external matrix $U$122 //! store result of \c mult_sym_t in external matrix \f$U\f$ 117 123 void mult_sym_t ( const mat &C, fsqmat &U) const; 118 124 void clear(); … … 140 146 double logdet() const {return log ( det ( M ) );}; 141 147 double qform (const vec &v ) const {return ( v* ( M*v ) );}; 148 double invqform (const vec &v ) const {return ( v* ( itpp::inv(M)*v ) );}; 142 149 vec sqrt_mult (const vec &v ) const {mat Ch=chol(M); return Ch*v;}; 143 150 … … 163 170 /*! \brief Matrix stored in LD form, (typically known as UD) 164 171 165 Matrix is decomposed as follows: \f[M = L'DL\f] where only $L$ and $D$ matrices are stored.172 Matrix is decomposed as follows: \f[M = L'DL\f] where only \f$L\f$ and \f$D\f$ matrices are stored. 166 173 All inplace operations modifies only these and the need to compose and decompose the matrix is avoided. 167 174 */ … … 194 201 double logdet() const; 195 202 double qform (const vec &v ) const; 203 double invqform (const vec &v ) const; 196 204 // sqmat& operator -= ( const sqmat & ld2 ); 197 205 void clear(); … … 205 213 virtual void inv ( ldmat &Inv ) const; 206 214 207 /*! \brief Symmetric multiplication of $U$ by a general matrix $C$, result of which is stored in the current class.215 /*! \brief Symmetric multiplication of \f$U\f$ by a general matrix \f$C\f$, result of which is stored in the current class. 208 216 @param C matrix to multiply with 209 217 @param U a space where the inverse is stored. … … 211 219 void mult_sym ( const mat &C, ldmat &U) const; 212 220 213 /*! \brief Symmetric multiplication of $U$ by a transpose of a general matrix $C$, result of which is stored in the current class.221 /*! \brief Symmetric multiplication of \f$U\f$ by a transpose of a general matrix \f$C\f$, result of which is stored in the current class. 214 222 @param C matrix to multiply with 215 223 @param U a space where the inverse is stored. … … 218 226 219 227 220 /*! \brief Transforms general $A'D0 A$ into pure $L'DL$221 222 The new decomposition fullfills: $A'*diag(D)*A = self.L'*diag(self.D)*self.L$228 /*! \brief Transforms general \f$A'D0 A\f$ into pure \f$L'DL\f$ 229 230 The new decomposition fullfills: \f$A'*diag(D)*A = self.L'*diag(self.D)*self.L\f$ 223 231 @param A general matrix 224 232 @param D0 general vector … … 245 253 246 254 protected: 247 //! Positive vector $D$255 //! Positive vector \f$D\f$ 248 256 vec D; 249 //! Lower-triangular matrix $L$257 //! Lower-triangular matrix \f$L\f$ 250 258 mat L; 251 259