Actual source code: lgc.c
2: #include <petscviewer.h>
3: #include <../src/sys/classes/draw/utils/lgimpl.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: {
30: *axis = lg->axis;
31: return(0);
32: }
34: /*@
35: PetscDrawLGGetDraw - Gets the draw context associated with a line graph.
37: Not Collective, if PetscDrawLG is parallel then PetscDraw is parallel
39: Input Parameter:
40: . lg - the line graph context
42: Output Parameter:
43: . draw - the draw context
45: Level: intermediate
47: .seealso: PetscDrawLGCreate(), PetscDraw
48: @*/
49: PetscErrorCode PetscDrawLGGetDraw(PetscDrawLG lg,PetscDraw *draw)
50: {
54: *draw = lg->win;
55: return(0);
56: }
59: /*@
60: PetscDrawLGSPDraw - Redraws a line graph.
62: Collective on PetscDrawLG
64: Input Parameter:
65: . lg - the line graph context
67: Level: intermediate
69: .seealso: PetscDrawLGDraw(), PetscDrawSPDraw()
71: Developer Notes:
72: This code cheats and uses the fact that the LG and SP structs are the same
74: @*/
75: PetscErrorCode PetscDrawLGSPDraw(PetscDrawLG lg,PetscDrawSP spin)
76: {
77: PetscDrawLG sp = (PetscDrawLG)spin;
78: PetscReal xmin,xmax,ymin,ymax;
80: PetscBool isnull;
81: PetscMPIInt rank;
82: PetscDraw draw;
87: PetscDrawIsNull(lg->win,&isnull);
88: if (isnull) return(0);
89: MPI_Comm_rank(PetscObjectComm((PetscObject)lg),&rank);
91: draw = lg->win;
92: PetscDrawCheckResizedWindow(draw);
93: PetscDrawClear(draw);
95: xmin = PetscMin(lg->xmin,sp->xmin); ymin = PetscMin(lg->ymin,sp->ymin);
96: xmax = PetscMax(lg->xmax,sp->xmax); ymax = PetscMax(lg->ymax,sp->ymax);
97: PetscDrawAxisSetLimits(lg->axis,xmin,xmax,ymin,ymax);
98: PetscDrawAxisDraw(lg->axis);
100: PetscDrawCollectiveBegin(draw);
101: if (!rank) {
102: int i,j,dim,nopts;
103: dim = lg->dim;
104: nopts = lg->nopts;
105: for (i=0; i<dim; i++) {
106: for (j=1; j<nopts; j++) {
107: 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);
108: if (lg->use_markers) {
109: PetscDrawMarker(draw,lg->x[j*dim+i],lg->y[j*dim+i],PETSC_DRAW_RED);
110: }
111: }
112: }
113: dim = sp->dim;
114: nopts = sp->nopts;
115: for (i=0; i<dim; i++) {
116: for (j=0; j<nopts; j++) {
117: PetscDrawMarker(draw,sp->x[j*dim+i],sp->y[j*dim+i],PETSC_DRAW_RED);
118: }
119: }
120: }
121: PetscDrawCollectiveEnd(draw);
123: PetscDrawFlush(draw);
124: PetscDrawPause(draw);
125: return(0);
126: }
129: /*@
130: PetscDrawLGCreate - Creates a line graph data structure.
132: Collective on PetscDraw
134: Input Parameters:
135: + draw - the window where the graph will be made.
136: - dim - the number of curves which will be drawn
138: Output Parameters:
139: . outlg - the line graph context
141: Level: intermediate
143: Notes:
144: 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
145: zeroth MPI process in the communicator. All MPI processes in the communicator must call PetscDrawLGDraw() to display the updated graph.
147: .seealso: PetscDrawLGDestroy(), PetscDrawLGAddPoint(), PetscDrawLGAddCommonPoint(), PetscDrawLGAddPoints(), PetscDrawLGDraw(), PetscDrawLGSave(),
148: PetscDrawLGView(), PetscDrawLGReset(), PetscDrawLGSetDimension(), PetscDrawLGGetDimension(), PetscDrawLGSetLegend(), PetscDrawLGGetAxis(),
149: PetscDrawLGGetDraw(), PetscDrawLGSetUseMarkers(), PetscDrawLGSetLimits(), PetscDrawLGSetColors(), PetscDrawLGSetOptionsPrefix(), PetscDrawLGSetFromOptions()
150: @*/
151: PetscErrorCode PetscDrawLGCreate(PetscDraw draw,PetscInt dim,PetscDrawLG *outlg)
152: {
153: PetscDrawLG lg;
161: PetscHeaderCreate(lg,PETSC_DRAWLG_CLASSID,"DrawLG","Line Graph","Draw",PetscObjectComm((PetscObject)draw),PetscDrawLGDestroy,NULL);
162: PetscLogObjectParent((PetscObject)draw,(PetscObject)lg);
163: PetscDrawLGSetOptionsPrefix(lg,((PetscObject)draw)->prefix);
165: PetscObjectReference((PetscObject)draw);
166: lg->win = draw;
168: lg->view = NULL;
169: lg->destroy = NULL;
170: lg->nopts = 0;
171: lg->dim = dim;
172: lg->xmin = 1.e20;
173: lg->ymin = 1.e20;
174: lg->xmax = -1.e20;
175: lg->ymax = -1.e20;
177: PetscMalloc2(dim*CHUNCKSIZE,&lg->x,dim*CHUNCKSIZE,&lg->y);
178: PetscLogObjectMemory((PetscObject)lg,2*dim*CHUNCKSIZE*sizeof(PetscReal));
180: lg->len = dim*CHUNCKSIZE;
181: lg->loc = 0;
182: lg->use_markers = PETSC_FALSE;
184: PetscDrawAxisCreate(draw,&lg->axis);
185: PetscLogObjectParent((PetscObject)lg,(PetscObject)lg->axis);
187: *outlg = lg;
188: return(0);
189: }
191: /*@
192: PetscDrawLGSetColors - Sets the color of each line graph drawn
194: Logically Collective on PetscDrawLG
196: Input Parameter:
197: + lg - the line graph context.
198: - colors - the colors
200: Level: intermediate
202: .seealso: PetscDrawLGCreate()
204: @*/
205: PetscErrorCode PetscDrawLGSetColors(PetscDrawLG lg,const int colors[])
206: {
213: PetscFree(lg->colors);
214: PetscMalloc1(lg->dim,&lg->colors);
215: PetscArraycpy(lg->colors,colors,lg->dim);
216: return(0);
217: }
219: /*@C
220: PetscDrawLGSetLegend - sets the names of each curve plotted
222: Logically Collective on PetscDrawLG
224: Input Parameter:
225: + lg - the line graph context.
226: - names - the names for each curve
228: Level: intermediate
230: Notes:
231: Call PetscDrawLGGetAxis() and then change properties of the PetscDrawAxis for detailed control of the plot
233: .seealso: PetscDrawLGGetAxis(), PetscDrawAxis, PetscDrawAxisSetColors(), PetscDrawAxisSetLabels(), PetscDrawAxisSetHoldLimits()
235: @*/
236: PetscErrorCode PetscDrawLGSetLegend(PetscDrawLG lg,const char *const *names)
237: {
239: PetscInt i;
245: if (lg->legend) {
246: for (i=0; i<lg->dim; i++) {
247: PetscFree(lg->legend[i]);
248: }
249: PetscFree(lg->legend);
250: }
251: if (names) {
252: PetscMalloc1(lg->dim,&lg->legend);
253: for (i=0; i<lg->dim; i++) {
254: PetscStrallocpy(names[i],&lg->legend[i]);
255: }
256: }
257: return(0);
258: }
260: /*@
261: PetscDrawLGGetDimension - Change the number of lines that are to be drawn.
263: Not Collective
265: Input Parameter:
266: . lg - the line graph context.
268: Output Parameter:
269: . dim - the number of curves.
271: Level: intermediate
273: .seealso: PetscDrawLGCreate(), PetscDrawLGSetDimension()
275: @*/
276: PetscErrorCode PetscDrawLGGetDimension(PetscDrawLG lg,PetscInt *dim)
277: {
281: *dim = lg->dim;
282: return(0);
283: }
285: /*@
286: PetscDrawLGSetDimension - Change the number of lines that are to be drawn.
288: Logically Collective on PetscDrawLG
290: Input Parameter:
291: + lg - the line graph context.
292: - dim - the number of curves.
294: Level: intermediate
296: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDimension()
297: @*/
298: PetscErrorCode PetscDrawLGSetDimension(PetscDrawLG lg,PetscInt dim)
299: {
301: PetscInt i;
306: if (lg->dim == dim) return(0);
308: PetscFree2(lg->x,lg->y);
309: if (lg->legend) {
310: for (i=0; i<lg->dim; i++) {
311: PetscFree(lg->legend[i]);
312: }
313: PetscFree(lg->legend);
314: }
315: PetscFree(lg->colors);
316: lg->dim = dim;
317: PetscMalloc2(dim*CHUNCKSIZE,&lg->x,dim*CHUNCKSIZE,&lg->y);
318: PetscLogObjectMemory((PetscObject)lg,2*dim*CHUNCKSIZE*sizeof(PetscReal));
319: lg->len = dim*CHUNCKSIZE;
320: return(0);
321: }
324: /*@
325: PetscDrawLGSetLimits - Sets the axis limits for a line graph. If more
326: points are added after this call, the limits will be adjusted to
327: include those additional points.
329: Logically Collective on PetscDrawLG
331: Input Parameters:
332: + xlg - the line graph context
333: - x_min,x_max,y_min,y_max - the limits
335: Level: intermediate
337: .seealso: PetscDrawLGCreate()
339: @*/
340: PetscErrorCode PetscDrawLGSetLimits(PetscDrawLG lg,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
341: {
345: (lg)->xmin = x_min;
346: (lg)->xmax = x_max;
347: (lg)->ymin = y_min;
348: (lg)->ymax = y_max;
349: return(0);
350: }
352: /*@
353: PetscDrawLGReset - Clears line graph to allow for reuse with new data.
355: Logically Collective on PetscDrawLG
357: Input Parameter:
358: . lg - the line graph context.
360: Level: intermediate
362: .seealso: PetscDrawLGCreate()
363: @*/
364: PetscErrorCode PetscDrawLGReset(PetscDrawLG lg)
365: {
368: lg->xmin = 1.e20;
369: lg->ymin = 1.e20;
370: lg->xmax = -1.e20;
371: lg->ymax = -1.e20;
372: lg->loc = 0;
373: lg->nopts = 0;
374: return(0);
375: }
377: /*@
378: PetscDrawLGDestroy - Frees all space taken up by line graph data structure.
380: Collective on PetscDrawLG
382: Input Parameter:
383: . lg - the line graph context
385: Level: intermediate
387: .seealso: PetscDrawLGCreate()
388: @*/
389: PetscErrorCode PetscDrawLGDestroy(PetscDrawLG *lg)
390: {
392: PetscInt i;
395: if (!*lg) return(0);
397: if (--((PetscObject)(*lg))->refct > 0) {*lg = NULL; return(0);}
399: if ((*lg)->legend) {
400: for (i=0; i<(*lg)->dim; i++) {
401: PetscFree((*lg)->legend[i]);
402: }
403: PetscFree((*lg)->legend);
404: }
405: PetscFree((*lg)->colors);
406: PetscFree2((*lg)->x,(*lg)->y);
407: PetscDrawAxisDestroy(&(*lg)->axis);
408: PetscDrawDestroy(&(*lg)->win);
409: PetscHeaderDestroy(lg);
410: return(0);
411: }
412: /*@
413: PetscDrawLGSetUseMarkers - Causes LG to draw a marker for each data-point.
415: Logically Collective on PetscDrawLG
417: Input Parameters:
418: + lg - the linegraph context
419: - flg - should mark each data point
421: Options Database:
422: . -lg_use_markers <true,false> - true means the graphPetscDrawLG draws a marker for each point
424: Level: intermediate
426: .seealso: PetscDrawLGCreate()
427: @*/
428: PetscErrorCode PetscDrawLGSetUseMarkers(PetscDrawLG lg,PetscBool flg)
429: {
433: lg->use_markers = flg;
434: return(0);
435: }
437: /*@
438: PetscDrawLGDraw - Redraws a line graph.
440: Collective on PetscDrawLG
442: Input Parameter:
443: . lg - the line graph context
445: Level: intermediate
447: .seealso: PetscDrawSPDraw(), PetscDrawLGSPDraw(), PetscDrawLGReset()
448: @*/
449: PetscErrorCode PetscDrawLGDraw(PetscDrawLG lg)
450: {
451: PetscReal xmin,xmax,ymin,ymax;
453: PetscMPIInt rank;
454: PetscDraw draw;
455: PetscBool isnull;
459: PetscDrawIsNull(lg->win,&isnull);
460: if (isnull) return(0);
461: MPI_Comm_rank(PetscObjectComm((PetscObject)lg),&rank);
463: draw = lg->win;
464: PetscDrawCheckResizedWindow(draw);
465: PetscDrawClear(draw);
467: xmin = lg->xmin; xmax = lg->xmax; ymin = lg->ymin; ymax = lg->ymax;
468: PetscDrawAxisSetLimits(lg->axis,xmin,xmax,ymin,ymax);
469: PetscDrawAxisDraw(lg->axis);
471: PetscDrawCollectiveBegin(draw);
472: if (!rank) {
473: int i,j,dim=lg->dim,nopts=lg->nopts,cl;
474: for (i=0; i<dim; i++) {
475: for (j=1; j<nopts; j++) {
476: cl = lg->colors ? lg->colors[i] : ((PETSC_DRAW_BLACK + i) % PETSC_DRAW_MAXCOLOR);
477: 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);
478: if (lg->use_markers) {PetscDrawMarker(draw,lg->x[j*dim+i],lg->y[j*dim+i],cl);}
479: }
480: }
481: }
482: if (!rank && lg->legend) {
483: int i,dim=lg->dim,cl;
484: PetscReal xl,yl,xr,yr,tw,th;
485: size_t slen,len=0;
486: PetscDrawAxisGetLimits(lg->axis,&xl,&xr,&yl,&yr);
487: PetscDrawStringGetSize(draw,&tw,&th);
488: for (i=0; i<dim; i++) {
489: PetscStrlen(lg->legend[i],&slen);
490: len = PetscMax(len,slen);
491: }
492: xr = xr - 1.5*tw; xl = xr - (len + 7)*tw;
493: yr = yr - 1.0*th; yl = yr - (dim + 1)*th;
494: PetscDrawLine(draw,xl,yl,xr,yl,PETSC_DRAW_BLACK);
495: PetscDrawLine(draw,xr,yl,xr,yr,PETSC_DRAW_BLACK);
496: PetscDrawLine(draw,xr,yr,xl,yr,PETSC_DRAW_BLACK);
497: PetscDrawLine(draw,xl,yr,xl,yl,PETSC_DRAW_BLACK);
498: for (i=0; i<dim; i++) {
499: cl = lg->colors ? lg->colors[i] : (PETSC_DRAW_BLACK + i);
500: PetscDrawLine(draw,xl + 1*tw,yr - (i + 1)*th,xl + 5*tw,yr - (i + 1)*th,cl);
501: PetscDrawString(draw,xl + 6*tw,yr - (i + 1.5)*th,PETSC_DRAW_BLACK,lg->legend[i]);
502: }
503: }
504: PetscDrawCollectiveEnd(draw);
506: PetscDrawFlush(draw);
507: PetscDrawPause(draw);
508: return(0);
509: }
511: /*@
512: PetscDrawLGSave - Saves a drawn image
514: Collective on PetscDrawLG
516: Input Parameter:
517: . lg - The line graph context
519: Level: intermediate
521: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDraw(), PetscDrawSetSave(), PetscDrawSave()
522: @*/
523: PetscErrorCode PetscDrawLGSave(PetscDrawLG lg)
524: {
529: PetscDrawSave(lg->win);
530: return(0);
531: }
533: /*@
534: PetscDrawLGView - Prints a line graph.
536: Collective on PetscDrawLG
538: Input Parameter:
539: . lg - the line graph context
541: Level: beginner
543: .seealso: PetscDrawLGCreate()
545: @*/
546: PetscErrorCode PetscDrawLGView(PetscDrawLG lg,PetscViewer viewer)
547: {
548: PetscReal xmin=lg->xmin, xmax=lg->xmax, ymin=lg->ymin, ymax=lg->ymax;
549: PetscInt i, j, dim = lg->dim, nopts = lg->nopts;
555: if (nopts < 1) return(0);
556: if (xmin > xmax || ymin > ymax) return(0);
558: if (!viewer){
559: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)lg),&viewer);
560: }
561: PetscObjectPrintClassNamePrefixType((PetscObject)lg,viewer);
562: for (i = 0; i < dim; i++) {
563: PetscViewerASCIIPrintf(viewer, "Line %D>\n", i);
564: for (j = 0; j < nopts; j++) {
565: PetscViewerASCIIPrintf(viewer, " X: %g Y: %g\n", (double)lg->x[j*dim+i], (double)lg->y[j*dim+i]);
566: }
567: }
568: return(0);
569: }
571: /*@C
572: PetscDrawLGSetOptionsPrefix - Sets the prefix used for searching for all
573: PetscDrawLG options in the database.
575: Logically Collective on PetscDrawLG
577: Input Parameter:
578: + lg - the line graph context
579: - prefix - the prefix to prepend to all option names
581: Level: advanced
583: .seealso: PetscDrawLGSetFromOptions(), PetscDrawLGCreate()
584: @*/
585: PetscErrorCode PetscDrawLGSetOptionsPrefix(PetscDrawLG lg,const char prefix[])
586: {
591: PetscObjectSetOptionsPrefix((PetscObject)lg,prefix);
592: return(0);
593: }
595: /*@
596: PetscDrawLGSetFromOptions - Sets options related to the PetscDrawLG
598: Collective on PetscDrawLG
600: Options Database:
602: Level: intermediate
604: .seealso: PetscDrawLGDestroy(), PetscDrawLGCreate()
605: @*/
606: PetscErrorCode PetscDrawLGSetFromOptions(PetscDrawLG lg)
607: {
608: PetscErrorCode ierr;
609: PetscBool usemarkers,set;
610: PetscDrawMarkerType markertype;
615: PetscDrawGetMarkerType(lg->win,&markertype);
616: PetscOptionsGetEnum(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_marker_type",PetscDrawMarkerTypes,(PetscEnum*)&markertype,&set);
617: if (set) {
618: PetscDrawLGSetUseMarkers(lg,PETSC_TRUE);
619: PetscDrawSetMarkerType(lg->win,markertype);
620: }
621: usemarkers = lg->use_markers;
622: PetscOptionsGetBool(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_use_markers",&usemarkers,&set);
623: if (set) {PetscDrawLGSetUseMarkers(lg,usemarkers);}
624: return(0);
625: }