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().nn";
5: #include petscmat.h
7: int main(int argc,char **argv)
8: {
9: Mat mat,submat,*submatrices;
10: int m = 10,n = 10,i = 4,tmp,ierr;
11: IS irkeep,ickeep;
12: PetscScalar value = 1.0;
13: PetscViewer sviewer;
15: PetscInitialize(&argc,&argv,(char *)0,help);
16: PetscViewerSetFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON);
17: PetscViewerSetFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_COMMON);
19: MatCreate(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,&mat);
20: MatSetFromOptions(mat);
21: for (i=0; i<m; i++) {
22: value = (PetscReal)i+1; tmp = i % 5;
23: MatSetValues(mat,1,&tmp,1,&i,&value,INSERT_VALUES);
24: }
25: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
26: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
27: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original matrixn");
28: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
30: /* Form submatrix with rows 2-4 and columns 4-8 */
31: ISCreateStride(PETSC_COMM_SELF,3,2,1,&irkeep);
32: ISCreateStride(PETSC_COMM_SELF,5,4,1,&ickeep);
33: MatGetSubMatrices(mat,1,&irkeep,&ickeep,MAT_INITIAL_MATRIX,&submatrices);
34: submat = *submatrices;
35: PetscFree(submatrices);
36: /*
37: sviewer will cause the submatrices (one per processor) to be printed in the correct order
38: */
39: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Submatricesn");
40: PetscViewerGetSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
41: MatView(submat,sviewer);
42: PetscViewerRestoreSingleton(PETSC_VIEWER_STDOUT_WORLD,&sviewer);
43: PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
45: /* Zero the original matrix */
46: PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original zeroed matrixn");
47: MatZeroEntries(mat);
48: MatView(mat,PETSC_VIEWER_STDOUT_WORLD);
50: ISDestroy(irkeep);
51: ISDestroy(ickeep);
52: MatDestroy(submat);
53: MatDestroy(mat);
54: PetscFinalize();
55: return 0;
56: }
57: