Actual source code: xmon.c


  2: #include <petsc/private/kspimpl.h>
  3: #include <petscdraw.h>

  5: 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)
  6: {
  7:   PetscDraw     draw;
  8:   PetscDrawAxis axis;
  9:   PetscDrawLG   lg;

 11:   PetscFunctionBegin;
 12:   PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw));
 13:   PetscCall(PetscDrawSetFromOptions(draw));
 14:   PetscCall(PetscDrawLGCreate(draw, l, &lg));
 15:   if (names) PetscCall(PetscDrawLGSetLegend(lg, names));
 16:   PetscCall(PetscDrawLGSetFromOptions(lg));
 17:   PetscCall(PetscDrawLGGetAxis(lg, &axis));
 18:   PetscCall(PetscDrawAxisSetLabels(axis, "Convergence", "Iteration", metric));
 19:   PetscCall(PetscDrawDestroy(&draw));
 20:   *lgctx = lg;
 21:   PetscFunctionReturn(PETSC_SUCCESS);
 22: }

 24: PetscErrorCode KSPMonitorLGRange(KSP ksp, PetscInt n, PetscReal rnorm, void *monctx)
 25: {
 26:   PetscDrawLG      lg;
 27:   PetscReal        x, y, per;
 28:   PetscViewer      v = (PetscViewer)monctx;
 29:   static PetscReal prev; /* should be in the context */
 30:   PetscDraw        draw;

 32:   PetscFunctionBegin;

 35:   PetscCall(KSPMonitorRange_Private(ksp, n, &per));
 36:   if (!n) prev = rnorm;

 38:   PetscCall(PetscViewerDrawGetDrawLG(v, 0, &lg));
 39:   if (!n) PetscCall(PetscDrawLGReset(lg));
 40:   PetscCall(PetscDrawLGGetDraw(lg, &draw));
 41:   PetscCall(PetscDrawSetTitle(draw, "Residual norm"));
 42:   x = (PetscReal)n;
 43:   if (rnorm > 0.0) y = PetscLog10Real(rnorm);
 44:   else y = -15.0;
 45:   PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
 46:   if (n < 20 || !(n % 5) || ksp->reason) {
 47:     PetscCall(PetscDrawLGDraw(lg));
 48:     PetscCall(PetscDrawLGSave(lg));
 49:   }

 51:   PetscCall(PetscViewerDrawGetDrawLG(v, 1, &lg));
 52:   if (!n) PetscCall(PetscDrawLGReset(lg));
 53:   PetscCall(PetscDrawLGGetDraw(lg, &draw));
 54:   PetscCall(PetscDrawSetTitle(draw, "% elements > .2*max element"));
 55:   x = (PetscReal)n;
 56:   y = 100.0 * per;
 57:   PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
 58:   if (n < 20 || !(n % 5) || ksp->reason) {
 59:     PetscCall(PetscDrawLGDraw(lg));
 60:     PetscCall(PetscDrawLGSave(lg));
 61:   }

 63:   PetscCall(PetscViewerDrawGetDrawLG(v, 2, &lg));
 64:   if (!n) PetscCall(PetscDrawLGReset(lg));
 65:   PetscCall(PetscDrawLGGetDraw(lg, &draw));
 66:   PetscCall(PetscDrawSetTitle(draw, "(norm-oldnorm)/oldnorm"));
 67:   x = (PetscReal)n;
 68:   y = (prev - rnorm) / prev;
 69:   PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
 70:   if (n < 20 || !(n % 5) || ksp->reason) {
 71:     PetscCall(PetscDrawLGDraw(lg));
 72:     PetscCall(PetscDrawLGSave(lg));
 73:   }

 75:   PetscCall(PetscViewerDrawGetDrawLG(v, 3, &lg));
 76:   if (!n) PetscCall(PetscDrawLGReset(lg));
 77:   PetscCall(PetscDrawLGGetDraw(lg, &draw));
 78:   PetscCall(PetscDrawSetTitle(draw, "(norm -oldnorm)/oldnorm*(% > .2 max)"));
 79:   x = (PetscReal)n;
 80:   y = (prev - rnorm) / (prev * per);
 81:   if (n > 5) PetscCall(PetscDrawLGAddPoint(lg, &x, &y));
 82:   if (n < 20 || !(n % 5) || ksp->reason) {
 83:     PetscCall(PetscDrawLGDraw(lg));
 84:     PetscCall(PetscDrawLGSave(lg));
 85:   }

 87:   prev = rnorm;
 88:   PetscFunctionReturn(PETSC_SUCCESS);
 89: }