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;

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

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


 37:   KSPMonitorRange_Private(ksp,n,&per);
 38:   if (!n) prev = rnorm;

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

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

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

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

 91:   prev = rnorm;
 92:   return(0);
 93: }