Actual source code: ex47.c
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";
6: #include <petscmat.h>
8: int main(int argc,char **args)
9: {
10: Mat A,B,C;
11: PetscViewer va,vb,vc;
12: Vec x,y;
13: PetscInt i,j,row,m,n,ncols1,ncols2,ct,m2,n2;
14: const PetscInt *cols1,*cols2;
15: char file[PETSC_MAX_PATH_LEN];
16: PetscBool tflg;
17: PetscScalar rval;
18: const PetscScalar *vals1,*vals2;
19: PetscReal norm1,norm2,rnorm;
20: PetscRandom r;
22: PetscInitialize(&argc,&args,(char*)0,help);
23: PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),NULL);
25: /* Load the matrix as AIJ format */
26: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&va);
27: MatCreate(PETSC_COMM_WORLD,&A);
28: MatSetType(A,MATSEQAIJ);
29: MatLoad(A,va);
30: PetscViewerDestroy(&va);
32: /* Load the matrix as BAIJ format */
33: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&vb);
34: MatCreate(PETSC_COMM_WORLD,&B);
35: MatSetType(B,MATSEQBAIJ);
36: MatLoad(B,vb);
37: PetscViewerDestroy(&vb);
39: /* Load the matrix as BAIJ format */
40: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&vc);
41: MatCreate(PETSC_COMM_WORLD,&C);
42: MatSetType(C,MATSEQBAIJ);
43: MatSetFromOptions(C);
44: MatLoad(C,vc);
45: PetscViewerDestroy(&vc);
47: MatGetSize(A,&m,&n);
48: MatGetSize(B,&m2,&n2);
51: /* Test MatEqual() */
52: MatEqual(B,C,&tflg);
55: /* Test MatGetDiagonal() */
56: VecCreateSeq(PETSC_COMM_SELF,m,&x);
57: VecCreateSeq(PETSC_COMM_SELF,m,&y);
59: MatGetDiagonal(A,x);
60: MatGetDiagonal(B,y);
62: VecEqual(x,y,&tflg);
65: /* Test MatDiagonalScale() */
66: PetscRandomCreate(PETSC_COMM_SELF,&r);
67: PetscRandomSetFromOptions(r);
68: VecSetRandom(x,r);
69: VecSetRandom(y,r);
71: MatDiagonalScale(A,x,y);
72: MatDiagonalScale(B,x,y);
73: MatMult(A,x,y);
74: VecNorm(y,NORM_2,&norm1);
75: MatMult(B,x,y);
76: VecNorm(y,NORM_2,&norm2);
77: rnorm = ((norm1-norm2)*100)/norm1;
80: /* Test MatGetRow()/ MatRestoreRow() */
81: for (ct=0; ct<100; ct++) {
82: PetscRandomGetValue(r,&rval);
83: row = (int)(rval*m);
84: MatGetRow(A,row,&ncols1,&cols1,&vals1);
85: MatGetRow(B,row,&ncols2,&cols2,&vals2);
87: for (i=0,j=0; i<ncols1 && j<ncols2; i++) {
88: while (cols2[j] != cols1[i]) j++;
90: }
93: MatRestoreRow(A,row,&ncols1,&cols1,&vals1);
94: MatRestoreRow(B,row,&ncols2,&cols2,&vals2);
95: }
97: MatDestroy(&A);
98: MatDestroy(&B);
99: MatDestroy(&C);
100: VecDestroy(&x);
101: VecDestroy(&y);
102: PetscRandomDestroy(&r);
103: PetscFinalize();
104: return 0;
105: }
107: /*TEST
109: build:
110: requires: !complex
112: test:
113: args: -f ${DATAFILESPATH}/matrices/cfd.1.10 -mat_block_size 5
114: requires: !complex double datafilespath !defined(PETSC_USE_64BIT_INDICES)
116: TEST*/