2: static char help[] = "Tests incorrect use of MatDiagonalSet() for SHELL matrices\n\n";

  4: #include <petscmat.h>

  6: typedef struct _n_User *User;
  7: struct _n_User {
  8:   Mat B;
  9: };

 11: static PetscErrorCode MatGetDiagonal_User(Mat A,Vec X)
 12: {
 13:   User           user;

 17:   MatShellGetContext(A,&user);
 18:   MatGetDiagonal(user->B,X);
 19:   return(0);
 20: }

 22: int main(int argc,char **args)
 23: {
 24:   const PetscScalar xvals[] = {11,13};
 25:   const PetscInt    inds[]  = {0,1};
 26:   PetscScalar       avals[] = {2,3,5,7};
 27:   Mat               A,S;
 28:   Vec               X,Y;
 29:   User              user;
 30:   PetscErrorCode    ierr;

 32:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 33:   MatCreateSeqAIJ(PETSC_COMM_WORLD,2,2,2,NULL,&A);
 34:   MatSetUp(A);
 35:   MatSetValues(A,2,inds,2,inds,avals,INSERT_VALUES);
 36:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 37:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 38:   VecCreateSeq(PETSC_COMM_WORLD,2,&X);
 39:   VecSetValues(X,2,inds,xvals,INSERT_VALUES);
 40:   VecAssemblyBegin(X);
 41:   VecAssemblyEnd(X);
 42:   VecDuplicate(X,&Y);

 44:   PetscNew(&user);
 45:   user->B = A;

 47:   MatCreateShell(PETSC_COMM_WORLD,2,2,2,2,user,&S);
 48:   MatShellSetOperation(S,MATOP_GET_DIAGONAL,(void (*)(void))MatGetDiagonal_User);
 49:   MatSetUp(S);

 51:   MatShift(S,42);
 52:   MatGetDiagonal(S,Y);
 53:   VecView(Y,PETSC_VIEWER_STDOUT_WORLD);
 54:   MatDiagonalSet(S,X,ADD_VALUES);
 55:   MatGetDiagonal(S,Y);
 56:   VecView(Y,PETSC_VIEWER_STDOUT_WORLD);
 57:   MatScale(S,42);
 58:   MatGetDiagonal(S,Y);
 59:   VecView(Y,PETSC_VIEWER_STDOUT_WORLD);

 61:   MatDestroy(&A);
 62:   MatDestroy(&S);
 63:   VecDestroy(&X);
 64:   VecDestroy(&Y);
 65:   PetscFree(user);
 66:   PetscFinalize();
 67:   return ierr;
 68: }

 70: /*TEST

 72:    test:
 73:       args: -malloc_dump

 75: TEST*/