[petsc-users] How to debug problem with MatSetValues

[Jed Brown]

On Mon, 05 Jul 2010 15:22:50 +0200, Thomas Witkowski <thomas.witkowski at tu-dresden.de> wrote:
> Hi,
> 
> I've some trouble with matrix values that are set by MatSetValues, but 
> are missing in the final matrix. I can reduce the problem to a 75x75 
> matrix that is created on four processors. I create it quite simple with:
> 
> MatCreateMPIAIJ(PETSC_COMM_WORLD, 12, 12, 75, 75, petscMatrix) on p0
> MatCreateMPIAIJ(PETSC_COMM_WORLD, 18, 18, 75, 75, petscMatrix) on p1
> MatCreateMPIAIJ(PETSC_COMM_WORLD, 18, 18, 75, 75, petscMatrix) on p2
> MatCreateMPIAIJ(PETSC_COMM_WORLD, 27, 27, 75, 75, petscMatrix) on p3
> 
> On all processors, the values are set with the following command:
> 
> MatSetValues(petscMatrix, 1, &rowIndex, cols.size(), &(cols[0]), 
> &(values[0]), ADD_VALUES);
> 
> where rowIndex is an integer, cols is of type std::vector<int> and 
> values if of type std::vector<double>. Before MatSetValues is called, I 
> run over the arrays and print all the entries that are added to the 
> matrix. Finally,
> 
> MatAssemblyBegin(petscMatrix, MAT_FINAL_ASSEMBLY);
> MatAssemblyEnd(petscMatrix, MAT_FINAL_ASSEMBLY);
> 
> are called. To check the matrix, I use the option -mat_view_matlab. 
> Okay, now some of the entries are missing and I've absolute no idea what 
> I did wrong (i.e. entry row 48-col 69, which is only once by rank 3).  

Are you checking (49,70) in Matlab (PETSc uses 0-based indexing)?  When
you say the value is missing, is it equal to zero, or not existent in
the sparsity pattern?

Can you confirm that MatSetValues is really being called with
rowIndex=48 and cols containing 69 (e.g. by setting a breakpoint there)?


Matt, std::vector is guaranteed to be contiguous in revisions of C++98
and subsequent draft standards.  I think all implementations observe
this.

  http://stackoverflow.com/questions/247738/is-it-safe-to-assume-that-stl-vector-storage-is-always-contiguous

Jed