Actual source code: ex177.c
petsc-3.12.2 2019-11-22
2: static char help[] = "Tests various routines in MatKAIJ format.\n";
4: #include <petscmat.h>
5: #define IMAX 15
7: int main(int argc,char **args)
8: {
9: Mat A,B,TA;
10: PetscScalar *S,*T;
11: PetscViewer fd;
12: char file[PETSC_MAX_PATH_LEN];
13: PetscInt m,n,M,N,p=1,q=1,i,j;
14: PetscMPIInt rank,size;
16: PetscBool flg;
18: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
19: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
20: MPI_Comm_size(PETSC_COMM_WORLD,&size);
22: /* Load aij matrix A */
23: PetscOptionsGetString(NULL,NULL,"-f",file,PETSC_MAX_PATH_LEN,&flg);
24: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_USER,"Must indicate binary file with the -f option");
25: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);
26: MatCreate(PETSC_COMM_WORLD,&A);
27: MatLoad(A,fd);
28: PetscViewerDestroy(&fd);
30: /* Get dof, then create S and T */
31: PetscOptionsGetInt(NULL,NULL,"-p",&p,NULL);
32: PetscOptionsGetInt(NULL,NULL,"-q",&q,NULL);
33: PetscMalloc2(p*q,&S,p*q,&T);
34: for (i=0; i<p*q; i++) S[i] = 0;
36: for (i=0; i<p; i++) {
37: for (j = 0; j<q; j++) {
38: /* set some random non-zero values */
39: S[i+p*j] = ((PetscReal) ((i+1)*(j+1))) / ((PetscReal) (p+q));
40: T[i+p*j] = ((PetscReal) ((p-i)+j)) / ((PetscReal) (p*q));
41: }
42: }
44: /* Test KAIJ when both S & T are not NULL */
46: /* create kaij matrix TA */
47: MatCreateKAIJ(A,p,q,S,T,&TA);
48: MatGetLocalSize(A,&m,&n);
49: MatGetSize(A,&M,&N);
51: MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B);
53: /* Test MatMult() */
54: MatMultEqual(TA,B,10,&flg);
55: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 1: MatMult() for KAIJ matrix");
56: /* Test MatMultAdd() */
57: MatMultAddEqual(TA,B,10,&flg);
58: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 1: MatMultAdd() for KAIJ matrix");
60: MatDestroy(&TA);
61: MatDestroy(&B);
63: /* Test KAIJ when S is NULL */
65: /* create kaij matrix TA */
66: MatCreateKAIJ(A,p,q,NULL,T,&TA);
67: MatGetLocalSize(A,&m,&n);
68: MatGetSize(A,&M,&N);
70: MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B);
72: /* Test MatMult() */
73: MatMultEqual(TA,B,10,&flg);
74: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 2: MatMult() for KAIJ matrix");
75: /* Test MatMultAdd() */
76: MatMultAddEqual(TA,B,10,&flg);
77: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 2: MatMultAdd() for KAIJ matrix");
79: MatDestroy(&TA);
80: MatDestroy(&B);
81:
82: /* Test KAIJ when T is NULL */
84: /* create kaij matrix TA */
85: MatCreateKAIJ(A,p,q,S,NULL,&TA);
86: MatGetLocalSize(A,&m,&n);
87: MatGetSize(A,&M,&N);
89: MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B);
91: /* Test MatMult() */
92: MatMultEqual(TA,B,10,&flg);
93: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 3: MatMult() for KAIJ matrix");
94: /* Test MatMultAdd() */
95: MatMultAddEqual(TA,B,10,&flg);
96: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 3: MatMultAdd() for KAIJ matrix");
98: MatDestroy(&TA);
99: MatDestroy(&B);
101: /* Test KAIJ when T is is an identity matrix */
103: if (p == q) {
104: for (i=0; i<p; i++) {
105: for (j=0; j<q; j++) {
106: if (i==j) T[i+j*p] = 1.0;
107: else T[i+j*p] = 0.0;
108: }
109: }
111: /* create kaij matrix TA */
112: MatCreateKAIJ(A,p,q,S,T,&TA);
113: MatGetLocalSize(A,&m,&n);
114: MatGetSize(A,&M,&N);
116: MatConvert(TA,MATAIJ,MAT_INITIAL_MATRIX,&B);
118: /* Test MatMult() */
119: MatMultEqual(TA,B,10,&flg);
120: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 4: MatMult() for KAIJ matrix");
121: /* Test MatMultAdd() */
122: MatMultAddEqual(TA,B,10,&flg);
123: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_CONV_FAILED,"Error in Test 4: MatMultAdd() for KAIJ matrix");
125: MatDestroy(&TA);
126: MatDestroy(&B);
127: }
129: /* Done with all tests */
131: PetscFree2(S,T);
132: MatDestroy(&A);
133: PetscFinalize();
134: return ierr;
135: }
137: /*TEST
139: build:
140: requires: !complex
142: test:
143: requires: datafilespath !complex double !define(PETSC_USE_64BIT_INDICES)
144: output_file: output/ex176.out
145: nsize: {{1 2 3 4}}
146: args: -f ${DATAFILESPATH}/matrices/small -p {{2 3 7}} -q {{3 7}} -viewer_binary_skip_info
148: TEST*/