1 2 static char help[] = "Tests MatShellTestMult()\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 MatMult_User(Mat A,Vec X,Vec Y) 12 { 13 User user; 14 15 PetscFunctionBegin; 16 CHKERRQ(MatShellGetContext(A,&user)); 17 CHKERRQ(MatMult(user->B,X,Y)); 18 PetscFunctionReturn(0); 19 } 20 21 static PetscErrorCode MatMultTranspose_User(Mat A,Vec X,Vec Y) 22 { 23 User user; 24 25 PetscFunctionBegin; 26 CHKERRQ(MatShellGetContext(A,&user)); 27 CHKERRQ(MatMultTranspose(user->B,X,Y)); 28 PetscFunctionReturn(0); 29 } 30 31 static PetscErrorCode MyFunction(void *ctx,Vec x,Vec y) 32 { 33 User user = (User) ctx; 34 35 PetscFunctionBegin; 36 CHKERRQ(MatMult(user->B,x,y)); 37 PetscFunctionReturn(0); 38 } 39 40 int main(int argc,char **args) 41 { 42 const PetscInt inds[] = {0,1}; 43 PetscScalar avals[] = {2,3,5,7}; 44 Mat S; 45 User user; 46 PetscErrorCode ierr; 47 Vec base; 48 49 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 50 CHKERRQ(PetscNew(&user)); 51 CHKERRQ(MatCreateSeqAIJ(PETSC_COMM_WORLD,2,2,2,NULL,&user->B)); 52 CHKERRQ(MatSetUp(user->B)); 53 CHKERRQ(MatSetValues(user->B,2,inds,2,inds,avals,INSERT_VALUES)); 54 CHKERRQ(MatAssemblyBegin(user->B,MAT_FINAL_ASSEMBLY)); 55 CHKERRQ(MatAssemblyEnd(user->B,MAT_FINAL_ASSEMBLY)); 56 CHKERRQ(MatCreateVecs(user->B,&base,NULL)); 57 CHKERRQ(MatCreateShell(PETSC_COMM_WORLD,2,2,2,2,user,&S)); 58 CHKERRQ(MatSetUp(S)); 59 CHKERRQ(MatShellSetOperation(S,MATOP_MULT,(void (*)(void))MatMult_User)); 60 CHKERRQ(MatShellSetOperation(S,MATOP_MULT_TRANSPOSE,(void (*)(void))MatMultTranspose_User)); 61 62 CHKERRQ(MatShellTestMult(S,MyFunction,base,user,NULL)); 63 CHKERRQ(MatShellTestMultTranspose(S,MyFunction,base,user,NULL)); 64 65 CHKERRQ(VecDestroy(&base)); 66 CHKERRQ(MatDestroy(&user->B)); 67 CHKERRQ(MatDestroy(&S)); 68 CHKERRQ(PetscFree(user)); 69 ierr = PetscFinalize(); 70 return ierr; 71 } 72 73 /*TEST 74 75 test: 76 args: -mat_shell_test_mult_view -mat_shell_test_mult_transpose_view 77 78 TEST*/ 79