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