1 2 static char help[] = "Demonstrates the use of fast Richardson for SOR. And tests\n\ 3 the MatSOR() routines.\n\n"; 4 5 #include <petscpc.h> 6 7 int main(int argc,char **args) 8 { 9 Mat mat; 10 Vec b,u; 11 PC pc; 12 PetscErrorCode ierr; 13 PetscInt n = 5,i,col[3]; 14 PetscScalar value[3]; 15 16 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 17 /* Create vectors */ 18 ierr = VecCreateSeq(PETSC_COMM_SELF,n,&b);CHKERRQ(ierr); 19 ierr = VecCreateSeq(PETSC_COMM_SELF,n,&u);CHKERRQ(ierr); 20 21 /* Create and assemble matrix */ 22 ierr = MatCreateSeqDense(PETSC_COMM_SELF,n,n,NULL,&mat);CHKERRQ(ierr); 23 value[0] = -1.0; value[1] = 2.0; value[2] = -1.0; 24 for (i=1; i<n-1; i++) { 25 col[0] = i-1; col[1] = i; col[2] = i+1; 26 ierr = MatSetValues(mat,1,&i,3,col,value,INSERT_VALUES);CHKERRQ(ierr); 27 } 28 i = n - 1; col[0] = n - 2; col[1] = n - 1; 29 ierr = MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 30 i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0; 31 ierr = MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 32 ierr = MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 33 ierr = MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 34 35 /* Create PC context and set up data structures */ 36 ierr = PCCreate(PETSC_COMM_WORLD,&pc);CHKERRQ(ierr); 37 ierr = PCSetType(pc,PCSOR);CHKERRQ(ierr); 38 ierr = PCSetFromOptions(pc);CHKERRQ(ierr); 39 ierr = PCSetOperators(pc,mat,mat);CHKERRQ(ierr); 40 ierr = PCSetUp(pc);CHKERRQ(ierr); 41 42 value[0] = 1.0; 43 for (i=0; i<n; i++) { 44 ierr = VecSet(u,0.0);CHKERRQ(ierr); 45 ierr = VecSetValues(u,1,&i,value,INSERT_VALUES);CHKERRQ(ierr); 46 ierr = VecAssemblyBegin(u);CHKERRQ(ierr); 47 ierr = VecAssemblyEnd(u);CHKERRQ(ierr); 48 ierr = PCApply(pc,u,b);CHKERRQ(ierr); 49 ierr = VecView(b,PETSC_VIEWER_STDOUT_SELF);CHKERRQ(ierr); 50 } 51 52 /* Free data structures */ 53 ierr = MatDestroy(&mat);CHKERRQ(ierr); 54 ierr = PCDestroy(&pc);CHKERRQ(ierr); 55 ierr = VecDestroy(&u);CHKERRQ(ierr); 56 ierr = VecDestroy(&b);CHKERRQ(ierr); 57 ierr = PetscFinalize(); 58 return ierr; 59 } 60 61 /*TEST 62 63 test: 64 65 TEST*/ 66 67