Actual source code: ex100.c
2: static char help[] = "Tests vatious routines in MatMAIJ format.\n";
4: #include petscmat.h
5: #define IMAX 15
8: int main(int argc,char **args)
9: {
10: Mat A,B,MA;
11: PetscViewer fd;
12: char file[PETSC_MAX_PATH_LEN];
13: PetscRandom rand;
14: Vec xx,yy,s1,s2;
15: PetscInt m,n,M,N,dof=1;
16: PetscMPIInt rank,size;
17: PetscErrorCode ierr;
18: PetscTruth flg;
20: PetscInitialize(&argc,&args,(char *)0,help);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
22: MPI_Comm_size(PETSC_COMM_WORLD,&size);
24: #if defined(PETSC_USE_COMPLEX)
25: SETERRQ(1,"This example does not work with complex numbers");
26: #else
28: /* Load aij matrix A */
29: PetscOptionsGetString(PETSC_NULL,"-f",file,PETSC_MAX_PATH_LEN-1,&flg);
30: if (!flg) SETERRQ(PETSC_ERR_USER,"Must indicate binary file with the -f option");
31: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
32: MatLoad(fd,MATAIJ,&A);
33: PetscViewerDestroy(fd);
35: /* Get dof, then create maij matrix MA */
36: PetscOptionsGetInt(PETSC_NULL,"-dof",&dof,PETSC_NULL);
37: MatCreateMAIJ(A,dof,&MA);
38: MatGetLocalSize(MA,&m,&n);
39: MatGetSize(MA,&M,&N);
41: if (size == 1){
42: MatConvert(MA,MATSEQAIJ,MAT_INITIAL_MATRIX,&B);
43: } else {
44: MatConvert(MA,MATMPIAIJ,MAT_INITIAL_MATRIX,&B);
45: }
46: PetscRandomCreate(PETSC_COMM_WORLD,RANDOM_DEFAULT,&rand);
47: VecCreate(PETSC_COMM_WORLD,&xx);
48: VecSetSizes(xx,m,PETSC_DECIDE);
49: VecSetFromOptions(xx);
50: VecDuplicate(xx,&s1);
51: VecDuplicate(xx,&s2);
52: VecDuplicate(xx,&yy);
54: /* Test MatMult() */
55: MatMultEqual(MA,B,10,&flg);
56: if (!flg){
57: SETERRQ(PETSC_ERR_CONV_FAILED,"Error: MatMul() for MAIJ matrix");
58: }
59: /* Test MatMultAdd() */
60: MatMultAddEqual(MA,B,10,&flg);
61: if (!flg){
62: SETERRQ(PETSC_ERR_CONV_FAILED,"Error: MatMulAdd() for MAIJ matrix");
63: }
65: /* Test MatMultTranspose() */
66: MatMultTransposeEqual(MA,B,10,&flg);
67: if (!flg){
68: SETERRQ(PETSC_ERR_CONV_FAILED,"Error: MatMulAdd() for MAIJ matrix");
69: }
70:
71: /* Test MatMultTransposeAdd() */
72: MatMultTransposeAddEqual(MA,B,10,&flg);
73: if (!flg){
74: SETERRQ(PETSC_ERR_CONV_FAILED,"Error: MatMulTransposeAdd() for MAIJ matrix");
75: }
77: MatDestroy(MA);
78: MatDestroy(A);
79: MatDestroy(B);
80: VecDestroy(xx);
81: VecDestroy(yy);
82: VecDestroy(s1);
83: VecDestroy(s2);
84: PetscRandomDestroy(rand);
85: PetscFinalize();
86: #endif
87: return 0;
88: }