| Elena, I have attached a modification to src/snes/tutorials/ex5.c that adds a monitor routine in the style I think you are suggesting. |
ex5.c
Description: Binary data
Below I cut and paste the beginning of the output from running the command ./ex5 -ksp_type cg -ksp_monitor_true_residual -ksp_norm_type unpreconditioned -pc_type jacobi -da_refine 3 0 KSP unpreconditioned resid norm 1.265943996096e+00 true resid norm 1.265943996096e+00 ||r(i)||/||b|| 1.000000000000e+00 My monitor 0 1.265943996096e+00 1 KSP unpreconditioned resid norm 1.030361071579e+00 true resid norm 1.030361071579e+00 ||r(i)||/||b|| 8.139073092933e-01 My monitor 1 1.030361071579e+00 2 KSP unpreconditioned resid norm 7.753237278390e-01 true resid norm 7.753237278390e-01 ||r(i)||/||b|| 6.124470989473e-01 My monitor 2 7.753237278390e-01 3 KSP unpreconditioned resid norm 6.674186105521e-01 true resid norm 6.674186105521e-01 ||r(i)||/||b|| 5.272102183115e-01 My monitor 3 6.674186105521e-01 4 KSP unpreconditioned resid norm 5.745948088398e-01 true resid norm 5.745948088398e-01 ||r(i)||/||b|| 4.538864362181e-01 My monitor 4 5.745948088398e-01 5 KSP unpreconditioned resid norm 5.103132053010e-01 true resid norm 5.103132053010e-01 ||r(i)||/||b|| 4.031088317292e-01 My monitor 5 5.103132053010e-01 6 KSP unpreconditioned resid norm 4.581737850155e-01 true resid norm 4.581737850155e-01 ||r(i)||/||b|| 3.619226335670e-01 My monitor 6 4.581737850155e-01 7 KSP unpreconditioned resid norm 4.202213342980e-01 true resid norm 4.202213342980e-01 ||r(i)||/||b|| 3.319430682509e-01 My monitor 7 4.202213342980e-01 8 KSP unpreconditioned resid norm 3.936600255123e-01 true resid norm 3.936600255123e-01 ||r(i)||/||b|| 3.109616434267e-01 My monitor 8 3.936600255123e-01 9 KSP unpreconditioned resid norm 3.811944420804e-01 true resid norm 3.811944420804e-01 ||r(i)||/||b|| 3.011147754212e-01 My monitor 9 3.811944420804e-01 10 KSP unpreconditioned resid norm 3.851182669108e-01 true resid norm 3.851182669108e-01 ||r(i)||/||b|| 3.042143002363e-01 My monitor 10 3.851182669108e-01 11 KSP unpreconditioned resid norm 4.107620195902e-01 true resid norm 4.107620195902e-01 ||r(i)||/||b|| 3.244709251411e-01 My monitor 11 4.107620195902e-01 12 KSP unpreconditioned resid norm 3.678610761984e-01 true resid norm 3.678610761984e-01 ||r(i)||/||b|| 2.905824249198e-01 My monitor 12 3.678610761984e-01 13 KSP unpreconditioned resid norm 3.891700469761e-01 true resid norm 3.891700469761e-01 ||r(i)||/||b|| 3.074149000083e-01 My monitor 13 3.891700469761e-01 14 KSP unpreconditioned resid norm 4.123002052123e-01 true resid norm 4.123002052123e-01 ||r(i)||/||b|| 3.256859754331e-01 My monitor 14 4.123002052123e-01 15 KSP unpreconditioned resid norm 4.456104079353e-01 true resid norm 4.456104079353e-01 ||r(i)||/||b|| 3.519985159765e-01 My monitor 15 4.456104079353e-01 16 KSP unpreconditioned resid norm 5.125721163597e-01 true resid norm 5.125721163597e-01 ||r(i)||/||b|| 4.048932005999e-01 My monitor 16 5.125721163597e-01 17 KSP unpreconditioned resid norm 4.475156370525e-01 true resid norm 4.475156370525e-01 ||r(i)||/||b|| 3.535035028662e-01 My monitor 17 4.475156370525e-01 18 KSP unpreconditioned resid norm 2.977755423590e-01 true resid norm 2.977755423590e-01 ||r(i)||/||b|| 2.352201545070e-01 My monitor 18 2.977755423590e-01 19 KSP unpreconditioned resid norm 2.317275576684e-01 true resid norm 2.317275576684e-01 ||r(i)||/||b|| 1.830472425186e-01 My monitor 19 2.317275576684e-01 20 KSP unpreconditioned resid norm 2.388542347249e-01 true resid norm 2.388542347249e-01 ||r(i)||/||b|| 1.886767783262e-01 My monitor 20 2.388542347249e-01 21 KSP unpreconditioned resid norm 1.722165062986e-01 true resid norm 1.722165062986e-01 ||r(i)||/||b|| 1.360380133953e-01 My monitor 21 1.722165062986e-01 22 KSP unpreconditioned resid norm 1.161869442046e-01 true resid norm 1.161869442046e-01 ||r(i)||/||b|| 9.177889745747e-02 My monitor 22 1.161869442046e-01 23 KSP unpreconditioned resid norm 6.594339583731e-02 true resid norm 6.594339583731e-02 ||r(i)||/||b|| 5.209029470549e-02 My monitor 23 6.594339583731e-02 24 KSP unpreconditioned resid norm 4.351679748574e-02 true resid norm 4.351679748574e-02 ||r(i)||/||b|| 3.437497837181e-02 My monitor 24 4.351679748573e-02 25 KSP unpreconditioned resid norm 3.847638846864e-02 true resid norm 3.847638846864e-02 ||r(i)||/||b|| 3.039343650847e-02 My monitor 25 3.847638846864e-02 26 KSP unpreconditioned resid norm 2.063424248358e-02 true resid norm 2.063424248358e-02 ||r(i)||/||b|| 1.629949077306e-02 My monitor 26 2.063424248358e-02 27 KSP unpreconditioned resid norm 1.402462240396e-02 true resid norm 1.402462240396e-02 ||r(i)||/||b|| 1.107839086659e-02 My monitor 27 1.402462240396e-02 28 KSP unpreconditioned resid norm 7.732817953098e-03 true resid norm 7.732817953098e-03 ||r(i)||/||b|| 6.108341267025e-03 My monitor 28 7.732817953099e-03 29 KSP unpreconditioned resid norm 5.109464751004e-03 true resid norm 5.109464751004e-03 ||r(i)||/||b|| 4.036090669698e-03 My monitor 29 5.109464751004e-03 30 KSP unpreconditioned resid norm 2.628714079103e-03 true resid norm 2.628714079103e-03 ||r(i)||/||b|| 2.076485284664e-03 My monitor 30 2.628714079103e-03 31 KSP unpreconditioned resid norm 1.211324322673e-03 true resid norm 1.211324322673e-03 ||r(i)||/||b|| 9.568545894671e-04 My monitor 31 1.211324322673e-03 At iteration 32 we see a slight difference in the reported norms 32 KSP unpreconditioned resid norm 5.638897702485e-04 true resid norm 5.638897702485e-04 ||r(i)||/||b|| 4.454302654678e-04 My monitor 32 5.638897702491e-04 Then they continue to be different with more and more digits 33 KSP unpreconditioned resid norm 2.557920120696e-04 true resid norm 2.557920120696e-04 ||r(i)||/||b|| 2.020563412429e-04 My monitor 33 2.557920120695e-04 34 KSP unpreconditioned resid norm 1.249567288159e-04 true resid norm 1.249567288159e-04 ||r(i)||/||b|| 9.870636394758e-05 My monitor 34 1.249567288156e-04 35 KSP unpreconditioned resid norm 6.554146400697e-05 true resid norm 6.554146400697e-05 ||r(i)||/||b|| 5.177279896194e-05 My monitor 35 6.554146400761e-05 36 KSP unpreconditioned resid norm 3.360138566154e-05 true resid norm 3.360138566154e-05 ||r(i)||/||b|| 2.654255303959e-05 My monitor 36 3.360138566057e-05 37 KSP unpreconditioned resid norm 1.963635751089e-05 true resid norm 1.963635751089e-05 ||r(i)||/||b|| 1.551123712537e-05 My monitor 37 1.963635751179e-05 38 KSP unpreconditioned resid norm 1.111922577034e-05 true resid norm 1.111922577034e-05 ||r(i)||/||b|| 8.783347292320e-06 My monitor 38 1.111922577016e-05 Is this the type of discrepancy you're seeing in your code, or are you seeing enormous differences right off the bat? The discrepancy shown above is normal. It arises because KSPSolve_CG() uses PetscCall(VecAXPY(X, a, P)); /* x <- x + ap */ PetscCall(VecAXPY(R, -a, W)); /* r <- r - aw */ to update the solution and the residual. Where W has been computed further up in the code as A*P. If I change KSPSolve_CG() to instead compute R = b - A X explicitly (attached) |
cg.c
Description: Binary data
then the output becomes 32 KSP unpreconditioned resid norm 5.638897702486e-04 true resid norm 5.638897702486e-04 ||r(i)||/||b|| 4.454302654678e-04 My monitor 32 5.638897702486e-04 33 KSP unpreconditioned resid norm 2.557920120698e-04 true resid norm 2.557920120698e-04 ||r(i)||/||b|| 2.020563412431e-04 My monitor 33 2.557920120698e-04 34 KSP unpreconditioned resid norm 1.249567288163e-04 true resid norm 1.249567288163e-04 ||r(i)||/||b|| 9.870636394784e-05 My monitor 34 1.249567288163e-04 35 KSP unpreconditioned resid norm 6.554146400720e-05 true resid norm 6.554146400720e-05 ||r(i)||/||b|| 5.177279896212e-05 My monitor 35 6.554146400720e-05 36 KSP unpreconditioned resid norm 3.360138566157e-05 true resid norm 3.360138566157e-05 ||r(i)||/||b|| 2.654255303962e-05 My monitor 36 3.360138566157e-05 37 KSP unpreconditioned resid norm 1.963635751099e-05 true resid norm 1.963635751099e-05 ||r(i)||/||b|| 1.551123712545e-05 My monitor 37 1.963635751099e-05 38 KSP unpreconditioned resid norm 1.111922577015e-05 true resid norm 1.111922577015e-05 ||r(i)||/||b|| 8.783347292168e-06 My monitor 38 1.111922577015e-05 Now the residual norm printed by -ksp_monitor and MyMonitor are the same to all digits for all iterations. KSPSolve_CG() uses R = R - a W = R - a (A * P) instead of R = B - A*X because it saves a matrix-vector multiply per iteration (generally, for CG the matrix-vector multiply dominates the solution time). One final note. How come "The convergence test is consistent with the 2-norm of KSPBuildResidual" but not computed explicitly from KSPBuildSolution()? That is because the KSPBuildResidual() routine cheats PETSC_INTERN PetscErrorCode KSPBuildResidual_CG(KSP ksp, Vec t, Vec v, Vec *V) { PetscFunctionBegin; PetscCall(VecCopy(ksp->work[0], v)); *V = v; It knows from the KSPSolve_CG() code where the computed residual is stored and gives you that (slightly incorrect :-) one. Now some Krylov methods do not explicitly store (or even compute) the residual vector so they explicitly compute it with PetscErrorCode KSPBuildResidualDefault(KSP ksp, Vec t, Vec v, Vec *V) { Mat Amat, Pmat; PetscFunctionBegin; if (!ksp->pc) PetscCall(KSPGetPC(ksp, &ksp->pc)); PetscCall(PCGetOperators(ksp->pc, &Amat, &Pmat)); PetscCall(KSPBuildSolution(ksp, t, NULL)); PetscCall(KSP_MatMult(ksp, Amat, t, v)); PetscCall(VecAYPX(v, -1.0, ksp->vec_rhs)); *V = v; PetscFunctionReturn(PETSC_SUCCESS); } Barry This is an interesting phenomenon that often is not discussed in elementary introductions to Krylov methods (or even in advanced discussions), so I think I will write an FAQ that explains it for petsc.org
|
