Actual source code: ex49.c
1: /*$Id: ex49.c,v 1.24 2001/08/07 21:30:08 bsmith Exp $*/
3: static char help[] = "Tests MatTranspose(), MatNorm(), MatValid(), and MatAXPY().\n\n";
5: #include petscmat.h
9: int main(int argc,char **argv)
10: {
11: Mat mat,tmat = 0;
12: int m = 4,n,i,j,ierr,size,rank;
13: int rstart,rend,rect = 0,nd,bs,*diag,*bdlen;
14: PetscTruth flg,isbdiag;
15: PetscScalar v,**diagv;
16: PetscReal normf,normi,norm1;
17: MatInfo info;
18:
19: PetscInitialize(&argc,&argv,(char*)0,help);
20: PetscOptionsGetInt(PETSC_NULL,"-m",&m,PETSC_NULL);
21: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
22: MPI_Comm_size(PETSC_COMM_WORLD,&size);
23: n = m;
24: PetscOptionsHasName(PETSC_NULL,"-rect1",&flg);
25: if (flg) {n += 2; rect = 1;}
26: PetscOptionsHasName(PETSC_NULL,"-rect2",&flg);
27: if (flg) {n -= 2; rect = 1;}
29: /* Create and assemble matrix */
30: MatCreate(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,&mat);
31: MatSetFromOptions(mat);
32: MatGetOwnershipRange(mat,&rstart,&rend);
33: for (i=rstart; i<rend; i++) {
34: for (j=0; j<n; j++) {
35: v=10*i+j;
36: MatSetValues(mat,1,&i,1,&j,&v,INSERT_VALUES);
37: }
38: }
39: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
40: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
42: /* Test whether matrix has been corrupted (just to demonstrate this
43: routine) not needed in most application codes. */
44: MatValid(mat,(PetscTruth*)&flg);
45: if (!flg) SETERRQ(1,"Corrupted matrix.");
47: /* Print info about original matrix */
48: MatGetInfo(mat,MAT_GLOBAL_SUM,&info);
49: PetscPrintf(PETSC_COMM_WORLD,"original matrix nonzeros = %d, allocated nonzeros = %d\n",
50: (int)info.nz_used,(int)info.nz_allocated);
51: MatNorm(mat,NORM_FROBENIUS,&normf);
52: MatNorm(mat,NORM_1,&norm1);
53: MatNorm(mat,NORM_INFINITY,&normi);
54: PetscPrintf(PETSC_COMM_WORLD,"original: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",
55: normf,norm1,normi);
56: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
58: PetscTypeCompare((PetscObject)mat,MATSEQBDIAG,&isbdiag);
59: if (!isbdiag) {
60: PetscTypeCompare((PetscObject)mat,MATMPIBDIAG,&isbdiag);
61: }
62: if (isbdiag) {
63: MatBDiagGetData(mat,&nd,&bs,&diag,&bdlen,&diagv);
64: PetscPrintf(PETSC_COMM_WORLD,"original matrix: block diag format: %d diagonals, block size = %d\n",nd,bs);
65: for (i=0; i<nd; i++) {
66: PetscPrintf(PETSC_COMM_WORLD," diag=%d, bdlength=%d\n",diag[i],bdlen[i]);
67: }
68: }
70: /* Form matrix transpose */
71: PetscOptionsHasName(PETSC_NULL,"-in_place",&flg);
72: if (!rect && flg) {
73: MatTranspose(mat,0); /* in-place transpose */
74: tmat = mat; mat = 0;
75: } else { /* out-of-place transpose */
76: MatTranspose(mat,&tmat);
77: }
79: /* Print info about transpose matrix */
80: MatGetInfo(tmat,MAT_GLOBAL_SUM,&info);
81: PetscPrintf(PETSC_COMM_WORLD,"transpose matrix nonzeros = %d, allocated nonzeros = %d\n",
82: (int)info.nz_used,(int)info.nz_allocated);
83: MatNorm(tmat,NORM_FROBENIUS,&normf);
84: MatNorm(tmat,NORM_1,&norm1);
85: MatNorm(tmat,NORM_INFINITY,&normi);
86: PetscPrintf(PETSC_COMM_WORLD,"transpose: Frobenious norm = %g, one norm = %g, infinity norm = %g\n",
87: normf,norm1,normi);
88: MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
90: if (isbdiag) {
91: MatBDiagGetData(tmat,&nd,&bs,&diag,&bdlen,&diagv);
92: PetscPrintf(PETSC_COMM_WORLD,"transpose matrix: block diag format: %d diagonals, block size = %d\n",nd,bs);
93: for (i=0; i<nd; i++) {
94: PetscPrintf(PETSC_COMM_WORLD," diag=%d, bdlength=%d\n",diag[i],bdlen[i]);
95: }
96: }
98: /* Test MatAXPY */
99: if (mat && !rect) {
100: PetscScalar alpha = 1.0;
101: PetscOptionsGetScalar(PETSC_NULL,"-alpha",&alpha,PETSC_NULL);
102: PetscPrintf(PETSC_COMM_WORLD,"matrix addition: B = B + alpha * A\n");
103: MatAXPY(&alpha,mat,tmat,DIFFERENT_NONZERO_PATTERN);
104: MatView(tmat,PETSC_VIEWER_STDOUT_WORLD);
105: }
107: /* Free data structures */
108: MatDestroy(tmat);
109: if (mat) {MatDestroy(mat);}
111: PetscFinalize();
112: return 0;
113: }
114: