1 2 static char help[] = "Tests the various routines in MatBAIJ format.\n\ 3 Input arguments are:\n\ 4 -f <input_file> : file to load. For example see $PETSC_DIR/share/petsc/datafiles/matrices\n\n"; 5 6 #include <petscmat.h> 7 8 int main(int argc,char **args) 9 { 10 Mat A,B,C; 11 PetscViewer va,vb,vc; 12 Vec x,y; 13 PetscErrorCode ierr; 14 PetscInt i,j,row,m,n,ncols1,ncols2,ct,m2,n2; 15 const PetscInt *cols1,*cols2; 16 char file[PETSC_MAX_PATH_LEN]; 17 PetscBool tflg; 18 PetscScalar rval; 19 const PetscScalar *vals1,*vals2; 20 PetscReal norm1,norm2,rnorm; 21 PetscRandom r; 22 23 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 24 ierr = PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),NULL);CHKERRQ(ierr); 25 26 /* Load the matrix as AIJ format */ 27 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&va);CHKERRQ(ierr); 28 ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 29 ierr = MatSetType(A,MATSEQAIJ);CHKERRQ(ierr); 30 ierr = MatLoad(A,va);CHKERRQ(ierr); 31 ierr = PetscViewerDestroy(&va);CHKERRQ(ierr); 32 33 /* Load the matrix as BAIJ format */ 34 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&vb);CHKERRQ(ierr); 35 ierr = MatCreate(PETSC_COMM_WORLD,&B);CHKERRQ(ierr); 36 ierr = MatSetType(B,MATSEQBAIJ);CHKERRQ(ierr); 37 ierr = MatLoad(B,vb);CHKERRQ(ierr); 38 ierr = PetscViewerDestroy(&vb);CHKERRQ(ierr); 39 40 /* Load the matrix as BAIJ format */ 41 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&vc);CHKERRQ(ierr); 42 ierr = MatCreate(PETSC_COMM_WORLD,&C);CHKERRQ(ierr); 43 ierr = MatSetType(C,MATSEQBAIJ);CHKERRQ(ierr); 44 ierr = MatSetFromOptions(C);CHKERRQ(ierr); 45 ierr = MatLoad(C,vc);CHKERRQ(ierr); 46 ierr = PetscViewerDestroy(&vc);CHKERRQ(ierr); 47 48 ierr = MatGetSize(A,&m,&n);CHKERRQ(ierr); 49 ierr = MatGetSize(B,&m2,&n2);CHKERRQ(ierr); 50 if (m!=m2) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ,"Matrices are of different size. Cannot run this example"); 51 52 /* Test MatEqual() */ 53 ierr = MatEqual(B,C,&tflg);CHKERRQ(ierr); 54 if (!tflg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatEqual() failed"); 55 56 /* Test MatGetDiagonal() */ 57 ierr = VecCreateSeq(PETSC_COMM_SELF,m,&x);CHKERRQ(ierr); 58 ierr = VecCreateSeq(PETSC_COMM_SELF,m,&y);CHKERRQ(ierr); 59 60 ierr = MatGetDiagonal(A,x);CHKERRQ(ierr); 61 ierr = MatGetDiagonal(B,y);CHKERRQ(ierr); 62 63 ierr = VecEqual(x,y,&tflg);CHKERRQ(ierr); 64 if (!tflg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatGetDiagonal() failed"); 65 66 /* Test MatDiagonalScale() */ 67 ierr = PetscRandomCreate(PETSC_COMM_SELF,&r);CHKERRQ(ierr); 68 ierr = PetscRandomSetFromOptions(r);CHKERRQ(ierr); 69 ierr = VecSetRandom(x,r);CHKERRQ(ierr); 70 ierr = VecSetRandom(y,r);CHKERRQ(ierr); 71 72 ierr = MatDiagonalScale(A,x,y);CHKERRQ(ierr); 73 ierr = MatDiagonalScale(B,x,y);CHKERRQ(ierr); 74 ierr = MatMult(A,x,y);CHKERRQ(ierr); 75 ierr = VecNorm(y,NORM_2,&norm1);CHKERRQ(ierr); 76 ierr = MatMult(B,x,y);CHKERRQ(ierr); 77 ierr = VecNorm(y,NORM_2,&norm2);CHKERRQ(ierr); 78 rnorm = ((norm1-norm2)*100)/norm1; 79 if (rnorm<-0.1 || rnorm>0.01) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatDiagonalScale() failed Norm1 %g Norm2 %g",(double)norm1,(double)norm2); 80 81 /* Test MatGetRow()/ MatRestoreRow() */ 82 for (ct=0; ct<100; ct++) { 83 ierr = PetscRandomGetValue(r,&rval);CHKERRQ(ierr); 84 row = (int)(rval*m); 85 ierr = MatGetRow(A,row,&ncols1,&cols1,&vals1);CHKERRQ(ierr); 86 ierr = MatGetRow(B,row,&ncols2,&cols2,&vals2);CHKERRQ(ierr); 87 88 for (i=0,j=0; i<ncols1 && j<ncols2; i++) { 89 while (cols2[j] != cols1[i]) j++; 90 if (vals1[i] != vals2[j]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatGetRow() failed - vals incorrect."); 91 } 92 if (i<ncols1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"MatGetRow() failed - cols incorrect"); 93 94 ierr = MatRestoreRow(A,row,&ncols1,&cols1,&vals1);CHKERRQ(ierr); 95 ierr = MatRestoreRow(B,row,&ncols2,&cols2,&vals2);CHKERRQ(ierr); 96 } 97 98 ierr = MatDestroy(&A);CHKERRQ(ierr); 99 ierr = MatDestroy(&B);CHKERRQ(ierr); 100 ierr = MatDestroy(&C);CHKERRQ(ierr); 101 ierr = VecDestroy(&x);CHKERRQ(ierr); 102 ierr = VecDestroy(&y);CHKERRQ(ierr); 103 ierr = PetscRandomDestroy(&r);CHKERRQ(ierr); 104 ierr = PetscFinalize(); 105 return ierr; 106 } 107 108 /*TEST 109 110 build: 111 requires: !complex 112 113 test: 114 args: -f ${DATAFILESPATH}/matrices/cfd.1.10 -mat_block_size 5 115 requires: !complex double datafilespath !define(PETSC_USE_64BIT_INDICES) 116 117 TEST*/ 118