[petsc-dev] Bug in MatShift_MPIAIJ ?

[Eric Chamberland]

Hi,

I think I have made simple (valid) example.

Please have a look at the attachment.

Compiled as is, the 2nd process builds a matrix with no local lines, 
copy it the MatShift it so it shows the bug (with debugging=yes and 
-on_error_attach_debugger ddd).  You can comment the line:

#define SET_2nd_PROC_TO_HAVE_NO_LOCAL_LINES

to see that it runs fine if the 2nd process have local lines...

The "magic" happens only if I call MatConvert on the "original" (C) 
matrix...

Hope this helps!

Thanks,

Eric

On 08/10/15 04:18 PM, Barry Smith wrote:
>
>> On Oct 8, 2015, at 2:50 PM, Eric Chamberland <Eric.Chamberland at giref.ulaval.ca> wrote:
>>
>> On 08/10/15 03:30 PM, Satish Balay wrote:
>>> The commented code is sequential code - and shouldn't make a
>>> difference.
>>
>> ...but it does!
>>
>>>
>>> Perhaps your application has other issues.
>>
>>>
>>> Can you verify if your code is valgrind clean?
>>> http://www.mcs.anl.gov/petsc/documentation/faq.html#valgrind
>>
>> Yes it is valgrind clean.  We run every night the same regression tests with valgrind and check morning results since jun 27 2012.  Valgrind run over ~2000 tests (in our suite) and takes almost 13 hours to complete. It is costly, but believe me, we have very useful results with the reports since we do this check.
>>
>> The point is that I do have "0 == aij->nz" because the matrix has 0 lines on one processor.  So why should it pass into  Mat*SetPreallocation?
>>
>> I understand it will speed-up things to preallocated the diagonal if you have not preallocated it, but the criterion of (0 ==  aij->nz) is not right in the case you have no lines on one processor!
>>
>> In other words, after the call to Mat*SetPreallocation, in that case, it must still have (0 == aij->nz) because there are no lines on the processor...
>
>    Sure, but why would that cause a hang or any other problem?  Calling  preallocation on a zero row matrix should be harmless. We need a test case that demonstrates the problem so we can reproduce the problem and determine the fundamental cause.
>
>    Barry
>
>>
>> Thanks!
>>
>> Eric
>>

-------------- next part --------------

static char help[] = "Show MatShift BUG happening after copying a matrix with no rows on a process";

#include <petscmat.h>

#undef __FUNCT__
#define __FUNCT__ "main"

/* DEFINE this to turn on/off the bug: */
#define SET_2nd_PROC_TO_HAVE_NO_LOCAL_LINES

int main(int argc,char **args)
{
  Mat               C;
  PetscInt          i,m = 3;
  PetscMPIInt       rank,size;
  PetscErrorCode    ierr;
  PetscScalar       v;


  PetscInitialize(&argc,&args,(char*)0,help);
  ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
  ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);

  if (2 != size) {
     printf("**** Relevant with 2 processes only.*****\n");
     ierr = PetscFinalize();
     return 1;
  }

  ierr = MatCreate(PETSC_COMM_WORLD,&C);CHKERRQ(ierr);
  PetscInt          locallines;
  PetscInt d_nnz[3] = {0,0,0};
  PetscInt o_nnz[3] = {0,0,0};

#ifdef SET_2nd_PROC_TO_HAVE_NO_LOCAL_LINES
  if (0 == rank) {
    locallines = m;
    d_nnz[0] = 1;
    d_nnz[1] = 1;
    d_nnz[2] = 1;
  }
  else {
   locallines = 0;
  }
#else
if (0 == rank) {
  locallines = m-1;
  d_nnz[0] = 1;
  d_nnz[1] = 1;
}
else {
 locallines = 1;
 d_nnz[0] = 1;
}
#endif

  ierr = MatSetSizes(C,locallines,locallines,m,m);CHKERRQ(ierr);

  ierr = MatSetFromOptions(C);CHKERRQ(ierr);

  ierr = MatXAIJSetPreallocation(C,1,d_nnz,o_nnz,NULL,NULL);CHKERRQ(ierr);

  v = 2;
  //Assembly on the diagonal:
  for (i=0; i<m; i++) {
     ierr = MatSetValues(C,1,&i,1,&i,&v,ADD_VALUES);CHKERRQ(ierr);
  }
  ierr = MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
  ierr = MatSetOption(C,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE);
  ierr = MatSetOption(C, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE);CHKERRQ(ierr);

  Mat lMatA;

  ierr = MatConvert(C,MATSAME, MAT_INITIAL_MATRIX, &lMatA);CHKERRQ(ierr);

  v = -1.0;
  ierr = MatShift(lMatA,v);
  ierr = PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO);CHKERRQ(ierr);
  ierr = MatView(lMatA,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
  ierr = MatView(lMatA,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);

  ierr = MatDestroy(&C);CHKERRQ(ierr);
  ierr = MatDestroy(&lMatA);CHKERRQ(ierr);
  ierr = PetscFinalize();
  return 0;
}