Actual source code: ex44.c
petsc-3.12.0 2019-09-29
2: static char help[] = "Solves a tridiagonal linear system. Designed to compare SOR for different Mat impls.\n\n";
4: #include <petscksp.h>
6: int main(int argc,char **args)
7: {
8: KSP ksp; /* linear solver context */
9: Mat A; /* linear system matrix */
10: Vec x,b; /* approx solution, RHS */
11: PetscInt Ii,Istart,Iend;
13: PetscScalar v[3] = {-1./2., 1., -1./2.};
14: PetscInt j[3];
15: PetscInt k=15;
16: PetscInt M,m=420;
18: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
19: PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
20: PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);
22: KSPCreate(PETSC_COMM_WORLD,&ksp);
23: KSPSetFromOptions(ksp);
24: KSPGetOperators(ksp,&A,NULL);
26: MatSetSizes(A,m,m,PETSC_DETERMINE,PETSC_DETERMINE);
27: MatSetFromOptions(A);
28: MatSetUp(A);
29: MatGetOwnershipRange(A,&Istart,&Iend);
30: MatGetSize(A,&M,NULL);
31: for (Ii=Istart; Ii<Iend; Ii++) {
32: j[0] = Ii - k;
33: j[1] = Ii;
34: j[2] = (Ii + k) < M ? (Ii + k) : -1;
35: MatSetValues(A,1,&Ii,3,j,v,INSERT_VALUES);
36: }
37: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
38: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
39: MatCreateVecs(A,&x,&b);
41: VecSetFromOptions(b);
42: VecSet(b,1.0);
43: VecSetFromOptions(x);
44: VecSet(x,2.0);
46: KSPSolve(ksp,b,x);
48: VecDestroy(&b);
49: VecDestroy(&x);
50: KSPDestroy(&ksp);
52: PetscFinalize();
53: return ierr;
54: }