1 2 static char help[] = "Tests incorrect use of MatDiagonalSet() for SHELL matrices\n\n"; 3 4 #include <petscmat.h> 5 6 typedef struct _n_User *User; 7 struct _n_User { 8 Mat B; 9 }; 10 11 static PetscErrorCode MatGetDiagonal_User(Mat A,Vec X) 12 { 13 User user; 14 15 PetscFunctionBegin; 16 CHKERRQ(MatShellGetContext(A,&user)); 17 CHKERRQ(MatGetDiagonal(user->B,X)); 18 PetscFunctionReturn(0); 19 } 20 21 int main(int argc,char **args) 22 { 23 const PetscScalar xvals[] = {11,13}; 24 const PetscInt inds[] = {0,1}; 25 PetscScalar avals[] = {2,3,5,7}; 26 Mat A,S; 27 Vec X,Y; 28 User user; 29 PetscErrorCode ierr; 30 31 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 32 CHKERRQ(MatCreateSeqAIJ(PETSC_COMM_WORLD,2,2,2,NULL,&A)); 33 CHKERRQ(MatSetUp(A)); 34 CHKERRQ(MatSetValues(A,2,inds,2,inds,avals,INSERT_VALUES)); 35 CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 36 CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 37 CHKERRQ(VecCreateSeq(PETSC_COMM_WORLD,2,&X)); 38 CHKERRQ(VecSetValues(X,2,inds,xvals,INSERT_VALUES)); 39 CHKERRQ(VecAssemblyBegin(X)); 40 CHKERRQ(VecAssemblyEnd(X)); 41 CHKERRQ(VecDuplicate(X,&Y)); 42 43 CHKERRQ(PetscNew(&user)); 44 user->B = A; 45 46 CHKERRQ(MatCreateShell(PETSC_COMM_WORLD,2,2,2,2,user,&S)); 47 CHKERRQ(MatShellSetOperation(S,MATOP_GET_DIAGONAL,(void (*)(void))MatGetDiagonal_User)); 48 CHKERRQ(MatSetUp(S)); 49 50 CHKERRQ(MatShift(S,42)); 51 CHKERRQ(MatGetDiagonal(S,Y)); 52 CHKERRQ(VecView(Y,PETSC_VIEWER_STDOUT_WORLD)); 53 CHKERRQ(MatDiagonalSet(S,X,ADD_VALUES)); 54 CHKERRQ(MatGetDiagonal(S,Y)); 55 CHKERRQ(VecView(Y,PETSC_VIEWER_STDOUT_WORLD)); 56 CHKERRQ(MatScale(S,42)); 57 CHKERRQ(MatGetDiagonal(S,Y)); 58 CHKERRQ(VecView(Y,PETSC_VIEWER_STDOUT_WORLD)); 59 60 CHKERRQ(MatDestroy(&A)); 61 CHKERRQ(MatDestroy(&S)); 62 CHKERRQ(VecDestroy(&X)); 63 CHKERRQ(VecDestroy(&Y)); 64 CHKERRQ(PetscFree(user)); 65 ierr = PetscFinalize(); 66 return ierr; 67 } 68 69 /*TEST 70 71 test: 72 args: -malloc_dump 73 74 TEST*/ 75