title: Cholesky with ScaLAPACK
body: <p>I am trying to do a Cholesky decomposition via <a href="https://software.intel.com/en-us/node/521465" rel="nofollow noreferrer">pdpotrf()</a> of MKL-Intel's library, which uses ScaLAPACK. I am reading the whole matrix in the master node and then distribute it like in this <a href="https://andyspiros.wordpress.com/2011/07/08/an-example-of-blacs-with-c/" rel="nofollow noreferrer">example</a>. Everything works fine when the dimension of the SPD matrix is even. However, when it's odd, <code>pdpotrf()</code> thinks that the matrix is not positive definite.</p>
  
  <p>Could it be because the submatrices are not SPD? I am working with this matrix: <img src="https://i.stack.imgur.com/NOzvb.gif" alt="enter image description here"></p>
  
  <p>and the submatrices are (with 4 processes and blocks of size 2x2):</p>
  
  <pre><code>A_loc on node 0
    4   1   2
    1 0.5   0
    2   0  16
  
  nrows = 3, ncols = 2
  A_loc on node 1
    2 0.5
    0   0
    0   0
  
  nrows = 2, ncols = 3
  A_loc on node 2
    2   0   0
  0.5   0   0
  
  nrows = 2, ncols = 2
  A_loc on node 3
    3   0
    0 0.625
  </code></pre>
  
  <p>Here, every submatrix is not SPD, however, the overall matrix is SPD (have checked with running with 1 process). What should I do? Or there is nothing I can do and <code>pdpotrf()</code> does not work with matrices of odd size?</p>
  
  <hr>
  
  <p>Here is how  I call the routine:</p>
  
  <pre><code>int iZERO = 0;
  int descA[9];
  // N, M dimensions of matrix. lda = N
  // Nb, Mb dimensions of block
  descinit_(descA, &amp;N, &amp;M, &amp;Nb, &amp;Mb, &amp;iZERO, &amp;iZERO, &amp;ctxt, &amp;lda, &amp;info);
  ...
  pdpotrf((char*)"L", &amp;ord, A_loc, &amp;IA, &amp;JA, descA, &amp;info);
  </code></pre>
  
  <p>I also tried this:</p>
  
  <pre><code>// nrows/ncols is the number of rows/columns a submatrix has
  descinit_(descA, &amp;N, &amp;M, &amp;nrows, &amp;ncols, &amp;iZERO, &amp;iZERO, &amp;ctxt, &amp;lda, &amp;info);
  </code></pre>
  
  <p>but I get an error:</p>
  
  <blockquote>
    <p>{    0,    0}:  On entry to {    0,    1}:  On entry to  PDPOTR{    1,
    0}:  On entry to  PDPOTRF parameter number  605 had an illegal value {
    1,    1}:  On entry to  PDPOTRF parameter number  605 had an illegal
    value  F parameter number  605 had an illegal value </p>
    
    <p>PDPOTRF parameter number  605 had an illegal value  info &lt; 0:  If the
    i-th argument is an array and the j-entry had an illegal value, then
    INFO = -(i*100+j), if the i-th argument is a scalar and had an illegal
    value, then INFO = -i.  info = -605</p>
  </blockquote>
  
  <p>From my <a href="https://stackoverflow.com/a/31137250/2411320">answer</a>, you can see what the arguments of the function mean.</p>
  
  <hr>
  
  <p>The code is based on this <a href="https://stackoverflow.com/questions/14147705/cholesky-decomposition-scalapack-error?rq=1">question</a>. Output:</p>
  
  <pre><code>gsamaras@pythagoras:~/konstantis/check_examples$ ../../mpich-install/bin/mpic++ -o test minor.cpp -I../../intel/mkl/include ../../intel/mkl/lib/intel64/libmkl_scalapack_lp64.a       -Wl,--start-group       ../../intel/mkl/lib/intel64/libmkl_intel_lp64.a ../../intel/mkl/lib/intel64/libmkl_core.a  ../../intel/mkl/lib/intel64/libmkl_sequential.a    -Wl,--end-group ../../intel/mkl/lib/intel64/libmkl_blacs_intelmpi_lp64.a -lpthread       -lm     -ldl
  gsamaras@pythagoras:~/konstantis/check_examples$ mpiexec -n 4 ./test
  Processes grid pattern:
  0 1
  2 3
  nrows = 3, ncols = 3
  A_loc on node 0
    4   1   2
    1 0.5   0
    2   0  16
  
  nrows = 3, ncols = 2
  A_loc on node 1
    2 0.5
    0   0
    0   0
  
  nrows = 2, ncols = 3
  A_loc on node 2
    2   0   0
  0.5   0   0
  
  nrows = 2, ncols = 2
  A_loc on node 3
    3   0
    0 0.625
  
  Description init sucesss!
  matrix is not positive definte
  Matrix A result:
    2   1   2 0.5   2
  0.5 0.5   0   0   0
    1   0   1   0 -0.25
  0.25  -1 -0.5 0.625   0
    1  -1  -2 -0.5  14
  </code></pre>