Actual source code: ex12.c
petsc-3.7.3 2016-08-01
2: static char help[] = "Tests the use of MatZeroRows() for parallel matrices.\n\
3: This example also tests the use of MatDuplicate() for both MPIAIJ and MPIBAIJ matrices";
5: #include <petscmat.h>
7: extern PetscErrorCode TestMatZeroRows_Basic(Mat,IS,PetscScalar);
8: extern PetscErrorCode TestMatZeroRows_with_no_allocation(Mat,IS,PetscScalar);
12: int main(int argc,char **args)
13: {
14: Mat A;
15: PetscInt i,j,m = 3,n,Ii,J,Imax;
16: PetscMPIInt rank,size;
18: PetscScalar v,diag=-4.0;
19: IS is;
21: PetscInitialize(&argc,&args,(char*)0,help);
22: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
23: MPI_Comm_size(PETSC_COMM_WORLD,&size);
24: n = 2*size;
26: /* create A Square matrix for the five point stencil,YET AGAIN*/
27: MatCreate(PETSC_COMM_WORLD,&A);
28: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);
29: MatSetFromOptions(A);
30: MatSetUp(A);
31: for (i=0; i<m; i++) {
32: for (j=2*rank; j<2*rank+2; j++) {
33: v = -1.0; Ii = j + n*i;
34: if (i>0) {J = Ii - n; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
35: if (i<m-1) {J = Ii + n; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
36: if (j>0) {J = Ii - 1; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
37: if (j<n-1) {J = Ii + 1; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
38: v = 4.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
39: }
40: }
41: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
42: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
44: /* Create AN IS required by MatZeroRows() */
45: Imax = n*rank; if (Imax>= n*m -m - 1) Imax = m*n - m - 1;
46: ISCreateStride(PETSC_COMM_SELF,m,Imax,1,&is);
48: TestMatZeroRows_Basic(A,is,0.0);
49: TestMatZeroRows_Basic(A,is,diag);
51: TestMatZeroRows_with_no_allocation(A,is,0.0);
52: TestMatZeroRows_with_no_allocation(A,is,diag);
54: MatDestroy(&A);
56: /* Now Create a rectangular matrix with five point stencil (app)
57: n+size is used so that this dimension is always divisible by size.
58: This way, we can always use bs = size for any number of procs */
59: MatCreate(PETSC_COMM_WORLD,&A);
60: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m*n,m*(n+size));
61: MatSetFromOptions(A);
62: MatSetUp(A);
63: for (i=0; i<m; i++) {
64: for (j=2*rank; j<2*rank+2; j++) {
65: v = -1.0; Ii = j + n*i;
66: if (i>0) {J = Ii - n; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
67: if (i<m-1) {J = Ii + n; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
68: if (j>0) {J = Ii - 1; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
69: if (j<n+size-1) {J = Ii + 1; MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);}
70: v = 4.0; MatSetValues(A,1,&Ii,1,&Ii,&v,INSERT_VALUES);
71: }
72: }
73: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
74: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
76: TestMatZeroRows_Basic(A,is,0.0);
77: TestMatZeroRows_Basic(A,is,diag);
79: MatDestroy(&A);
80: ISDestroy(&is);
81: PetscFinalize();
82: return 0;
83: }
87: PetscErrorCode TestMatZeroRows_Basic(Mat A,IS is,PetscScalar diag)
88: {
89: Mat B;
91: PetscBool keepnonzeropattern;
93: /* Now copy A into B, and test it with MatZeroRows() */
94: MatDuplicate(A,MAT_COPY_VALUES,&B);
96: PetscOptionsHasName(NULL,NULL,"-keep_nonzero_pattern",&keepnonzeropattern);
97: if (keepnonzeropattern) {
98: MatSetOption(B,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE);
99: }
101: MatZeroRowsIS(B,is,diag,0,0);
102: MatView(B,PETSC_VIEWER_STDOUT_WORLD);
103: MatDestroy(&B);
104: return 0;
105: }
109: PetscErrorCode TestMatZeroRows_with_no_allocation(Mat A,IS is,PetscScalar diag)
110: {
111: Mat B;
114: /* Now copy A into B, and test it with MatZeroRows() */
115: MatDuplicate(A,MAT_COPY_VALUES,&B);
116: /* Set this flag after assembly. This way, it affects only MatZeroRows() */
117: MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_TRUE);
119: MatZeroRowsIS(B,is,diag,0,0);
120: MatView(B,PETSC_VIEWER_STDOUT_WORLD);
121: MatDestroy(&B);
122: return 0;
123: }