Actual source code: ex3.c
2: static char help[] = "Demonstrates the use of fast Richardson for SOR. And\n\
3: also tests the MatSOR() routines. Input parameters are:\n\
4: -n <n> : problem dimension\n\n";
6: #include <petscksp.h>
7: #include <petscpc.h>
9: int main(int argc,char **args)
10: {
11: Mat mat; /* matrix */
12: Vec b,ustar,u; /* vectors (RHS, exact solution, approx solution) */
13: PC pc; /* PC context */
14: KSP ksp; /* KSP context */
15: PetscInt n = 10,i,its,col[3];
16: PetscScalar value[3];
17: KSPType kspname;
18: PCType pcname;
20: PetscInitialize(&argc,&args,(char*)0,help);
21: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
23: /* Create and initialize vectors */
24: VecCreateSeq(PETSC_COMM_SELF,n,&b);
25: VecCreateSeq(PETSC_COMM_SELF,n,&ustar);
26: VecCreateSeq(PETSC_COMM_SELF,n,&u);
27: VecSet(ustar,1.0);
28: VecSet(u,0.0);
30: /* Create and assemble matrix */
31: MatCreate(PETSC_COMM_SELF,&mat);
32: MatSetType(mat,MATSEQAIJ);
33: MatSetSizes(mat,n,n,n,n);
34: MatSetFromOptions(mat);
35: MatSeqAIJSetPreallocation(mat,3,NULL);
36: MatSeqBAIJSetPreallocation(mat,1,3,NULL);
37: MatSeqSBAIJSetPreallocation(mat,1,3,NULL);
38: MatSeqSELLSetPreallocation(mat,3,NULL);
39: value[0] = -1.0; value[1] = 2.0; value[2] = -1.0;
40: for (i=1; i<n-1; i++) {
41: col[0] = i-1; col[1] = i; col[2] = i+1;
42: MatSetValues(mat,1,&i,3,col,value,INSERT_VALUES);
43: }
44: i = n - 1; col[0] = n - 2; col[1] = n - 1;
45: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
46: i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
47: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
48: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
49: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
51: /* Compute right-hand-side vector */
52: MatMult(mat,ustar,b);
54: /* Create PC context and set up data structures */
55: PCCreate(PETSC_COMM_WORLD,&pc);
56: PCSetType(pc,PCNONE);
57: PCSetFromOptions(pc);
58: PCSetOperators(pc,mat,mat);
59: PCSetUp(pc);
61: /* Create KSP context and set up data structures */
62: KSPCreate(PETSC_COMM_WORLD,&ksp);
63: KSPSetType(ksp,KSPRICHARDSON);
64: KSPSetFromOptions(ksp);
65: PCSetOperators(pc,mat,mat);
66: KSPSetPC(ksp,pc);
67: KSPSetUp(ksp);
69: /* Solve the problem */
70: KSPGetType(ksp,&kspname);
71: PCGetType(pc,&pcname);
72: PetscPrintf(PETSC_COMM_SELF,"Running %s with %s preconditioning\n",kspname,pcname);
73: KSPSolve(ksp,b,u);
74: KSPGetIterationNumber(ksp,&its);
75: PetscPrintf(PETSC_COMM_WORLD,"Number of iterations %D\n",its);
77: /* Free data structures */
78: KSPDestroy(&ksp);
79: VecDestroy(&u);
80: VecDestroy(&ustar);
81: VecDestroy(&b);
82: MatDestroy(&mat);
83: PCDestroy(&pc);
84: PetscFinalize();
85: return 0;
86: }
88: /*TEST
90: testset:
91: args: -ksp_type gmres -ksp_monitor_short -pc_type sor -pc_sor_symmetric
92: output_file: output/ex3_1.out
93: test:
94: suffix: sor_aij
95: test:
96: suffix: sor_seqbaij
97: args: -mat_type seqbaij
98: test:
99: suffix: sor_seqsbaij
100: args: -mat_type seqbaij
101: test:
102: suffix: sor_seqsell
103: args: -mat_type seqsell
105: TEST*/