1 static char help[] = "Tests MatComputeOperator() and MatComputeOperatorTranspose()\n\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc,char **argv) 6 { 7 Mat A,Ae,Aet; 8 char filename[PETSC_MAX_PATH_LEN]; 9 char expltype[128],*etype = NULL; 10 PetscInt bs = 1; 11 PetscBool flg, check = PETSC_TRUE; 12 13 PetscCall(PetscInitialize(&argc,&argv,(char*) 0,help)); 14 15 PetscCall(PetscOptionsGetString(NULL,NULL,"-expl_type",expltype,sizeof(expltype),&flg)); 16 if (flg) { 17 PetscCall(PetscStrallocpy(expltype,&etype)); 18 } 19 PetscCall(PetscOptionsGetString(NULL,NULL,"-f",filename,sizeof(filename),&flg)); 20 PetscCall(PetscOptionsGetInt(NULL,NULL,"-bs",&bs,NULL)); 21 if (!flg) { 22 PetscInt M = 13,N = 6; 23 24 PetscCall(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 25 PetscCall(PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL)); 26 PetscCall(MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&A)); 27 PetscCall(MatSetBlockSize(A,bs)); 28 PetscCall(MatSetRandom(A,NULL)); 29 } else { 30 PetscViewer viewer; 31 32 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer)); 33 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 34 PetscCall(MatSetBlockSize(A,bs)); 35 PetscCall(MatSetFromOptions(A)); 36 PetscCall(MatLoad(A,viewer)); 37 PetscCall(PetscViewerDestroy(&viewer)); 38 } 39 PetscCall(PetscObjectSetName((PetscObject)A,"Matrix")); 40 PetscCall(MatViewFromOptions(A,NULL,"-view_expl")); 41 42 PetscCall(MatComputeOperator(A,etype,&Ae)); 43 PetscCall(PetscObjectSetName((PetscObject)Ae,"Explicit matrix")); 44 PetscCall(MatViewFromOptions(Ae,NULL,"-view_expl")); 45 46 PetscCall(PetscOptionsGetBool(NULL,NULL,"-check",&check,NULL)); 47 if (check) { 48 Mat A2; 49 PetscReal err,tol = PETSC_SMALL; 50 51 PetscCall(PetscOptionsGetReal(NULL,NULL,"-tol",&tol,NULL)); 52 PetscCall(MatConvert(A,etype,MAT_INITIAL_MATRIX,&A2)); 53 PetscCall(MatAXPY(A2,-1.0,Ae,DIFFERENT_NONZERO_PATTERN)); 54 PetscCall(MatNorm(A2,NORM_FROBENIUS,&err)); 55 if (err > tol) { 56 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Error %g > %g (type %s)\n",(double)err,(double)tol,etype)); 57 } 58 PetscCall(MatDestroy(&A2)); 59 } 60 61 PetscCall(MatComputeOperatorTranspose(A,etype,&Aet)); 62 PetscCall(PetscObjectSetName((PetscObject)Aet,"Explicit matrix transpose")); 63 PetscCall(MatViewFromOptions(Aet,NULL,"-view_expl")); 64 65 PetscCall(PetscFree(etype)); 66 PetscCall(MatDestroy(&Ae)); 67 PetscCall(MatDestroy(&Aet)); 68 PetscCall(MatDestroy(&A)); 69 PetscCall(PetscFinalize()); 70 return 0; 71 } 72 73 /*TEST 74 75 test: 76 output_file: output/ex222_null.out 77 78 testset: 79 suffix: matexpl_rect 80 output_file: output/ex222_null.out 81 nsize: {{1 3}} 82 args: -expl_type {{dense aij baij}} 83 84 testset: 85 suffix: matexpl_square 86 output_file: output/ex222_null.out 87 nsize: {{1 3}} 88 args: -bs {{1 2 3}} -M 36 -N 36 -expl_type {{dense aij baij sbaij}} 89 90 TEST*/ 91