Actual source code: ex44.c


  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;
 12:   PetscScalar    v[3] = {-1./2., 1., -1./2.};
 13:   PetscInt       j[3];
 14:   PetscInt       k=15;
 15:   PetscInt       M,m=420;

 17:   PetscInitialize(&argc,&args,(char*)0,help);
 18:   PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
 19:   PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL);

 21:   KSPCreate(PETSC_COMM_WORLD,&ksp);
 22:   KSPSetFromOptions(ksp);
 23:   KSPGetOperators(ksp,&A,NULL);

 25:   MatSetSizes(A,m,m,PETSC_DETERMINE,PETSC_DETERMINE);
 26:   MatSetFromOptions(A);
 27:   MatSetUp(A);
 28:   MatGetOwnershipRange(A,&Istart,&Iend);
 29:   MatGetSize(A,&M,NULL);
 30:   for (Ii=Istart; Ii<Iend; Ii++) {
 31:     j[0] = Ii - k;
 32:     j[1] = Ii;
 33:     j[2] = (Ii + k) < M ? (Ii + k) : -1;
 34:     MatSetValues(A,1,&Ii,3,j,v,INSERT_VALUES);
 35:   }
 36:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 37:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 38:   MatCreateVecs(A,&x,&b);

 40:   VecSetFromOptions(b);
 41:   VecSet(b,1.0);
 42:   VecSetFromOptions(x);
 43:   VecSet(x,2.0);

 45:   KSPSolve(ksp,b,x);

 47:   VecDestroy(&b);
 48:   VecDestroy(&x);
 49:   KSPDestroy(&ksp);

 51:   PetscFinalize();
 52:   return 0;
 53: }