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