Actual source code: ex100.c
petsc-3.3-p6 2013-02-11
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: PetscInt m,n,M,N,dof=1;
14: PetscMPIInt rank,size;
15: PetscErrorCode ierr;
16: PetscBool flg;
18: PetscInitialize(&argc,&args,(char *)0,help);
19: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
22: #if defined(PETSC_USE_COMPLEX)
23: SETERRQ(PETSC_COMM_WORLD,1,"This example does not work with complex numbers");
24: #else
26: /* Load aij matrix A */
27: PetscOptionsGetString(PETSC_NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);
28: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option");
29: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
30: MatCreate(PETSC_COMM_WORLD,&A);
31: MatLoad(A,fd);
32: PetscViewerDestroy(&fd);
34: /* Get dof, then create maij matrix MA */
35: PetscOptionsGetInt(PETSC_NULL,"-dof",&dof,PETSC_NULL);
36: MatCreateMAIJ(A,dof,&MA);
37: MatGetLocalSize(MA,&m,&n);
38: MatGetSize(MA,&M,&N);
39:
40: if (size == 1){
41: MatConvert(MA,MATSEQAIJ,MAT_INITIAL_MATRIX,&B);
42: } else {
43: MatConvert(MA,MATMPIAIJ,MAT_INITIAL_MATRIX,&B);
44: }
46: /* Test MatMult() */
47: MatMultEqual(MA,B,10,&flg);
48: if (!flg){
49: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error: MatMul() for MAIJ matrix");
50: }
51: /* Test MatMultAdd() */
52: MatMultAddEqual(MA,B,10,&flg);
53: if (!flg){
54: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error: MatMulAdd() for MAIJ matrix");
55: }
57: /* Test MatMultTranspose() */
58: MatMultTransposeEqual(MA,B,10,&flg);
59: if (!flg){
60: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error: MatMulAdd() for MAIJ matrix");
61: }
62:
63: /* Test MatMultTransposeAdd() */
64: MatMultTransposeAddEqual(MA,B,10,&flg);
65: if (!flg){
66: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error: MatMulTransposeAdd() for MAIJ matrix");
67: }
69: MatDestroy(&MA);
70: MatDestroy(&A);
71: MatDestroy(&B);
72: PetscFinalize();
73: #endif
74: return 0;
75: }