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 7*d71ae5a4SJacob Faibussowitsch int main(int argc, char **args) 8*d71ae5a4SJacob Faibussowitsch { 9c4762a1bSJed Brown Mat mat; 10c4762a1bSJed Brown Vec b, u; 11c4762a1bSJed Brown PC pc; 12c4762a1bSJed Brown PetscInt n = 5, i, col[3]; 13c4762a1bSJed Brown PetscScalar value[3]; 14c4762a1bSJed Brown 15327415f7SBarry Smith PetscFunctionBeginUser; 169566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 17c4762a1bSJed Brown /* Create vectors */ 189566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &b)); 199566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &u)); 20c4762a1bSJed Brown 21c4762a1bSJed Brown /* Create and assemble matrix */ 229566063dSJacob Faibussowitsch PetscCall(MatCreateSeqDense(PETSC_COMM_SELF, n, n, NULL, &mat)); 239371c9d4SSatish Balay value[0] = -1.0; 249371c9d4SSatish Balay value[1] = 2.0; 259371c9d4SSatish Balay value[2] = -1.0; 26c4762a1bSJed Brown for (i = 1; i < n - 1; i++) { 279371c9d4SSatish Balay col[0] = i - 1; 289371c9d4SSatish Balay col[1] = i; 299371c9d4SSatish Balay col[2] = i + 1; 309566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 1, &i, 3, col, value, INSERT_VALUES)); 31c4762a1bSJed Brown } 329371c9d4SSatish Balay i = n - 1; 339371c9d4SSatish Balay col[0] = n - 2; 349371c9d4SSatish Balay col[1] = n - 1; 359566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 1, &i, 2, col, value, INSERT_VALUES)); 369371c9d4SSatish Balay i = 0; 379371c9d4SSatish Balay col[0] = 0; 389371c9d4SSatish Balay col[1] = 1; 399371c9d4SSatish Balay value[0] = 2.0; 409371c9d4SSatish Balay value[1] = -1.0; 419566063dSJacob Faibussowitsch PetscCall(MatSetValues(mat, 1, &i, 2, col, value, INSERT_VALUES)); 429566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(mat, MAT_FINAL_ASSEMBLY)); 439566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(mat, MAT_FINAL_ASSEMBLY)); 44c4762a1bSJed Brown 45c4762a1bSJed Brown /* Create PC context and set up data structures */ 469566063dSJacob Faibussowitsch PetscCall(PCCreate(PETSC_COMM_WORLD, &pc)); 479566063dSJacob Faibussowitsch PetscCall(PCSetType(pc, PCSOR)); 489566063dSJacob Faibussowitsch PetscCall(PCSetFromOptions(pc)); 499566063dSJacob Faibussowitsch PetscCall(PCSetOperators(pc, mat, mat)); 509566063dSJacob Faibussowitsch PetscCall(PCSetUp(pc)); 51c4762a1bSJed Brown 52c4762a1bSJed Brown value[0] = 1.0; 53c4762a1bSJed Brown for (i = 0; i < n; i++) { 549566063dSJacob Faibussowitsch PetscCall(VecSet(u, 0.0)); 559566063dSJacob Faibussowitsch PetscCall(VecSetValues(u, 1, &i, value, INSERT_VALUES)); 569566063dSJacob Faibussowitsch PetscCall(VecAssemblyBegin(u)); 579566063dSJacob Faibussowitsch PetscCall(VecAssemblyEnd(u)); 589566063dSJacob Faibussowitsch PetscCall(PCApply(pc, u, b)); 599566063dSJacob Faibussowitsch PetscCall(VecView(b, PETSC_VIEWER_STDOUT_SELF)); 60c4762a1bSJed Brown } 61c4762a1bSJed Brown 62c4762a1bSJed Brown /* Free data structures */ 639566063dSJacob Faibussowitsch PetscCall(MatDestroy(&mat)); 649566063dSJacob Faibussowitsch PetscCall(PCDestroy(&pc)); 659566063dSJacob Faibussowitsch PetscCall(VecDestroy(&u)); 669566063dSJacob Faibussowitsch PetscCall(VecDestroy(&b)); 679566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 68b122ec5aSJacob Faibussowitsch return 0; 69c4762a1bSJed Brown } 70c4762a1bSJed Brown 71c4762a1bSJed Brown /*TEST 72c4762a1bSJed Brown 73c4762a1bSJed Brown test: 74c4762a1bSJed Brown 75c4762a1bSJed Brown TEST*/ 76