Number of non-zeros less than nRows

[Pietro_B_]

Hello. I have two very related questions that I hope you can answer. I am writing code that uses the MKL to get the LDL' decomposition of a matrix that is real and symmetric but, in general, indefinite.

First question: given the type of the problem I'm solving, the matrix D might (and should) be composed of 1x1 and 2x2 blocks. Can dss_factor_real(), when the opt parameter is set to MKL_DSS_INDEFINITE, generate such a matrix?

Second question: from the dss_sym_c.c example it appears that, when passing the option MKL_DSS_SYMMETRIC to dss_define_structure(), one does not need to pass the indices and values of the lower triangle of the matrix, but only those of the diagonal and the upper triangle. This said, the LDL' decomposition for the matrix (using Matlab notation)

A=[0 1; 1 0]

should be L=I and D=A. If I call dss_define_structure(), dss_reorder(), and then dss_factor_real() on A, I need to define it with only one nonzero, given that the other one is below the diagonal. As a result, the parameters passed to dss_define_structure() dss_factor_real() should be as follows (see mkl/examples/solverc/source/dss_sym_c.c), where NROWS=NCOLS=2 and NNONZEROS=1:

[cpp]
  _INTEGER_t rowIndex[NROWS + 1] = { 1, 2, 2};
  _INTEGER_t columns[NNONZEROS] = { 2};
  _DOUBLE_PRECISION_t values[NNONZEROS] = { 3};
  _DOUBLE_PRECISION_t rhs[NROWS] = { 1, 6};

[/cpp]

However, when I do so, the program fails after printing "MKL-DSS-DSS-Error, Number of non-zeros less than nRows".

I can avoid this error by passing both off-diagonal values of the matrix, i.e.,

[cpp]
  _INTEGER_t rowIndex[NROWS + 1] = { 1, 2, 3};
  _INTEGER_t columns[NNONZEROS] = { 2, 1};
  _DOUBLE_PRECISION_t values[NNONZEROS] = { 3, 3};
  _DOUBLE_PRECISION_t rhs[NRHS*NROWS] = { 1, 6};[/cpp]

However this seems a little inconsistent with the standard passing method, and I am afraid that it might either cause errors or be inefficient since twice as many numbers must then be passed. Do I have to pass all values of the matrix?

Thanks,

Pietro

[mecej4]

Petro B. wrote:
given the type of the problem I'm solving, the matrix D might (and should) be composed of 1x1 and 2x2 blocks.

The conventioin is that the D in an L.D.L^T decomposition is a diagonal matrix, not a block-diagonal matrix. Often, such routines take D as a 1-D array argument since that representation suffices for a diagonal matrix.

You should read the documentation on the DSS routine that you are calling carefully, and establish whether D is required to be diagonal or it can be block-diagonal.

[Alexander_K_Intel2]

Hi Petro,

We know about this problem, to resolve it just put zero diagonal elements in structure. Which version of MKL do you use currently?

Thanks,

Alexander Kalinkin

[Pietro_B_]

Thank you Alexander. Putting zeros on the diagonal and having three "nonzeros" passed to dss_factor_real (one for the real nonzero off-diagonal and the two zeros on the diagonal) solved the second problem. I am not sure about the version I'm using, it is the MKL that came with the Composer XE 2011 SP1 (Update 13), hence it should be MKL version 10.3.12 according to http://software.intel.com/en-us/articles/which-version-of-the-intel-ipp-intel-mkl-and-intel-tbb-libraries-are-included-in-the-intel

About the first problem, @mecej4: in the reference manual 

http://software.intel.com/sites/products/documentation/doclib/mkl_sa/11/mklman/mklman.pdf

I could not find hints as to what structure the D matrix has in return from dss_factor_real. However, for A=[0 1; 1 0], any D that would allow the L matrix to be lower triangular is a non-diagonal 2x2 matrix D with off-diagonal elements equal to one. Hence dss_factor_real should either return a 2x2 matrix or an error/warning message.

Thanks,

Pietro