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: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
12: PetscDrawSetFromOptions(draw);
13: PetscDrawLGCreate(draw,l,&lg);
14: if (names) PetscDrawLGSetLegend(lg,names);
15: PetscDrawLGSetFromOptions(lg);
16: PetscDrawLGGetAxis(lg,&axis);
17: PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric);
18: PetscDrawDestroy(&draw);
19: *lgctx = lg;
20: return 0;
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;
33: KSPMonitorRange_Private(ksp,n,&per);
34: if (!n) prev = rnorm;
36: PetscViewerDrawGetDrawLG(v,0,&lg);
37: if (!n) PetscDrawLGReset(lg);
38: PetscDrawLGGetDraw(lg,&draw);
39: PetscDrawSetTitle(draw,"Residual norm");
40: x = (PetscReal) n;
41: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
42: else y = -15.0;
43: PetscDrawLGAddPoint(lg,&x,&y);
44: if (n < 20 || !(n % 5) || ksp->reason) {
45: PetscDrawLGDraw(lg);
46: PetscDrawLGSave(lg);
47: }
49: PetscViewerDrawGetDrawLG(v,1,&lg);
50: if (!n) PetscDrawLGReset(lg);
51: PetscDrawLGGetDraw(lg,&draw);
52: PetscDrawSetTitle(draw,"% elemts > .2*max elemt");
53: x = (PetscReal) n;
54: y = 100.0*per;
55: PetscDrawLGAddPoint(lg,&x,&y);
56: if (n < 20 || !(n % 5) || ksp->reason) {
57: PetscDrawLGDraw(lg);
58: PetscDrawLGSave(lg);
59: }
61: PetscViewerDrawGetDrawLG(v,2,&lg);
62: if (!n) PetscDrawLGReset(lg);
63: PetscDrawLGGetDraw(lg,&draw);
64: PetscDrawSetTitle(draw,"(norm-oldnorm)/oldnorm");
65: x = (PetscReal) n;
66: y = (prev - rnorm)/prev;
67: PetscDrawLGAddPoint(lg,&x,&y);
68: if (n < 20 || !(n % 5) || ksp->reason) {
69: PetscDrawLGDraw(lg);
70: PetscDrawLGSave(lg);
71: }
73: PetscViewerDrawGetDrawLG(v,3,&lg);
74: if (!n) PetscDrawLGReset(lg);
75: PetscDrawLGGetDraw(lg,&draw);
76: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
77: x = (PetscReal) n;
78: y = (prev - rnorm)/(prev*per);
79: if (n > 5) {
80: PetscDrawLGAddPoint(lg,&x,&y);
81: }
82: if (n < 20 || !(n % 5) || ksp->reason) {
83: PetscDrawLGDraw(lg);
84: PetscDrawLGSave(lg);
85: }
87: prev = rnorm;
88: return 0;
89: }