1 static char help[] = "Tests solving linear system on 0 by 0 matrix, and KSPLSQR convergence test handling.\n\n";
2
3 #include <petscksp.h>
4
GetConvergenceTestName(KSPConvergenceTestFn * converged,char name[],size_t n)5 static PetscErrorCode GetConvergenceTestName(KSPConvergenceTestFn *converged, char name[], size_t n)
6 {
7 PetscFunctionBegin;
8 if (converged == KSPConvergedDefault) {
9 PetscCall(PetscStrncpy(name, "default", n));
10 } else if (converged == KSPConvergedSkip) {
11 PetscCall(PetscStrncpy(name, "skip", n));
12 } else if (converged == KSPLSQRConvergedDefault) {
13 PetscCall(PetscStrncpy(name, "lsqr", n));
14 } else {
15 PetscCall(PetscStrncpy(name, "other", n));
16 }
17 PetscFunctionReturn(PETSC_SUCCESS);
18 }
19
main(int argc,char ** args)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;
28
29 PetscFunctionBeginUser;
30 PetscCall(PetscInitialize(&argc, &args, NULL, help));
31
32 /* create stiffness matrix */
33 PetscCall(MatCreate(PETSC_COMM_WORLD, &C));
34 PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, N, N));
35 PetscCall(MatSetFromOptions(C));
36 PetscCall(MatSetUp(C));
37 PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
38 PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));
39
40 /* create right-hand side and solution */
41 PetscCall(VecCreate(PETSC_COMM_WORLD, &u));
42 PetscCall(VecSetSizes(u, PETSC_DECIDE, N));
43 PetscCall(VecSetFromOptions(u));
44 PetscCall(VecDuplicate(u, &b));
45 PetscCall(VecDuplicate(u, &x));
46 PetscCall(VecSet(u, 0.0));
47 PetscCall(VecSet(b, 0.0));
48
49 PetscCall(VecAssemblyBegin(b));
50 PetscCall(VecAssemblyEnd(b));
51
52 /* solve linear system */
53 PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
54 PetscCall(KSPSetOperators(ksp, C, C));
55 PetscCall(KSPSetFromOptions(ksp));
56 PetscCall(KSPSolve(ksp, b, u));
57
58 /* test proper handling of convergence test by KSPLSQR */
59 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_lsqr", &flg, NULL));
60 if (flg) {
61 char *type;
62 char convtestname[16];
63 PetscBool islsqr;
64 KSPConvergenceTestFn *converged, *converged1;
65 PetscCtxDestroyFn *destroy, *destroy1;
66 void *ctx, *ctx1;
67
68 {
69 const char *typeP;
70 PetscCall(KSPGetType(ksp, &typeP));
71 PetscCall(PetscStrallocpy(typeP, &type));
72 }
73 PetscCall(PetscStrcmp(type, KSPLSQR, &islsqr));
74 PetscCall(KSPGetConvergenceTest(ksp, &converged, &ctx, &destroy));
75 PetscCall(GetConvergenceTestName(converged, convtestname, 16));
76 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "convergence test: %s\n", convtestname));
77 PetscCall(KSPSetType(ksp, KSPLSQR));
78 PetscCall(KSPGetConvergenceTest(ksp, &converged1, &ctx1, &destroy1));
79 PetscCheck(converged1 == KSPLSQRConvergedDefault, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test should be KSPLSQRConvergedDefault");
80 PetscCheck(destroy1 == KSPConvergedDefaultDestroy, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test destroy function should be KSPConvergedDefaultDestroy");
81 if (islsqr) {
82 PetscCheck(converged1 == converged, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test should be kept");
83 PetscCheck(destroy1 == destroy, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test destroy function should be kept");
84 PetscCheck(ctx1 == ctx, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test context should be kept");
85 }
86 PetscCall(GetConvergenceTestName(converged1, convtestname, 16));
87 PetscCall(KSPViewFromOptions(ksp, NULL, "-ksp1_view"));
88 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "convergence test: %s\n", convtestname));
89 PetscCall(KSPSetType(ksp, type));
90 PetscCall(KSPGetConvergenceTest(ksp, &converged1, &ctx1, &destroy1));
91 PetscCheck(converged1 == converged, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test not reverted properly");
92 PetscCheck(destroy1 == destroy, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test destroy function not reverted properly");
93 PetscCheck(ctx1 == ctx, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "convergence test context not reverted properly");
94 PetscCall(GetConvergenceTestName(converged1, convtestname, 16));
95 PetscCall(KSPViewFromOptions(ksp, NULL, "-ksp2_view"));
96 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "convergence test: %s\n", convtestname));
97 PetscCall(PetscFree(type));
98 }
99
100 PetscCall(MatMult(C, u, x));
101 PetscCall(VecAXPY(x, -1.0, b));
102 PetscCall(VecNorm(x, NORM_2, &norm));
103
104 PetscCall(KSPDestroy(&ksp));
105 PetscCall(VecDestroy(&u));
106 PetscCall(VecDestroy(&x));
107 PetscCall(VecDestroy(&b));
108 PetscCall(MatDestroy(&C));
109 PetscCall(PetscFinalize());
110 return 0;
111 }
112
113 /*TEST
114
115 test:
116 args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
117
118 test:
119 suffix: 2
120 nsize: 2
121 args: -pc_type jacobi -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
122
123 test:
124 suffix: 3
125 args: -pc_type sor -pc_sor_symmetric -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
126
127 test:
128 suffix: 5
129 args: -pc_type eisenstat -ksp_monitor_short -ksp_gmres_cgs_refinement_type refine_always
130
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
143
144 TEST*/
145