title: Cholesky decomposition of large sparse matrices in Java
body: <p>I want to do Cholesky decomposition of large sparse matrices in Java. Currently I'm using the <a href="https://sites.google.com/site/piotrwendykier/software/parallelcolt" rel="nofollow noreferrer">Parallel Colt</a> library class <a href="http://www.mas.ncl.ac.uk/~ndjw1/share/ParallelColt/doc/cern/colt/matrix/tdouble/algo/decomposition/SparseDoubleCholeskyDecomposition.html" rel="nofollow noreferrer">SparseDoubleCholeskyDecomposition</a> but it's much slower than using <a href="https://stackoverflow.com/questions/22479258/cholesky-decomposition-with-openmp">my code I wrote in C for dense matrices which I use in java with JNI</a>.</p>
  
  <p>For example for 5570x5570 matrices with a non-zero density of 0.25% with SparseDoubleCholeskyDecomposition takes <strong>26.6 seconds</strong> to factor and my own code for the same matrix using dense storage only takes <strong>1.12 seconds</strong>. However, if I set the density to 0.025% then the colt library only takes 0.13 seconds.</p>
  
  <p>I have also written my own sparse matrix Cholesky decomposition in C using <a href="https://en.wikipedia.org/wiki/Sparse_matrix#Compressed_row_Storage_.28CRS_or_CSR.29" rel="nofollow noreferrer">compressed row storage</a> and OpenMP and it's quite a bit faster than SparseDoubleCholeskyDecomposition as well but still slower than my algorithm using dense storage. I could probably optimize it further but in any case the Cholesky factors are dense so it neither solves the speed nor the storage problem.</p>
  
  <p>But I want times in the milliseconds and I eventually need to scale to over 10000x10000 matrices so even my own dense code for dense matrices will become too slow and use too much memory.</p>
  
  <p><a href="https://stackoverflow.com/questions/29603982/cholesky-decomposition-of-sparse-matricies-using-permutation-matrices">I have reason to believe the Cholesky decomposition of sparse matrices can be done much faster and use less memory</a>. Maybe I need a <a href="https://stackoverflow.com/questions/17046585/cholmod-in-java">better Java BLAS library</a>? Is there a library for Java which does Cholesky decomposition for sparse matrices efficiently? Maybe I'm not using Parallel Colt optimally?</p>
  
  <pre><code>import cern.colt.matrix.tdouble.algo.decomposition.SparseDoubleCholeskyDecomposition;
  import cern.colt.matrix.tdouble.impl.*;
  import java.util.Random;
  
  public class Cholesky {
      static SparseRCDoubleMatrix2D create_sparse(double[] a, int n, double r) {
          SparseRCDoubleMatrix2D result = new SparseRCDoubleMatrix2D(n, n);
          Random rand = new Random();
          for(int i=0; i&lt;n; i++) {
              for(int j=i; j&lt;n; j++) {
                  double c = rand.nextDouble();
                  double element = c&lt;r ? rand.nextDouble() : 0;
                  a[i*n+j] = element;
                  a[j*n+i] = element;
                  if (element != 0) {
                      result.setQuick(i, j, element);
                      result.setQuick(j, i, element);
                  }
              }
          }
          for(int i=0; i&lt;n; i++) {
              a[i * n + i] += n;
              result.setQuick(i,i, result.getQuick(i,i) + n);
          }
          return result;
      }
  
      public static void main(String[] args) {
          int n = 5570;
          //int n = 2048;
          double[] a = new double[n*n];
          SparseRCDoubleMatrix2D sparseMatrix = create_sparse(a, n, 0.0025);
  
          long startTime, endTime, duration;
  
          startTime = System.nanoTime();
          SparseDoubleCholeskyDecomposition sparseDoubleCholeskyDecomposition = new SparseDoubleCholeskyDecomposition(sparseMatrix, 0);
          endTime = System.nanoTime();
          duration = (endTime - startTime);
          System.out.printf("colt time construct %.2f s\n", 1E-9*duration);
      }
  }
  </code></pre>