Actual source code: xmon.c
1: #include <petsc/private/kspimpl.h>
2: #include <petscdraw.h>
4: PetscErrorCode KSPMonitorLGCreate(MPI_Comm comm, const char host[], const char label[], const char metric[], PetscInt l, const char *names[], int x, int y, int m, int n, PetscDrawLG *lgctx)
5: {
6: PetscDraw draw;
7: PetscDrawAxis axis;
8: PetscDrawLG lg;
10: PetscFunctionBegin;
11: PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
12: PetscCall(PetscDrawSetFromOptions(draw));
13: PetscCall(PetscDrawLGCreate(draw, l, &lg));
14: if (names) PetscCall(PetscDrawLGSetLegend(lg, names));
15: PetscCall(PetscDrawLGSetFromOptions(lg));
16: PetscCall(PetscDrawLGGetAxis(lg, &axis));
17: PetscCall(PetscDrawAxisSetLabels(axis, "Convergence", "Iteration", metric));
18: PetscCall(PetscDrawDestroy(&draw));
19: *lgctx = lg;
20: PetscFunctionReturn(PETSC_SUCCESS);
21: }
23: PetscErrorCode KSPMonitorLGRange(KSP ksp, PetscInt n, PetscReal rnorm, void *monctx)
24: {
25: PetscDrawLG lg;
26: PetscReal x, y, per;
27: PetscViewer v = (PetscViewer)monctx;
28: static PetscReal prev; /* should be in the context */
29: PetscDraw draw;
31: PetscFunctionBegin;
34: PetscCall(KSPMonitorRange_Private(ksp, n, &per));
35: if (!n) prev = rnorm;
37: PetscCall(PetscViewerDrawGetDrawLG(v, 0, &lg));
38: if (!n) PetscCall(PetscDrawLGReset(lg));
39: PetscCall(PetscDrawLGGetDraw(lg, &draw));
40: PetscCall(PetscDrawSetTitle(draw, "Residual norm"));
41: x = (PetscReal)n;
42: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
43: else y = -15.0;
44: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
45: if (n < 20 || !(n % 5) || ksp->reason) {
46: PetscCall(PetscDrawLGDraw(lg));
47: PetscCall(PetscDrawLGSave(lg));
48: }
50: PetscCall(PetscViewerDrawGetDrawLG(v, 1, &lg));
51: if (!n) PetscCall(PetscDrawLGReset(lg));
52: PetscCall(PetscDrawLGGetDraw(lg, &draw));
53: PetscCall(PetscDrawSetTitle(draw, "% elements > .2*max element"));
54: x = (PetscReal)n;
55: y = 100.0 * per;
56: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
57: if (n < 20 || !(n % 5) || ksp->reason) {
58: PetscCall(PetscDrawLGDraw(lg));
59: PetscCall(PetscDrawLGSave(lg));
60: }
62: PetscCall(PetscViewerDrawGetDrawLG(v, 2, &lg));
63: if (!n) PetscCall(PetscDrawLGReset(lg));
64: PetscCall(PetscDrawLGGetDraw(lg, &draw));
65: PetscCall(PetscDrawSetTitle(draw, "(norm-oldnorm)/oldnorm"));
66: x = (PetscReal)n;
67: y = (prev - rnorm) / prev;
68: PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
69: if (n < 20 || !(n % 5) || ksp->reason) {
70: PetscCall(PetscDrawLGDraw(lg));
71: PetscCall(PetscDrawLGSave(lg));
72: }
74: PetscCall(PetscViewerDrawGetDrawLG(v, 3, &lg));
75: if (!n) PetscCall(PetscDrawLGReset(lg));
76: PetscCall(PetscDrawLGGetDraw(lg, &draw));
77: PetscCall(PetscDrawSetTitle(draw, "(norm -oldnorm)/oldnorm*(% > .2 max)"));
78: x = (PetscReal)n;
79: y = (prev - rnorm) / (prev * per);
80: if (n > 5) PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
81: if (n < 20 || !(n % 5) || ksp->reason) {
82: PetscCall(PetscDrawLGDraw(lg));
83: PetscCall(PetscDrawLGSave(lg));
84: }
86: prev = rnorm;
87: PetscFunctionReturn(PETSC_SUCCESS);
88: }