Actual source code: xmon.c
petsc-3.9.4 2018-09-11
2: #include <petsc/private/kspimpl.h>
3: #include <petscdraw.h>
5: /*@C
6: KSPMonitorLGResidualNormCreate - Creates a line graph context for use with
7: KSP to monitor convergence of preconditioned residual norms.
9: Collective on KSP
11: Input Parameters:
12: + comm - communicator context
13: . host - the X display to open, or null for the local machine
14: . label - the title to put in the title bar
15: . x, y - the screen coordinates of the upper left coordinate of
16: the window
17: - m, n - the screen width and height in pixels
19: Output Parameter:
20: . lgctx - the drawing context
22: Options Database Key:
23: . -ksp_monitor_lg_residualnorm - Sets line graph monitor
25: Notes:
26: Use PetscDrawLGDestroy() to destroy this line graph.
28: Level: intermediate
30: .keywords: KSP, monitor, line graph, residual, create
32: .seealso: KSPMonitorSet(), KSPMonitorLGTrueResidualCreate()
33: @*/
34: PetscErrorCode KSPMonitorLGResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
35: {
36: PetscDraw draw;
38: PetscDrawAxis axis;
39: PetscDrawLG lg;
42: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
43: PetscDrawSetFromOptions(draw);
44: PetscDrawLGCreate(draw,1,&lg);
45: PetscDrawLGSetFromOptions(lg);
46: PetscDrawLGGetAxis(lg,&axis);
47: PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");
48: PetscDrawDestroy(&draw);
49: *lgctx = lg;
50: return(0);
51: }
53: PetscErrorCode KSPMonitorLGResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *ctx)
54: {
55: PetscDrawLG lg = (PetscDrawLG) ctx;
56: PetscReal x,y;
61: if (!n) {PetscDrawLGReset(lg);}
62: x = (PetscReal) n;
63: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
64: else y = -15.0;
65: PetscDrawLGAddPoint(lg,&x,&y);
66: if (n <= 20 || !(n % 5) || ksp->reason) {
67: PetscDrawLGDraw(lg);
68: PetscDrawLGSave(lg);
69: }
70: return(0);
71: }
73: extern PetscErrorCode KSPMonitorRange_Private(KSP,PetscInt,PetscReal*);
74: PetscErrorCode KSPMonitorLGRange(KSP ksp,PetscInt n,PetscReal rnorm,void *monctx)
75: {
76: PetscDrawLG lg;
77: PetscErrorCode ierr;
78: PetscReal x,y,per;
79: PetscViewer v = (PetscViewer)monctx;
80: static PetscReal prev; /* should be in the context */
81: PetscDraw draw;
86: KSPMonitorRange_Private(ksp,n,&per);
87: if (!n) prev = rnorm;
89: PetscViewerDrawGetDrawLG(v,0,&lg);
90: if (!n) {PetscDrawLGReset(lg);}
91: PetscDrawLGGetDraw(lg,&draw);
92: PetscDrawSetTitle(draw,"Residual norm");
93: x = (PetscReal) n;
94: if (rnorm > 0.0) y = PetscLog10Real(rnorm);
95: else y = -15.0;
96: PetscDrawLGAddPoint(lg,&x,&y);
97: if (n < 20 || !(n % 5) || ksp->reason) {
98: PetscDrawLGDraw(lg);
99: PetscDrawLGSave(lg);
100: }
102: PetscViewerDrawGetDrawLG(v,1,&lg);
103: if (!n) {PetscDrawLGReset(lg);}
104: PetscDrawLGGetDraw(lg,&draw);
105: PetscDrawSetTitle(draw,"% elemts > .2*max elemt");
106: x = (PetscReal) n;
107: y = 100.0*per;
108: PetscDrawLGAddPoint(lg,&x,&y);
109: if (n < 20 || !(n % 5) || ksp->reason) {
110: PetscDrawLGDraw(lg);
111: PetscDrawLGSave(lg);
112: }
114: PetscViewerDrawGetDrawLG(v,2,&lg);
115: if (!n) {PetscDrawLGReset(lg);}
116: PetscDrawLGGetDraw(lg,&draw);
117: PetscDrawSetTitle(draw,"(norm-oldnorm)/oldnorm");
118: x = (PetscReal) n;
119: y = (prev - rnorm)/prev;
120: PetscDrawLGAddPoint(lg,&x,&y);
121: if (n < 20 || !(n % 5) || ksp->reason) {
122: PetscDrawLGDraw(lg);
123: PetscDrawLGSave(lg);
124: }
126: PetscViewerDrawGetDrawLG(v,3,&lg);
127: if (!n) {PetscDrawLGReset(lg);}
128: PetscDrawLGGetDraw(lg,&draw);
129: PetscDrawSetTitle(draw,"(norm -oldnorm)/oldnorm*(% > .2 max)");
130: x = (PetscReal) n;
131: y = (prev - rnorm)/(prev*per);
132: if (n > 5) {
133: PetscDrawLGAddPoint(lg,&x,&y);
134: }
135: if (n < 20 || !(n % 5) || ksp->reason) {
136: PetscDrawLGDraw(lg);
137: PetscDrawLGSave(lg);
138: }
140: prev = rnorm;
141: return(0);
142: }
144: /*@C
145: KSPMonitorLGTrueResidualNormCreate - Creates a line graph context for use with
146: KSP to monitor convergence of true residual norms (as opposed to
147: preconditioned residual norms).
149: Collective on KSP
151: Input Parameters:
152: + comm - communicator context
153: . host - the X display to open, or null for the local machine
154: . label - the title to put in the title bar
155: . x, y - the screen coordinates of the upper left coordinate of
156: the window
157: - m, n - the screen width and height in pixels
159: Output Parameter:
160: . lgctx - the drawing context
162: Options Database Key:
163: . -ksp_monitor_lg_true_residualnorm - Sets true line graph monitor
165: Notes:
166: Use PetscDrawLGDestroy() to destroy this line graph.
168: Level: intermediate
170: .keywords: KSP, monitor, line graph, residual, create, true
172: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGResidualNormCreate()
173: @*/
174: PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm comm,const char host[],const char label[],int x,int y,int m,int n,PetscDrawLG *lgctx)
175: {
176: PetscDraw draw;
178: PetscDrawAxis axis;
179: PetscDrawLG lg;
180: const char *names[] = {"Preconditioned","True"};
183: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
184: PetscDrawSetFromOptions(draw);
185: PetscDrawLGCreate(draw,2,&lg);
186: PetscDrawLGSetLegend(lg,names);
187: PetscDrawLGSetFromOptions(lg);
188: PetscDrawLGGetAxis(lg,&axis);
189: PetscDrawAxisSetLabels(axis,"Convergence","Iteration","Residual Norm");
190: PetscDrawDestroy(&draw);
191: *lgctx = lg;
192: return(0);
193: }
195: PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *ctx)
196: {
197: PetscDrawLG lg = (PetscDrawLG) ctx;
198: PetscReal x[2],y[2],scnorm;
199: Vec resid,work;
205: if (!n) {PetscDrawLGReset(lg);}
206: x[0] = x[1] = (PetscReal) n;
207: if (rnorm > 0.0) y[0] = PetscLog10Real(rnorm);
208: else y[0] = -15.0;
209: VecDuplicate(ksp->vec_rhs,&work);
210: KSPBuildResidual(ksp,NULL,work,&resid);
211: VecNorm(resid,NORM_2,&scnorm);
212: VecDestroy(&work);
213: if (scnorm > 0.0) y[1] = PetscLog10Real(scnorm);
214: else y[1] = -15.0;
215: PetscDrawLGAddPoint(lg,x,y);
216: if (n <= 20 || !(n % 5)) {
217: PetscDrawLGDraw(lg);
218: PetscDrawLGSave(lg);
219: }
220: return(0);
221: }