Actual source code: ex1.c
petsc-3.12.2 2019-11-22
2: static char help[] = "Tests solving linear system on 0 by 0 matrix.\n\n";
4: #include <petscksp.h>
6: int main(int argc,char **args)
7: {
8: Mat C;
10: PetscInt N = 0;
11: Vec u,b,x;
12: KSP ksp;
13: PetscReal norm;
15: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
17: /* create stiffness matrix */
18: MatCreate(PETSC_COMM_WORLD,&C);
19: MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,N,N);
20: MatSetFromOptions(C);
21: MatSetUp(C);
22: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
23: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
25: /* create right hand side and solution */
26: VecCreate(PETSC_COMM_WORLD,&u);
27: VecSetSizes(u,PETSC_DECIDE,N);
28: VecSetFromOptions(u);
29: VecDuplicate(u,&b);
30: VecDuplicate(u,&x);
31: VecSet(u,0.0);
32: VecSet(b,0.0);
34: VecAssemblyBegin(b);
35: VecAssemblyEnd(b);
37: /* solve linear system */
38: KSPCreate(PETSC_COMM_WORLD,&ksp);
39: KSPSetOperators(ksp,C,C);
40: KSPSetFromOptions(ksp);
41: KSPSolve(ksp,b,u);
43: MatMult(C,u,x);
44: VecAXPY(x,-1.0,b);
45: VecNorm(x,NORM_2,&norm);
47: KSPDestroy(&ksp);
48: VecDestroy(&u);
49: VecDestroy(&x);
50: VecDestroy(&b);
51: MatDestroy(&C);
52: PetscFinalize();
53: return ierr;
54: }
56: /*TEST
58: test:
59: args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
61: test:
62: suffix: 2
63: nsize: 2
64: args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
66: test:
67: suffix: 3
68: args: -pc_type sor -pc_sor_symmetric -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
70: test:
71: suffix: 5
72: args: -pc_type eisenstat -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
74: TEST*/