root/library/bdm/math/chmat.cpp @ 434

#include "chmat.h"

//using std::endl;


void chmat::opupdt ( const vec &v, double w ) {
//TODO see cholupdt in lhotse
        mat Z;
        mat R;
        mat V(1,v.length());
        V.set_row(0,v*sqrt(w));
        Z = concat_vertical ( Ch,V );
        qr ( Z,R );
        Ch = R ( 0, Ch.rows()-1, 0, Ch.cols()-1 );
};
mat chmat::to_mat() const {mat F=Ch.T() *Ch;return F;};
void chmat::mult_sym ( const mat &C ) {
        it_error ( "not implemented" );
};
void chmat::mult_sym ( const mat &C , chmat &U ) const {
        mat Z=C*Ch;
        U.Ch= chol(Z*Z.T());
};
void chmat::mult_sym_t ( const mat &C ) {
        it_error ( "not implemented" );
};
void chmat::mult_sym_t ( const mat &C, chmat &U ) const {
        it_error ( "not implemented" );
};
double chmat::logdet() const {
        double ldet=0.0; int i;
        //sum of logs of (possibly negative!) diagonal entries
        for ( i=0;i<Ch.rows();i++ ) {ldet+=log ( std::fabs ( Ch ( i,i ) ) );}
        return 2*ldet; //compensate for Ch being sqrt()
};
//TODO can be done more efficiently using BLAS, see triangular matrices
vec chmat::sqrt_mult ( const vec &v ) const {vec pom; pom = Ch*v; return pom;};
double chmat::qform ( const vec &v ) const {vec pom; pom = Ch*v; return pom*pom;};
double chmat::invqform ( const vec &v ) const {
        vec pom(v.length());
        forward_substitution(Ch.T(),v,pom);
        return pom*pom;
};
void chmat::clear() {Ch.clear();};