Examples of com.nr.la.Cholesky

com.nr.la.Cholesky
Cholesky decomposition of a matrix A. - PTC Copyright (C) Numerical Recipes Software 1986-2007 Java translation Copyright (C) Huang Wen Hui 2012 @author hwh

    


    // Test cholesky
    System.out.println("Testing cholesky");
    double[][] aposdef = matmul(a,transpose(a));
    Cholesky ach  = new Cholesky(aposdef);
    ach.solve(r,y);
    sbeps = 5.e-15;
    localflag = maxel(vecsub(matmul(aposdef,y),r)) > sbeps;
    globalflag = globalflag || localflag;
    if (localflag) {
      fail("*** cholesky: Error in solve() method");
      
    }


    ranvec(y);
    double[] yy = new double[y.length];
    ach.elmult(y,r);
    ach.elsolve(r,yy);
    sbeps = 5.e-15;
    localflag = maxel(vecsub(y,yy)) > sbeps;
    globalflag = globalflag || localflag;
    if (localflag) {
      fail("*** cholesky: Error in method elmult() or elsolve()");
      
    }


    double[][] ainv = buildMatrix(a);
    ach.inverse(ainv);
    sbeps = 5.e-15;
    localflag = maxel(matsub(matmul(ainv,aposdef),ident(aposdef.length,1.))) > sbeps;
    globalflag = globalflag || localflag;
    if (localflag) {
      fail("*** cholesky: Error in method inverse()");

View Full Code Here

    double tmp,sum,max,oldloglike;
    double[] u = new double[mm],v = new double[mm];
    oldloglike = loglike;
    for (k=0;k<kk;k++) {  // outer loop for computing the p_{nk}'s
 //  Decompose \sigma_{k} in the outer loop     
      Cholesky choltmp = new Cholesky(sig[k]);  
      lndets[k] = choltmp.logdet();
      for (n=0;n<nn;n++) {
        for (m=0;m<mm;m++) u[m] = data[n][m]-means[k][m];
        choltmp.elsolve(u,v);
        for (sum=0.,m=0; m<mm; m++) sum += SQR(v[m]);
        resp[n][k] = -0.5*(sum + lndets[k]) + log(frac[k]);
      }
    }
    loglike = 0;

View Full Code Here

  public Multinormaldev(final long j, final double[] mmean, final double[][] vvar) {
    super(j);
    mm = mmean.length;
    mean = buildVector(mmean);
    var = buildMatrix(vvar);
    chol = new Cholesky(var);
    spt = new double[mm];
    pt = new double[mm];
    if (var[0].length != mm || var.length != mm) throw new IllegalArgumentException("bad sizes");
  }

View Full Code Here

    double[][] aa = NRUtil.buildMatrix(3,3, a);
    double[][] inv = new double[3][3];
    double[] b = new double[]{2.0, 4.0,1.0};
    double[] x = new double[3];
    
    Cholesky ch = new Cholesky(aa);
    ch.solve(b, x);
    System.out.println(NRUtil.toString(x));
    ch.elsolve(b, x);
    System.out.println(NRUtil.toString(x));
    System.out.println(ch.logdet());
    
    ch.inverse(inv);
    
    System.out.println(NRUtil.toString(inv));
  }

View Full Code Here

TOP

Related Classes of com.nr.la.Cholesky

com.nr.ci.Gaumixmod

com.nr.example.Cholesky_Ex

com.nr.ran.Multinormaldev

com.nr.test.test_chapter2.Test_Cholesky

All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.