Actual source code: xmon.c
1: #include <petsc/private/kspimpl.h>
2: #include <petscdraw.h>
4: PetscErrorCode KSPMonitorLGRange(KSP ksp, PetscInt n, PetscReal rnorm, void *monctx)
5: {
6: PetscDrawLG lg;
7: PetscReal x, y, per;
8: PetscViewer v = (PetscViewer)monctx;
9: static PetscReal prev; /* should be in the context */
10: PetscDraw draw;
12: PetscFunctionBegin;
15: PetscCall(KSPMonitorRange_Private(ksp, n, &per));
16: if (!n) prev = rnorm;
18: PetscCall(PetscViewerDrawGetDrawLG(v, 0, &lg));
19: if (!n) PetscCall(PetscDrawLGReset(lg));
20: PetscCall(PetscDrawLGGetDraw(lg, &draw));
21: PetscCall(PetscDrawSetTitle(draw, "Residual norm"));
22: x = (PetscReal)n;
23: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
24: else y = -15.0;
25: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
26: if (n < 20 || !(n % 5) || ksp->reason) {
27: PetscCall(PetscDrawLGDraw(lg));
28: PetscCall(PetscDrawLGSave(lg));
29: }
31: PetscCall(PetscViewerDrawGetDrawLG(v, 1, &lg));
32: if (!n) PetscCall(PetscDrawLGReset(lg));
33: PetscCall(PetscDrawLGGetDraw(lg, &draw));
34: PetscCall(PetscDrawSetTitle(draw, "% elements > .2*max element"));
35: x = (PetscReal)n;
36: y = 100.0 * per;
37: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
38: if (n < 20 || !(n % 5) || ksp->reason) {
39: PetscCall(PetscDrawLGDraw(lg));
40: PetscCall(PetscDrawLGSave(lg));
41: }
43: PetscCall(PetscViewerDrawGetDrawLG(v, 2, &lg));
44: if (!n) PetscCall(PetscDrawLGReset(lg));
45: PetscCall(PetscDrawLGGetDraw(lg, &draw));
46: PetscCall(PetscDrawSetTitle(draw, "(norm - oldnorm)/oldnorm"));
47: x = (PetscReal)n;
48: y = (prev - rnorm) / prev;
49: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
50: if (n < 20 || !(n % 5) || ksp->reason) {
51: PetscCall(PetscDrawLGDraw(lg));
52: PetscCall(PetscDrawLGSave(lg));
53: }
55: PetscCall(PetscViewerDrawGetDrawLG(v, 3, &lg));
56: if (!n) PetscCall(PetscDrawLGReset(lg));
57: PetscCall(PetscDrawLGGetDraw(lg, &draw));
58: PetscCall(PetscDrawSetTitle(draw, "(norm - oldnorm)/oldnorm*(% > .2 max)"));
59: x = (PetscReal)n;
60: y = (prev - rnorm) / (prev * per);
61: if (n > 5) PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
62: if (n < 20 || !(n % 5) || ksp->reason) {
63: PetscCall(PetscDrawLGDraw(lg));
64: PetscCall(PetscDrawLGSave(lg));
65: }
66: prev = rnorm;
67: PetscFunctionReturn(PETSC_SUCCESS);
68: }