Actual source code: ex4.c
1: /*$Id: ex4.c,v 1.23 2001/08/07 21:30:08 bsmith Exp $*/
3: static char help[] = "Creates a matrix, inserts some values, and tests MatGetSubMatrices() and MatZeroEntries().\n\n";
5: #include petscmat.h
9: int main(int argc,char **argv)
10: {
11: Mat mat,submat,*submatrices;
12: int m = 10,n = 10,i = 4,tmp,ierr;
13: IS irkeep,ickeep;
14: PetscScalar value = 1.0;
15: PetscViewer sviewer;
17: PetscInitialize(&argc,&argv,(char *)0,help);
18: PetscViewerSetFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON);
19: PetscViewerSetFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_COMMON);
21: MatCreate(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,&mat);
22: MatSetFromOptions(mat);
23: for (i=0; i<m; i++) {
24: value = (PetscReal)i+1; tmp = i % 5;
25: MatSetValues(mat,1,&tmp,1,&i,&value,INSERT_VALUES);
26: }
27: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
28: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
29: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original matrix\n");
30: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
32: /* Form submatrix with rows 2-4 and columns 4-8 */
33: ISCreateStride(PETSC_COMM_SELF,3,2,1,&irkeep);
34: ISCreateStride(PETSC_COMM_SELF,5,4,1,&ickeep);
35: MatGetSubMatrices(mat,1,&irkeep,&ickeep,MAT_INITIAL_MATRIX,&submatrices);
36: submat = *submatrices;
37: PetscFree(submatrices);
38: /*
39: sviewer will cause the submatrices (one per processor) to be printed in the correct order
40: */
41: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Submatrices\n");
42: PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
43: MatView(submat,sviewer);
44: PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
45: PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
47: /* Zero the original matrix */
48: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original zeroed matrix\n");
49: MatZeroEntries(mat);
50: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
52: ISDestroy(irkeep);
53: ISDestroy(ickeep);
54: MatDestroy(submat);
55: MatDestroy(mat);
56: PetscFinalize();
57: return 0;
58: }
59: