Actual source code: ex2.c
1: /*$Id: ex2.c,v 1.53 2001/08/07 21:30:37 bsmith Exp $*/
3: static char help[] = "Tests PC and KSP on a tridiagonal matrix. Note that most\n\
4: users should employ the SLES interface instead of using PC directly.\n\n";
6: #include petscksp.h
7: #include petscpc.h
8: #include petsc.h
9: #include <stdio.h>
13: int main(int argc,char **args)
14: {
15: Mat mat; /* matrix */
16: Vec b,ustar,u; /* vectors (RHS, exact solution, approx solution) */
17: PC pc; /* PC context */
18: KSP ksp; /* KSP context */
19: int ierr,n = 10,i,its,col[3];
20: PetscScalar value[3],mone = -1.0,one = 1.0,zero = 0.0;
21: PCType pcname;
22: KSPType kspname;
23: PetscReal norm;
25: PetscInitialize(&argc,&args,(char *)0,help);
27: /* Create and initialize vectors */
28: VecCreateSeq(PETSC_COMM_SELF,n,&b);
29: VecCreateSeq(PETSC_COMM_SELF,n,&ustar);
30: VecCreateSeq(PETSC_COMM_SELF,n,&u);
31: VecSet(&one,ustar);
32: VecSet(&zero,u);
34: /* Create and assemble matrix */
35: MatCreateSeqAIJ(PETSC_COMM_SELF,n,n,3,PETSC_NULL,&mat);
36: value[0] = -1.0; value[1] = 2.0; value[2] = -1.0;
37: for (i=1; i<n-1; i++) {
38: col[0] = i-1; col[1] = i; col[2] = i+1;
39: MatSetValues(mat,1,&i,3,col,value,INSERT_VALUES);
40: }
41: i = n - 1; col[0] = n - 2; col[1] = n - 1;
42: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
43: i = 0; col[0] = 0; col[1] = 1; value[0] = 2.0; value[1] = -1.0;
44: MatSetValues(mat,1,&i,2,col,value,INSERT_VALUES);
45: MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
46: MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
48: /* Compute right-hand-side vector */
49: MatMult(mat,ustar,b);
51: /* Create PC context and set up data structures */
52: PCCreate(PETSC_COMM_WORLD,&pc);
53: PCSetType(pc,PCNONE);
54: PCSetFromOptions(pc);
55: PCSetOperators(pc,mat,mat,DIFFERENT_NONZERO_PATTERN);
56: PCSetVector(pc,u);
57: PCSetUp(pc);
59: /* Create KSP context and set up data structures */
60: KSPCreate(PETSC_COMM_WORLD,&ksp);
61: KSPSetType(ksp,KSPRICHARDSON);
62: KSPSetFromOptions(ksp);
63: KSPSetSolution(ksp,u);
64: KSPSetRhs(ksp,b);
65: PCSetOperators(pc,mat,mat,DIFFERENT_NONZERO_PATTERN);
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,&its);
74: VecAXPY(&mone,ustar,u);
75: VecNorm(u,NORM_2,&norm);
76: PetscPrintf(PETSC_COMM_SELF,"2 norm of error %A Number of iterations %d\n",norm,its);
78: /* Free data structures */
79: KSPDestroy(ksp);
80: VecDestroy(u);
81: VecDestroy(ustar);
82: VecDestroy(b);
83: MatDestroy(mat);
84: PCDestroy(pc);
86: PetscFinalize();
87: return 0;
88: }
89: