Actual source code: ex1.c
2: static char help[] = "Tests solving linear system on 0 by 0 matrix, and KSPLSQR convergence test handling.\n\n";
4: #include <petscksp.h>
6: static PetscErrorCode GetConvergenceTestName(PetscErrorCode (*converged)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),char name[],size_t n)
7: {
8: if (converged == KSPConvergedDefault) {
9: PetscStrncpy(name,"default",n);
10: } else if (converged == KSPConvergedSkip) {
11: PetscStrncpy(name,"skip",n);
12: } else if (converged == KSPLSQRConvergedDefault) {
13: PetscStrncpy(name,"lsqr",n);
14: } else {
15: PetscStrncpy(name,"other",n);
16: }
17: return 0;
18: }
20: int main(int argc,char **args)
21: {
22: Mat C;
23: PetscInt N = 0;
24: Vec u,b,x;
25: KSP ksp;
26: PetscReal norm;
27: PetscBool flg=PETSC_FALSE;
29: PetscInitialize(&argc,&args,(char*)0,help);
31: /* create stiffness matrix */
32: MatCreate(PETSC_COMM_WORLD,&C);
33: MatSetSizes(C,PETSC_DECIDE,PETSC_DECIDE,N,N);
34: MatSetFromOptions(C);
35: MatSetUp(C);
36: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
37: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
39: /* create right hand side and solution */
40: VecCreate(PETSC_COMM_WORLD,&u);
41: VecSetSizes(u,PETSC_DECIDE,N);
42: VecSetFromOptions(u);
43: VecDuplicate(u,&b);
44: VecDuplicate(u,&x);
45: VecSet(u,0.0);
46: VecSet(b,0.0);
48: VecAssemblyBegin(b);
49: VecAssemblyEnd(b);
51: /* solve linear system */
52: KSPCreate(PETSC_COMM_WORLD,&ksp);
53: KSPSetOperators(ksp,C,C);
54: KSPSetFromOptions(ksp);
55: KSPSolve(ksp,b,u);
57: /* test proper handling of convergence test by KSPLSQR */
58: PetscOptionsGetBool(NULL,NULL,"-test_lsqr",&flg,NULL);
59: if (flg) {
60: char *type;
61: char convtestname[16];
62: PetscBool islsqr;
63: PetscErrorCode (*converged)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*);
64: PetscErrorCode (*converged1)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*);
65: PetscErrorCode (*destroy)(void*),(*destroy1)(void*);
66: void *ctx,*ctx1;
68: {
69: const char *typeP;
70: KSPGetType(ksp,&typeP);
71: PetscStrallocpy(typeP,&type);
72: }
73: PetscStrcmp(type,KSPLSQR,&islsqr);
74: KSPGetConvergenceTest(ksp,&converged,&ctx,&destroy);
75: GetConvergenceTestName(converged,convtestname,16);
76: PetscPrintf(PETSC_COMM_WORLD,"convergence test: %s\n",convtestname);
77: KSPSetType(ksp,KSPLSQR);
78: KSPGetConvergenceTest(ksp,&converged1,&ctx1,&destroy1);
81: if (islsqr) {
85: }
86: GetConvergenceTestName(converged1,convtestname,16);
87: KSPViewFromOptions(ksp,NULL,"-ksp1_view");
88: PetscPrintf(PETSC_COMM_WORLD,"convergence test: %s\n",convtestname);
89: KSPSetType(ksp,type);
90: KSPGetConvergenceTest(ksp,&converged1,&ctx1,&destroy1);
94: GetConvergenceTestName(converged1,convtestname,16);
95: KSPViewFromOptions(ksp,NULL,"-ksp2_view");
96: PetscPrintf(PETSC_COMM_WORLD,"convergence test: %s\n",convtestname);
97: PetscFree(type);
98: }
100: MatMult(C,u,x);
101: VecAXPY(x,-1.0,b);
102: VecNorm(x,NORM_2,&norm);
104: KSPDestroy(&ksp);
105: VecDestroy(&u);
106: VecDestroy(&x);
107: VecDestroy(&b);
108: MatDestroy(&C);
109: PetscFinalize();
110: return 0;
111: }
113: /*TEST
115: test:
116: args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
118: test:
119: suffix: 2
120: nsize: 2
121: args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
123: test:
124: suffix: 3
125: args: -pc_type sor -pc_sor_symmetric -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
127: test:
128: suffix: 5
129: args: -pc_type eisenstat -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
131: testset:
132: args: -test_lsqr -ksp{,1,2}_view -pc_type jacobi
133: filter: grep -E "(^ type:|preconditioning|norm type|convergence test:)"
134: test:
135: suffix: lsqr_0
136: args: -ksp_convergence_test {{default skip}separate output}
137: test:
138: suffix: lsqr_1
139: args: -ksp_type cg -ksp_convergence_test {{default skip}separate output}
140: test:
141: suffix: lsqr_2
142: args: -ksp_type lsqr
144: TEST*/