Actual source code: lgc.c
2: #include <petscviewer.h>
3: #include <petsc/private/drawimpl.h>
4: PetscClassId PETSC_DRAWLG_CLASSID = 0;
6: /*@
7: PetscDrawLGGetAxis - Gets the axis context associated with a line graph.
8: This is useful if one wants to change some axis property, such as
9: labels, color, etc. The axis context should not be destroyed by the
10: application code.
12: Not Collective, if PetscDrawLG is parallel then PetscDrawAxis is parallel
14: Input Parameter:
15: . lg - the line graph context
17: Output Parameter:
18: . axis - the axis context
20: Level: advanced
22: .seealso: PetscDrawLGCreate(), PetscDrawAxis
24: @*/
25: PetscErrorCode PetscDrawLGGetAxis(PetscDrawLG lg,PetscDrawAxis *axis)
26: {
29: *axis = lg->axis;
30: return 0;
31: }
33: /*@
34: PetscDrawLGGetDraw - Gets the draw context associated with a line graph.
36: Not Collective, if PetscDrawLG is parallel then PetscDraw is parallel
38: Input Parameter:
39: . lg - the line graph context
41: Output Parameter:
42: . draw - the draw context
44: Level: intermediate
46: .seealso: PetscDrawLGCreate(), PetscDraw
47: @*/
48: PetscErrorCode PetscDrawLGGetDraw(PetscDrawLG lg,PetscDraw *draw)
49: {
52: *draw = lg->win;
53: return 0;
54: }
56: /*@
57: PetscDrawLGSPDraw - Redraws a line graph.
59: Collective on PetscDrawLG
61: Input Parameter:
62: . lg - the line graph context
64: Level: intermediate
66: .seealso: PetscDrawLGDraw(), PetscDrawSPDraw()
68: Developer Notes:
69: This code cheats and uses the fact that the LG and SP structs are the same
71: @*/
72: PetscErrorCode PetscDrawLGSPDraw(PetscDrawLG lg,PetscDrawSP spin)
73: {
74: PetscDrawLG sp = (PetscDrawLG)spin;
75: PetscReal xmin,xmax,ymin,ymax;
76: PetscBool isnull;
77: PetscMPIInt rank;
78: PetscDraw draw;
83: PetscDrawIsNull(lg->win,&isnull);
84: if (isnull) return 0;
85: MPI_Comm_rank(PetscObjectComm((PetscObject)lg),&rank);
87: draw = lg->win;
88: PetscDrawCheckResizedWindow(draw);
89: PetscDrawClear(draw);
91: xmin = PetscMin(lg->xmin,sp->xmin); ymin = PetscMin(lg->ymin,sp->ymin);
92: xmax = PetscMax(lg->xmax,sp->xmax); ymax = PetscMax(lg->ymax,sp->ymax);
93: PetscDrawAxisSetLimits(lg->axis,xmin,xmax,ymin,ymax);
94: PetscDrawAxisDraw(lg->axis);
96: PetscDrawCollectiveBegin(draw);
97: if (rank == 0) {
98: int i,j,dim,nopts;
99: dim = lg->dim;
100: nopts = lg->nopts;
101: for (i=0; i<dim; i++) {
102: for (j=1; j<nopts; j++) {
103: PetscDrawLine(draw,lg->x[(j-1)*dim+i],lg->y[(j-1)*dim+i],lg->x[j*dim+i],lg->y[j*dim+i],PETSC_DRAW_BLACK+i);
104: if (lg->use_markers) {
105: PetscDrawMarker(draw,lg->x[j*dim+i],lg->y[j*dim+i],PETSC_DRAW_RED);
106: }
107: }
108: }
109: dim = sp->dim;
110: nopts = sp->nopts;
111: for (i=0; i<dim; i++) {
112: for (j=0; j<nopts; j++) {
113: PetscDrawMarker(draw,sp->x[j*dim+i],sp->y[j*dim+i],PETSC_DRAW_RED);
114: }
115: }
116: }
117: PetscDrawCollectiveEnd(draw);
119: PetscDrawFlush(draw);
120: PetscDrawPause(draw);
121: return 0;
122: }
124: /*@
125: PetscDrawLGCreate - Creates a line graph data structure.
127: Collective on PetscDraw
129: Input Parameters:
130: + draw - the window where the graph will be made.
131: - dim - the number of curves which will be drawn
133: Output Parameters:
134: . outlg - the line graph context
136: Level: intermediate
138: Notes:
139: The MPI communicator that owns the PetscDraw owns this PetscDrawLG, but the calls to set options and add points are ignored on all processes except the
140: zeroth MPI process in the communicator. All MPI processes in the communicator must call PetscDrawLGDraw() to display the updated graph.
142: .seealso: PetscDrawLGDestroy(), PetscDrawLGAddPoint(), PetscDrawLGAddCommonPoint(), PetscDrawLGAddPoints(), PetscDrawLGDraw(), PetscDrawLGSave(),
143: PetscDrawLGView(), PetscDrawLGReset(), PetscDrawLGSetDimension(), PetscDrawLGGetDimension(), PetscDrawLGSetLegend(), PetscDrawLGGetAxis(),
144: PetscDrawLGGetDraw(), PetscDrawLGSetUseMarkers(), PetscDrawLGSetLimits(), PetscDrawLGSetColors(), PetscDrawLGSetOptionsPrefix(), PetscDrawLGSetFromOptions()
145: @*/
146: PetscErrorCode PetscDrawLGCreate(PetscDraw draw,PetscInt dim,PetscDrawLG *outlg)
147: {
148: PetscDrawLG lg;
154: PetscHeaderCreate(lg,PETSC_DRAWLG_CLASSID,"DrawLG","Line Graph","Draw",PetscObjectComm((PetscObject)draw),PetscDrawLGDestroy,NULL);
155: PetscLogObjectParent((PetscObject)draw,(PetscObject)lg);
156: PetscDrawLGSetOptionsPrefix(lg,((PetscObject)draw)->prefix);
158: PetscObjectReference((PetscObject)draw);
159: lg->win = draw;
161: lg->view = NULL;
162: lg->destroy = NULL;
163: lg->nopts = 0;
164: lg->dim = dim;
165: lg->xmin = 1.e20;
166: lg->ymin = 1.e20;
167: lg->xmax = -1.e20;
168: lg->ymax = -1.e20;
170: PetscMalloc2(dim*PETSC_DRAW_LG_CHUNK_SIZE,&lg->x,dim*PETSC_DRAW_LG_CHUNK_SIZE,&lg->y);
171: PetscLogObjectMemory((PetscObject)lg,2*dim*PETSC_DRAW_LG_CHUNK_SIZE*sizeof(PetscReal));
173: lg->len = dim*PETSC_DRAW_LG_CHUNK_SIZE;
174: lg->loc = 0;
175: lg->use_markers = PETSC_FALSE;
177: PetscDrawAxisCreate(draw,&lg->axis);
178: PetscLogObjectParent((PetscObject)lg,(PetscObject)lg->axis);
180: *outlg = lg;
181: return 0;
182: }
184: /*@
185: PetscDrawLGSetColors - Sets the color of each line graph drawn
187: Logically Collective on PetscDrawLG
189: Input Parameters:
190: + lg - the line graph context.
191: - colors - the colors
193: Level: intermediate
195: .seealso: PetscDrawLGCreate()
197: @*/
198: PetscErrorCode PetscDrawLGSetColors(PetscDrawLG lg,const int colors[])
199: {
203: PetscFree(lg->colors);
204: PetscMalloc1(lg->dim,&lg->colors);
205: PetscArraycpy(lg->colors,colors,lg->dim);
206: return 0;
207: }
209: /*@C
210: PetscDrawLGSetLegend - sets the names of each curve plotted
212: Logically Collective on PetscDrawLG
214: Input Parameters:
215: + lg - the line graph context.
216: - names - the names for each curve
218: Level: intermediate
220: Notes:
221: Call PetscDrawLGGetAxis() and then change properties of the PetscDrawAxis for detailed control of the plot
223: .seealso: PetscDrawLGGetAxis(), PetscDrawAxis, PetscDrawAxisSetColors(), PetscDrawAxisSetLabels(), PetscDrawAxisSetHoldLimits()
225: @*/
226: PetscErrorCode PetscDrawLGSetLegend(PetscDrawLG lg,const char *const *names)
227: {
228: PetscInt i;
233: if (lg->legend) {
234: for (i=0; i<lg->dim; i++) {
235: PetscFree(lg->legend[i]);
236: }
237: PetscFree(lg->legend);
238: }
239: if (names) {
240: PetscMalloc1(lg->dim,&lg->legend);
241: for (i=0; i<lg->dim; i++) {
242: PetscStrallocpy(names[i],&lg->legend[i]);
243: }
244: }
245: return 0;
246: }
248: /*@
249: PetscDrawLGGetDimension - Change the number of lines that are to be drawn.
251: Not Collective
253: Input Parameter:
254: . lg - the line graph context.
256: Output Parameter:
257: . dim - the number of curves.
259: Level: intermediate
261: .seealso: PetscDrawLGCreate(), PetscDrawLGSetDimension()
263: @*/
264: PetscErrorCode PetscDrawLGGetDimension(PetscDrawLG lg,PetscInt *dim)
265: {
268: *dim = lg->dim;
269: return 0;
270: }
272: /*@
273: PetscDrawLGSetDimension - Change the number of lines that are to be drawn.
275: Logically Collective on PetscDrawLG
277: Input Parameters:
278: + lg - the line graph context.
279: - dim - the number of curves.
281: Level: intermediate
283: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDimension()
284: @*/
285: PetscErrorCode PetscDrawLGSetDimension(PetscDrawLG lg,PetscInt dim)
286: {
287: PetscInt i;
291: if (lg->dim == dim) return 0;
293: PetscFree2(lg->x,lg->y);
294: if (lg->legend) {
295: for (i=0; i<lg->dim; i++) {
296: PetscFree(lg->legend[i]);
297: }
298: PetscFree(lg->legend);
299: }
300: PetscFree(lg->colors);
301: lg->dim = dim;
302: PetscMalloc2(dim*PETSC_DRAW_LG_CHUNK_SIZE,&lg->x,dim*PETSC_DRAW_LG_CHUNK_SIZE,&lg->y);
303: PetscLogObjectMemory((PetscObject)lg,2*dim*PETSC_DRAW_LG_CHUNK_SIZE*sizeof(PetscReal));
304: lg->len = dim*PETSC_DRAW_LG_CHUNK_SIZE;
305: return 0;
306: }
308: /*@
309: PetscDrawLGSetLimits - Sets the axis limits for a line graph. If more
310: points are added after this call, the limits will be adjusted to
311: include those additional points.
313: Logically Collective on PetscDrawLG
315: Input Parameters:
316: + xlg - the line graph context
317: - x_min,x_max,y_min,y_max - the limits
319: Level: intermediate
321: .seealso: PetscDrawLGCreate()
323: @*/
324: PetscErrorCode PetscDrawLGSetLimits(PetscDrawLG lg,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
325: {
328: (lg)->xmin = x_min;
329: (lg)->xmax = x_max;
330: (lg)->ymin = y_min;
331: (lg)->ymax = y_max;
332: return 0;
333: }
335: /*@
336: PetscDrawLGReset - Clears line graph to allow for reuse with new data.
338: Logically Collective on PetscDrawLG
340: Input Parameter:
341: . lg - the line graph context.
343: Level: intermediate
345: .seealso: PetscDrawLGCreate()
346: @*/
347: PetscErrorCode PetscDrawLGReset(PetscDrawLG lg)
348: {
350: lg->xmin = 1.e20;
351: lg->ymin = 1.e20;
352: lg->xmax = -1.e20;
353: lg->ymax = -1.e20;
354: lg->loc = 0;
355: lg->nopts = 0;
356: return 0;
357: }
359: /*@
360: PetscDrawLGDestroy - Frees all space taken up by line graph data structure.
362: Collective on PetscDrawLG
364: Input Parameter:
365: . lg - the line graph context
367: Level: intermediate
369: .seealso: PetscDrawLGCreate()
370: @*/
371: PetscErrorCode PetscDrawLGDestroy(PetscDrawLG *lg)
372: {
373: PetscInt i;
375: if (!*lg) return 0;
377: if (--((PetscObject)(*lg))->refct > 0) {*lg = NULL; return 0;}
379: if ((*lg)->legend) {
380: for (i=0; i<(*lg)->dim; i++) {
381: PetscFree((*lg)->legend[i]);
382: }
383: PetscFree((*lg)->legend);
384: }
385: PetscFree((*lg)->colors);
386: PetscFree2((*lg)->x,(*lg)->y);
387: PetscDrawAxisDestroy(&(*lg)->axis);
388: PetscDrawDestroy(&(*lg)->win);
389: PetscHeaderDestroy(lg);
390: return 0;
391: }
392: /*@
393: PetscDrawLGSetUseMarkers - Causes LG to draw a marker for each data-point.
395: Logically Collective on PetscDrawLG
397: Input Parameters:
398: + lg - the linegraph context
399: - flg - should mark each data point
401: Options Database:
402: . -lg_use_markers <true,false> - true means the graphPetscDrawLG draws a marker for each point
404: Level: intermediate
406: .seealso: PetscDrawLGCreate()
407: @*/
408: PetscErrorCode PetscDrawLGSetUseMarkers(PetscDrawLG lg,PetscBool flg)
409: {
412: lg->use_markers = flg;
413: return 0;
414: }
416: /*@
417: PetscDrawLGDraw - Redraws a line graph.
419: Collective on PetscDrawLG
421: Input Parameter:
422: . lg - the line graph context
424: Level: intermediate
426: .seealso: PetscDrawSPDraw(), PetscDrawLGSPDraw(), PetscDrawLGReset()
427: @*/
428: PetscErrorCode PetscDrawLGDraw(PetscDrawLG lg)
429: {
430: PetscReal xmin,xmax,ymin,ymax;
431: PetscMPIInt rank;
432: PetscDraw draw;
433: PetscBool isnull;
437: PetscDrawIsNull(lg->win,&isnull);
438: if (isnull) return 0;
439: MPI_Comm_rank(PetscObjectComm((PetscObject)lg),&rank);
441: draw = lg->win;
442: PetscDrawCheckResizedWindow(draw);
443: PetscDrawClear(draw);
445: xmin = lg->xmin; xmax = lg->xmax; ymin = lg->ymin; ymax = lg->ymax;
446: PetscDrawAxisSetLimits(lg->axis,xmin,xmax,ymin,ymax);
447: PetscDrawAxisDraw(lg->axis);
449: PetscDrawCollectiveBegin(draw);
450: if (rank == 0) {
451: int i,j,dim=lg->dim,nopts=lg->nopts,cl;
452: for (i=0; i<dim; i++) {
453: for (j=1; j<nopts; j++) {
454: cl = lg->colors ? lg->colors[i] : ((PETSC_DRAW_BLACK + i) % PETSC_DRAW_MAXCOLOR);
455: PetscDrawLine(draw,lg->x[(j-1)*dim+i],lg->y[(j-1)*dim+i],lg->x[j*dim+i],lg->y[j*dim+i],cl);
456: if (lg->use_markers) PetscDrawMarker(draw,lg->x[j*dim+i],lg->y[j*dim+i],cl);
457: }
458: }
459: }
460: if (rank == 0 && lg->legend) {
461: PetscBool right = PETSC_FALSE;
462: int i,dim=lg->dim,cl;
463: PetscReal xl,yl,xr,yr,tw,th;
464: size_t slen,len=0;
465: PetscDrawAxisGetLimits(lg->axis,&xl,&xr,&yl,&yr);
466: PetscDrawStringGetSize(draw,&tw,&th);
467: for (i=0; i<dim; i++) {
468: PetscStrlen(lg->legend[i],&slen);
469: len = PetscMax(len,slen);
470: }
471: if (right) {
472: xr = xr - 1.5*tw; xl = xr - (len + 7)*tw;
473: } else {
474: xl = xl + 1.5*tw; xr = xl + (len + 7)*tw;
475: }
476: yr = yr - 1.0*th; yl = yr - (dim + 1)*th;
477: PetscDrawLine(draw,xl,yl,xr,yl,PETSC_DRAW_BLACK);
478: PetscDrawLine(draw,xr,yl,xr,yr,PETSC_DRAW_BLACK);
479: PetscDrawLine(draw,xr,yr,xl,yr,PETSC_DRAW_BLACK);
480: PetscDrawLine(draw,xl,yr,xl,yl,PETSC_DRAW_BLACK);
481: for (i=0; i<dim; i++) {
482: cl = lg->colors ? lg->colors[i] : (PETSC_DRAW_BLACK + i);
483: PetscDrawLine(draw,xl + 1*tw,yr - (i + 1)*th,xl + 5*tw,yr - (i + 1)*th,cl);
484: PetscDrawString(draw,xl + 6*tw,yr - (i + 1.5)*th,PETSC_DRAW_BLACK,lg->legend[i]);
485: }
486: }
487: PetscDrawCollectiveEnd(draw);
489: PetscDrawFlush(draw);
490: PetscDrawPause(draw);
491: return 0;
492: }
494: /*@
495: PetscDrawLGSave - Saves a drawn image
497: Collective on PetscDrawLG
499: Input Parameter:
500: . lg - The line graph context
502: Level: intermediate
504: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDraw(), PetscDrawSetSave(), PetscDrawSave()
505: @*/
506: PetscErrorCode PetscDrawLGSave(PetscDrawLG lg)
507: {
509: PetscDrawSave(lg->win);
510: return 0;
511: }
513: /*@
514: PetscDrawLGView - Prints a line graph.
516: Collective on PetscDrawLG
518: Input Parameter:
519: . lg - the line graph context
521: Level: beginner
523: .seealso: PetscDrawLGCreate()
525: @*/
526: PetscErrorCode PetscDrawLGView(PetscDrawLG lg,PetscViewer viewer)
527: {
528: PetscReal xmin=lg->xmin, xmax=lg->xmax, ymin=lg->ymin, ymax=lg->ymax;
529: PetscInt i, j, dim = lg->dim, nopts = lg->nopts;
533: if (nopts < 1) return 0;
534: if (xmin > xmax || ymin > ymax) return 0;
536: if (!viewer) {
537: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)lg),&viewer);
538: }
539: PetscObjectPrintClassNamePrefixType((PetscObject)lg,viewer);
540: for (i = 0; i < dim; i++) {
541: PetscViewerASCIIPrintf(viewer, "Line %" PetscInt_FMT ">\n", i);
542: for (j = 0; j < nopts; j++) {
543: PetscViewerASCIIPrintf(viewer, " X: %g Y: %g\n", (double)lg->x[j*dim+i], (double)lg->y[j*dim+i]);
544: }
545: }
546: return 0;
547: }
549: /*@C
550: PetscDrawLGSetOptionsPrefix - Sets the prefix used for searching for all
551: PetscDrawLG options in the database.
553: Logically Collective on PetscDrawLG
555: Input Parameters:
556: + lg - the line graph context
557: - prefix - the prefix to prepend to all option names
559: Level: advanced
561: .seealso: PetscDrawLGSetFromOptions(), PetscDrawLGCreate()
562: @*/
563: PetscErrorCode PetscDrawLGSetOptionsPrefix(PetscDrawLG lg,const char prefix[])
564: {
566: PetscObjectSetOptionsPrefix((PetscObject)lg,prefix);
567: return 0;
568: }
570: /*@
571: PetscDrawLGSetFromOptions - Sets options related to the PetscDrawLG
573: Collective on PetscDrawLG
575: Options Database:
577: Level: intermediate
579: .seealso: PetscDrawLGDestroy(), PetscDrawLGCreate()
580: @*/
581: PetscErrorCode PetscDrawLGSetFromOptions(PetscDrawLG lg)
582: {
583: PetscBool usemarkers,set;
584: PetscDrawMarkerType markertype;
588: PetscDrawGetMarkerType(lg->win,&markertype);
589: PetscOptionsGetEnum(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_marker_type",PetscDrawMarkerTypes,(PetscEnum*)&markertype,&set);
590: if (set) {
591: PetscDrawLGSetUseMarkers(lg,PETSC_TRUE);
592: PetscDrawSetMarkerType(lg->win,markertype);
593: }
594: usemarkers = lg->use_markers;
595: PetscOptionsGetBool(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_use_markers",&usemarkers,&set);
596: if (set) PetscDrawLGSetUseMarkers(lg,usemarkers);
597: return 0;
598: }