Actual source code: ex3.c
petsc-3.6.1 2015-08-06
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>
11: int main(int argc,char **args)
12: {
13: Mat mat; /* matrix */
14: Vec b,ustar,u; /* vectors (RHS, exact solution, approx solution) */
15: PC pc; /* PC context */
16: KSP ksp; /* KSP context */
18: PetscInt n = 10,i,its,col[3];
19: PetscScalar value[3];
20: KSPType kspname;
21: PCType pcname;
23: PetscInitialize(&argc,&args,(char*)0,help);
24: PetscOptionsGetInt(NULL,"-n",&n,NULL);
26: /* Create and initialize vectors */
27: VecCreateSeq(PETSC_COMM_SELF,n,&b);
28: VecCreateSeq(PETSC_COMM_SELF,n,&ustar);
29: VecCreateSeq(PETSC_COMM_SELF,n,&u);
30: VecSet(ustar,1.0);
31: VecSet(u,0.0);
33: /* Create and assemble matrix */
34: MatCreateSeqAIJ(PETSC_COMM_SELF,n,n,3,NULL,&mat);
35: value[0] = -1.0; value[1] = 2.0; value[2] = -1.0;
36: for (i=1; i<n-1; i++) {
37: col[0] = i-1; col[1] = i; col[2] = i+1;
38: MatSetValues(mat,1,&i,3,col,value,INSERT_VALUES);
39: }
40: i = n - 1; col[0] = n - 2; col[1] = n - 1;
41: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
42: i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
43: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
44: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
45: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
47: /* Compute right-hand-side vector */
48: MatMult(mat,ustar,b);
50: /* Create PC context and set up data structures */
51: PCCreate(PETSC_COMM_WORLD,&pc);
52: PCSetType(pc,PCNONE);
53: PCSetFromOptions(pc);
54: PCSetOperators(pc,mat,mat);
55: PCSetUp(pc);
57: /* Create KSP context and set up data structures */
58: KSPCreate(PETSC_COMM_WORLD,&ksp);
59: KSPSetType(ksp,KSPRICHARDSON);
60: KSPSetFromOptions(ksp);
61: PCSetOperators(pc,mat,mat);
62: KSPSetPC(ksp,pc);
63: KSPSetUp(ksp);
65: /* Solve the problem */
66: KSPGetType(ksp,&kspname);
67: PCGetType(pc,&pcname);
68: PetscPrintf(PETSC_COMM_SELF,"Running %s with %s preconditioning\n",kspname,pcname);
69: KSPSolve(ksp,b,u);
70: KSPGetIterationNumber(ksp,&its);
71: PetscPrintf(PETSC_COMM_WORLD,"Number of iterations %D\n",its);
73: /* Free data structures */
74: KSPDestroy(&ksp);
75: VecDestroy(&u);
76: VecDestroy(&ustar);
77: VecDestroy(&b);
78: MatDestroy(&mat);
79: PCDestroy(&pc);
80: PetscFinalize();
81: return 0;
82: }