Actual source code: ex2.c
petsc-3.7.3 2016-08-01
2: static char help[] = "Tests MatTranspose(), MatNorm(), MatAXPY() and MatAYPX().\n\n";
4: #include <petscmat.h>
8: int main(int argc,char **argv)
9: {
10: Mat mat,tmat = 0;
11: PetscInt m = 7,n,i,j,rstart,rend,rect = 0;
13: PetscMPIInt size,rank;
14: PetscBool flg;
15: PetscScalar v, alpha;
16: PetscReal normf,normi,norm1;
18: PetscInitialize(&argc,&argv,(char*)0,help);
19: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON);
20: PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
22: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: n = m;
24: PetscOptionsHasName(NULL,NULL,"-rectA",&flg);
25: if (flg) {n += 2; rect = 1;}
26: PetscOptionsHasName(NULL,NULL,"-rectB",&flg);
27: if (flg) {n -= 2; rect = 1;}
29: /* ------- Assemble matrix, test MatValid() --------- */
30: MatCreate(PETSC_COMM_WORLD,&mat);
31: MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n);
32: MatSetFromOptions(mat);
33: MatSetUp(mat);
34: MatGetOwnershipRange(mat,&rstart,&rend);
35: for (i=rstart; i<rend; i++) {
36: for (j=0; j<n; j++) {
37: v = 10.0*i+j;
38: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
39: }
40: }
41: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
42: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
44: /* ----------------- Test MatNorm() ----------------- */
45: MatNorm(mat,NORM_FROBENIUS,&normf);
46: MatNorm(mat,NORM_1,&norm1);
47: MatNorm(mat,NORM_INFINITY,&normi);
48: PetscPrintf(PETSC_COMM_WORLD,"original A: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi);
49: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
51: /* --------------- Test MatTranspose() -------------- */
52: PetscOptionsHasName(NULL,NULL,"-in_place",&flg);
53: if (!rect && flg) {
54: MatTranspose(mat,MAT_REUSE_MATRIX,&mat); /* in-place transpose */
55: tmat = mat; mat = 0;
56: } else { /* out-of-place transpose */
57: MatTranspose(mat,MAT_INITIAL_MATRIX,&tmat);
58: }
60: /* ----------------- Test MatNorm() ----------------- */
61: /* Print info about transpose matrix */
62: MatNorm(tmat,NORM_FROBENIUS,&normf);
63: MatNorm(tmat,NORM_1,&norm1);
64: MatNorm(tmat,NORM_INFINITY,&normi);
65: PetscPrintf(PETSC_COMM_WORLD,"B = A^T: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",(double)normf,(double)norm1,(double)normi);
66: MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
68: /* ----------------- Test MatAXPY(), MatAYPX() ----------------- */
69: if (mat && !rect) {
70: alpha = 1.0;
71: PetscOptionsGetScalar(NULL,NULL,"-alpha",&alpha,NULL);
72: PetscPrintf(PETSC_COMM_WORLD,"MatAXPY: B = B + alpha * A\n");
73: MatAXPY(tmat,alpha,mat,DIFFERENT_NONZERO_PATTERN);
74: MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
76: PetscPrintf(PETSC_COMM_WORLD,"MatAYPX: B = alpha*B + A\n");
77: MatAYPX(tmat,alpha,mat,DIFFERENT_NONZERO_PATTERN);
78: MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
79: }
81: {
82: Mat C;
83: alpha = 1.0;
84: PetscPrintf(PETSC_COMM_WORLD,"MatAXPY: C = C + alpha * A, C=A, SAME_NONZERO_PATTERN\n");
85: MatDuplicate(mat,MAT_COPY_VALUES,&C);
86: MatAXPY(C,alpha,mat,SAME_NONZERO_PATTERN);
87: MatView(C,PETSC_VIEWER_STDOUT_WORLD);
88: MatDestroy(&C);
89: }
91: {
92: Mat matB;
93: /* get matB that has nonzeros of mat in all even numbers of row and col */
94: MatCreate(PETSC_COMM_WORLD,&matB);
95: MatSetSizes(matB,PETSC_DECIDE,PETSC_DECIDE,m,n);
96: MatSetFromOptions(matB);
97: MatSetUp(matB);
98: MatGetOwnershipRange(matB,&rstart,&rend);
99: if (rstart % 2 != 0) rstart++;
100: for (i=rstart; i<rend; i += 2) {
101: for (j=0; j<n; j += 2) {
102: v = 10.0*i+j;
103: MatSetValues(matB,1,&i,1,&j,&v,INSERT_VALUES);
104: }
105: }
106: MatAssemblyBegin(matB,MAT_FINAL_ASSEMBLY);
107: MatAssemblyEnd(matB,MAT_FINAL_ASSEMBLY);
108: PetscPrintf(PETSC_COMM_WORLD," A: original matrix:\n");
109: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
110: PetscPrintf(PETSC_COMM_WORLD," B(a subset of A):\n");
111: MatView(matB,PETSC_VIEWER_STDOUT_WORLD);
112: PetscPrintf(PETSC_COMM_WORLD,"MatAXPY: B = B + alpha * A, SUBSET_NONZERO_PATTERN\n");
113: MatAXPY(mat,alpha,matB,SUBSET_NONZERO_PATTERN);
114: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
115: MatDestroy(&matB);
116: }
118: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD);
119: /* Free data structures */
120: MatDestroy(&mat);
121: MatDestroy(&tmat);
122: PetscFinalize();
123: return 0;
124: }