1 2 static char help[] = "Test MatAXPY(), and illustrate how to reduce number of mallocs used during MatSetValues() calls \n\ 3 Matrix C is copied from ~petsc/src/ksp/ksp/tutorials/ex5.c\n\n"; 4 /* 5 Example: ./ex132 -mat_view ascii::ascii_info 6 */ 7 8 #include <petscmat.h> 9 10 int main(int argc,char **args) 11 { 12 Mat C,C1,C2; /* matrix */ 13 PetscScalar v; 14 PetscInt Ii,J,Istart,Iend; 15 PetscErrorCode ierr; 16 PetscInt i,j,m = 3,n = 2; 17 PetscMPIInt size,rank; 18 PetscBool mat_nonsymmetric = PETSC_FALSE; 19 MatInfo info; 20 21 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 22 ierr = PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);CHKERRQ(ierr); 23 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 24 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 25 n = 2*size; 26 27 /* Set flag if we are doing a nonsymmetric problem; the default is symmetric. */ 28 ierr = PetscOptionsGetBool(NULL,NULL,"-mat_nonsym",&mat_nonsymmetric,NULL);CHKERRQ(ierr); 29 30 ierr = MatCreate(PETSC_COMM_WORLD,&C);CHKERRQ(ierr); 31 ierr = MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr); 32 ierr = MatSetFromOptions(C);CHKERRQ(ierr); 33 ierr = MatSeqAIJSetPreallocation(C,5,NULL);CHKERRQ(ierr); 34 35 ierr = MatGetOwnershipRange(C,&Istart,&Iend);CHKERRQ(ierr); 36 for (Ii=Istart; Ii<Iend; Ii++) { 37 v = -1.0; i = Ii/n; j = Ii - i*n; 38 if (i>0) {J = Ii - n; ierr = MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);} 39 if (i<m-1) {J = Ii + n; ierr = MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);} 40 if (j>0) {J = Ii - 1; ierr = MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);} 41 if (j<n-1) {J = Ii + 1; ierr = MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);} 42 v = 4.0; ierr = MatSetValues(C,1,&Ii,1,&Ii,&v,ADD_VALUES);CHKERRQ(ierr); 43 } 44 45 /* Make the matrix nonsymmetric if desired */ 46 if (mat_nonsymmetric) { 47 for (Ii=Istart; Ii<Iend; Ii++) { 48 v = -1.5; i = Ii/n; 49 if (i>1) {J = Ii-n-1; ierr = MatSetValues(C,1,&Ii,1,&J,&v,ADD_VALUES);CHKERRQ(ierr);} 50 } 51 } else { 52 ierr = MatSetOption(C,MAT_SYMMETRIC,PETSC_TRUE);CHKERRQ(ierr); 53 ierr = MatSetOption(C,MAT_SYMMETRY_ETERNAL,PETSC_TRUE);CHKERRQ(ierr); 54 } 55 ierr = MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 56 ierr = MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 57 58 /* First, create C1 = 2.0*C1 + C, C1 has less non-zeros than C */ 59 ierr = PetscPrintf(PETSC_COMM_WORLD, "\ncreate C1 = 2.0*C1 + C, C1 has less non-zeros than C \n");CHKERRQ(ierr); 60 ierr = MatCreate(PETSC_COMM_WORLD,&C1);CHKERRQ(ierr); 61 ierr = MatSetSizes(C1,PETSC_DECIDE,PETSC_DECIDE,m*n,m*n);CHKERRQ(ierr); 62 ierr = MatSetFromOptions(C1);CHKERRQ(ierr); 63 ierr = MatSeqAIJSetPreallocation(C1,1,NULL);CHKERRQ(ierr); 64 for (Ii=Istart; Ii<Iend; Ii++) { 65 ierr = MatSetValues(C1,1,&Ii,1,&Ii,&v,ADD_VALUES);CHKERRQ(ierr); 66 } 67 ierr = MatAssemblyBegin(C1,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 68 ierr = MatAssemblyEnd(C1,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 69 ierr = PetscPrintf(PETSC_COMM_WORLD," MatAXPY(C1,2.0,C,DIFFERENT_NONZERO_PATTERN)...\n");CHKERRQ(ierr); 70 ierr = MatAXPY(C1,2.0,C,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 71 ierr = MatGetInfo(C1,MAT_GLOBAL_SUM,&info);CHKERRQ(ierr); 72 ierr = PetscPrintf(PETSC_COMM_WORLD," C1: nz_allocated = %g; nz_used = %g; nz_unneeded = %g\n",info.nz_allocated,info.nz_used, info.nz_unneeded);CHKERRQ(ierr); 73 74 /* Secondly, create C2 = 2.0*C2 + C, C2 has non-zero pattern of C2 + C */ 75 ierr = PetscPrintf(PETSC_COMM_WORLD, "\ncreate C2 = 2.0*C2 + C, C2 has non-zero pattern of C2 + C \n");CHKERRQ(ierr); 76 ierr = MatDuplicate(C,MAT_DO_NOT_COPY_VALUES,&C2);CHKERRQ(ierr); 77 78 for (Ii=Istart; Ii<Iend; Ii++) { 79 v = 1.0; 80 ierr = MatSetValues(C2,1,&Ii,1,&Ii,&v,ADD_VALUES);CHKERRQ(ierr); 81 } 82 ierr = MatAssemblyBegin(C2,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 83 ierr = MatAssemblyEnd(C2,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 84 ierr = PetscPrintf(PETSC_COMM_WORLD," MatAXPY(C2,2.0,C,SUBSET_NONZERO_PATTERN)...\n");CHKERRQ(ierr); 85 ierr = MatAXPY(C2,2.0,C,SUBSET_NONZERO_PATTERN);CHKERRQ(ierr); 86 ierr = MatGetInfo(C2,MAT_GLOBAL_SUM,&info);CHKERRQ(ierr); 87 ierr = PetscPrintf(PETSC_COMM_WORLD," C2: nz_allocated = %g; nz_used = %g; nz_unneeded = %g\n",info.nz_allocated,info.nz_used, info.nz_unneeded);CHKERRQ(ierr); 88 89 ierr = MatDestroy(&C1);CHKERRQ(ierr); 90 ierr = MatDestroy(&C2);CHKERRQ(ierr); 91 ierr = MatDestroy(&C);CHKERRQ(ierr); 92 93 ierr = PetscFinalize(); 94 return ierr; 95 } 96 97 /*TEST 98 99 test: 100 101 TEST*/ 102