Actual source code: lgc.c
petsc-3.10.5 2019-03-28
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: Concepts: line graph^creating
149: .seealso: PetscDrawLGDestroy(), PetscDrawLGAddPoint(), PetscDrawLGAddCommonPoint(), PetscDrawLGAddPoints(), PetscDrawLGDraw(), PetscDrawLGSave(),
150: PetscDrawLGView(), PetscDrawLGReset(), PetscDrawLGSetDimension(), PetscDrawLGGetDimension(), PetscDrawLGSetLegend(), PetscDrawLGGetAxis(),
151: PetscDrawLGGetDraw(), PetscDrawLGSetUseMarkers(), PetscDrawLGSetLimits(), PetscDrawLGSetColors(), PetscDrawLGSetOptionsPrefix(), PetscDrawLGSetFromOptions()
152: @*/
153: PetscErrorCode PetscDrawLGCreate(PetscDraw draw,PetscInt dim,PetscDrawLG *outlg)
154: {
155: PetscDrawLG lg;
163: PetscHeaderCreate(lg,PETSC_DRAWLG_CLASSID,"DrawLG","Line Graph","Draw",PetscObjectComm((PetscObject)draw),PetscDrawLGDestroy,NULL);
164: PetscLogObjectParent((PetscObject)draw,(PetscObject)lg);
165: PetscDrawLGSetOptionsPrefix(lg,((PetscObject)draw)->prefix);
167: PetscObjectReference((PetscObject)draw);
168: lg->win = draw;
170: lg->view = NULL;
171: lg->destroy = NULL;
172: lg->nopts = 0;
173: lg->dim = dim;
174: lg->xmin = 1.e20;
175: lg->ymin = 1.e20;
176: lg->xmax = -1.e20;
177: lg->ymax = -1.e20;
179: PetscMalloc2(dim*CHUNCKSIZE,&lg->x,dim*CHUNCKSIZE,&lg->y);
180: PetscLogObjectMemory((PetscObject)lg,2*dim*CHUNCKSIZE*sizeof(PetscReal));
182: lg->len = dim*CHUNCKSIZE;
183: lg->loc = 0;
184: lg->use_markers = PETSC_FALSE;
186: PetscDrawAxisCreate(draw,&lg->axis);
187: PetscLogObjectParent((PetscObject)lg,(PetscObject)lg->axis);
189: *outlg = lg;
190: return(0);
191: }
193: /*@
194: PetscDrawLGSetColors - Sets the color of each line graph drawn
196: Logically Collective on PetscDrawLG
198: Input Parameter:
199: + lg - the line graph context.
200: - colors - the colors
202: Level: intermediate
204: Concepts: line graph^setting number of lines
206: .seealso: PetscDrawLGCreate()
208: @*/
209: PetscErrorCode PetscDrawLGSetColors(PetscDrawLG lg,const int colors[])
210: {
217: PetscFree(lg->colors);
218: PetscMalloc1(lg->dim,&lg->colors);
219: PetscMemcpy(lg->colors,colors,lg->dim*sizeof(int));
220: return(0);
221: }
223: /*@C
224: PetscDrawLGSetLegend - sets the names of each curve plotted
226: Logically Collective on PetscDrawLG
228: Input Parameter:
229: + lg - the line graph context.
230: - names - the names for each curve
232: Level: intermediate
234: Notes:
235: Call PetscDrawLGGetAxis() and then change properties of the PetscDrawAxis for detailed control of the plot
237: Concepts: line graph^setting number of lines
239: .seealso: PetscDrawLGGetAxis(), PetscDrawAxis, PetscDrawAxisSetColors(), PetscDrawAxisSetLabels(), PetscDrawAxisSetHoldLimits()
241: @*/
242: PetscErrorCode PetscDrawLGSetLegend(PetscDrawLG lg,const char *const *names)
243: {
245: PetscInt i;
251: if (lg->legend) {
252: for (i=0; i<lg->dim; i++) {
253: PetscFree(lg->legend[i]);
254: }
255: PetscFree(lg->legend);
256: }
257: if (names) {
258: PetscMalloc1(lg->dim,&lg->legend);
259: for (i=0; i<lg->dim; i++) {
260: PetscStrallocpy(names[i],&lg->legend[i]);
261: }
262: }
263: return(0);
264: }
266: /*@
267: PetscDrawLGGetDimension - Change the number of lines that are to be drawn.
269: Not Collective
271: Input Parameter:
272: . lg - the line graph context.
274: Output Parameter:
275: . dim - the number of curves.
277: Level: intermediate
279: Concepts: line graph^setting number of lines
281: .seealso: PetscDrawLGCreate(), PetscDrawLGSetDimension()
283: @*/
284: PetscErrorCode PetscDrawLGGetDimension(PetscDrawLG lg,PetscInt *dim)
285: {
289: *dim = lg->dim;
290: return(0);
291: }
293: /*@
294: PetscDrawLGSetDimension - Change the number of lines that are to be drawn.
296: Logically Collective on PetscDrawLG
298: Input Parameter:
299: + lg - the line graph context.
300: - dim - the number of curves.
302: Level: intermediate
304: Concepts: line graph^setting number of lines
306: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDimension()
307: @*/
308: PetscErrorCode PetscDrawLGSetDimension(PetscDrawLG lg,PetscInt dim)
309: {
311: PetscInt i;
316: if (lg->dim == dim) return(0);
318: PetscFree2(lg->x,lg->y);
319: if (lg->legend) {
320: for (i=0; i<lg->dim; i++) {
321: PetscFree(lg->legend[i]);
322: }
323: PetscFree(lg->legend);
324: }
325: PetscFree(lg->colors);
326: lg->dim = dim;
327: PetscMalloc2(dim*CHUNCKSIZE,&lg->x,dim*CHUNCKSIZE,&lg->y);
328: PetscLogObjectMemory((PetscObject)lg,2*dim*CHUNCKSIZE*sizeof(PetscReal));
329: lg->len = dim*CHUNCKSIZE;
330: return(0);
331: }
334: /*@
335: PetscDrawLGSetLimits - Sets the axis limits for a line graph. If more
336: points are added after this call, the limits will be adjusted to
337: include those additional points.
339: Logically Collective on PetscDrawLG
341: Input Parameters:
342: + xlg - the line graph context
343: - x_min,x_max,y_min,y_max - the limits
345: Level: intermediate
347: Concepts: line graph^setting axis
349: .seealso: PetscDrawLGCreate()
351: @*/
352: PetscErrorCode PetscDrawLGSetLimits(PetscDrawLG lg,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max)
353: {
357: (lg)->xmin = x_min;
358: (lg)->xmax = x_max;
359: (lg)->ymin = y_min;
360: (lg)->ymax = y_max;
361: return(0);
362: }
364: /*@
365: PetscDrawLGReset - Clears line graph to allow for reuse with new data.
367: Logically Collective on PetscDrawLG
369: Input Parameter:
370: . lg - the line graph context.
372: Level: intermediate
374: Concepts: line graph^restarting
376: .seealso: PetscDrawLGCreate()
378: @*/
379: PetscErrorCode PetscDrawLGReset(PetscDrawLG lg)
380: {
383: lg->xmin = 1.e20;
384: lg->ymin = 1.e20;
385: lg->xmax = -1.e20;
386: lg->ymax = -1.e20;
387: lg->loc = 0;
388: lg->nopts = 0;
389: return(0);
390: }
392: /*@
393: PetscDrawLGDestroy - Frees all space taken up by line graph data structure.
395: Collective on PetscDrawLG
397: Input Parameter:
398: . lg - the line graph context
400: Level: intermediate
402: .seealso: PetscDrawLGCreate()
403: @*/
404: PetscErrorCode PetscDrawLGDestroy(PetscDrawLG *lg)
405: {
407: PetscInt i;
410: if (!*lg) return(0);
412: if (--((PetscObject)(*lg))->refct > 0) {*lg = NULL; return(0);}
414: if ((*lg)->legend) {
415: for (i=0; i<(*lg)->dim; i++) {
416: PetscFree((*lg)->legend[i]);
417: }
418: PetscFree((*lg)->legend);
419: }
420: PetscFree((*lg)->colors);
421: PetscFree2((*lg)->x,(*lg)->y);
422: PetscDrawAxisDestroy(&(*lg)->axis);
423: PetscDrawDestroy(&(*lg)->win);
424: PetscHeaderDestroy(lg);
425: return(0);
426: }
427: /*@
428: PetscDrawLGSetUseMarkers - Causes LG to draw a marker for each data-point.
430: Logically Collective on PetscDrawLG
432: Input Parameters:
433: + lg - the linegraph context
434: - flg - should mark each data point
436: Options Database:
437: . -lg_use_markers <true,false>
439: Level: intermediate
441: Concepts: line graph^showing points
443: .seealso: PetscDrawLGCreate()
445: @*/
446: PetscErrorCode PetscDrawLGSetUseMarkers(PetscDrawLG lg,PetscBool flg)
447: {
451: lg->use_markers = flg;
452: return(0);
453: }
455: /*@
456: PetscDrawLGDraw - Redraws a line graph.
458: Collective on PetscDrawLG
460: Input Parameter:
461: . lg - the line graph context
463: Level: intermediate
465: .seealso: PetscDrawSPDraw(), PetscDrawLGSPDraw(), PetscDrawLGReset()
467: @*/
468: PetscErrorCode PetscDrawLGDraw(PetscDrawLG lg)
469: {
470: PetscReal xmin,xmax,ymin,ymax;
472: PetscMPIInt rank;
473: PetscDraw draw;
474: PetscBool isnull;
478: PetscDrawIsNull(lg->win,&isnull);
479: if (isnull) return(0);
480: MPI_Comm_rank(PetscObjectComm((PetscObject)lg),&rank);
482: draw = lg->win;
483: PetscDrawCheckResizedWindow(draw);
484: PetscDrawClear(draw);
486: xmin = lg->xmin; xmax = lg->xmax; ymin = lg->ymin; ymax = lg->ymax;
487: PetscDrawAxisSetLimits(lg->axis,xmin,xmax,ymin,ymax);
488: PetscDrawAxisDraw(lg->axis);
490: PetscDrawCollectiveBegin(draw);
491: if (!rank) {
492: int i,j,dim=lg->dim,nopts=lg->nopts,cl;
493: for (i=0; i<dim; i++) {
494: for (j=1; j<nopts; j++) {
495: cl = lg->colors ? lg->colors[i] : ((PETSC_DRAW_BLACK + i) % PETSC_DRAW_MAXCOLOR);
496: 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);
497: if (lg->use_markers) {PetscDrawMarker(draw,lg->x[j*dim+i],lg->y[j*dim+i],cl);}
498: }
499: }
500: }
501: if (!rank && lg->legend) {
502: int i,dim=lg->dim,cl;
503: PetscReal xl,yl,xr,yr,tw,th;
504: size_t slen,len=0;
505: PetscDrawAxisGetLimits(lg->axis,&xl,&xr,&yl,&yr);
506: PetscDrawStringGetSize(draw,&tw,&th);
507: for (i=0; i<dim; i++) {
508: PetscStrlen(lg->legend[i],&slen);
509: len = PetscMax(len,slen);
510: }
511: xr = xr - 1.5*tw; xl = xr - (len + 7)*tw;
512: yr = yr - 1.0*th; yl = yr - (dim + 1)*th;
513: PetscDrawLine(draw,xl,yl,xr,yl,PETSC_DRAW_BLACK);
514: PetscDrawLine(draw,xr,yl,xr,yr,PETSC_DRAW_BLACK);
515: PetscDrawLine(draw,xr,yr,xl,yr,PETSC_DRAW_BLACK);
516: PetscDrawLine(draw,xl,yr,xl,yl,PETSC_DRAW_BLACK);
517: for (i=0; i<dim; i++) {
518: cl = lg->colors ? lg->colors[i] : (PETSC_DRAW_BLACK + i);
519: PetscDrawLine(draw,xl + 1*tw,yr - (i + 1)*th,xl + 5*tw,yr - (i + 1)*th,cl);
520: PetscDrawString(draw,xl + 6*tw,yr - (i + 1.5)*th,PETSC_DRAW_BLACK,lg->legend[i]);
521: }
522: }
523: PetscDrawCollectiveEnd(draw);
525: PetscDrawFlush(draw);
526: PetscDrawPause(draw);
527: return(0);
528: }
530: /*@
531: PetscDrawLGSave - Saves a drawn image
533: Collective on PetscDrawLG
535: Input Parameter:
536: . lg - The line graph context
538: Level: intermediate
540: Concepts: line graph^saving
542: .seealso: PetscDrawLGCreate(), PetscDrawLGGetDraw(), PetscDrawSetSave(), PetscDrawSave()
543: @*/
544: PetscErrorCode PetscDrawLGSave(PetscDrawLG lg)
545: {
550: PetscDrawSave(lg->win);
551: return(0);
552: }
554: /*@
555: PetscDrawLGView - Prints a line graph.
557: Collective on PetscDrawLG
559: Input Parameter:
560: . lg - the line graph context
562: Level: beginner
564: .seealso: PetscDrawLGCreate()
566: .keywords: draw, line, graph
567: @*/
568: PetscErrorCode PetscDrawLGView(PetscDrawLG lg,PetscViewer viewer)
569: {
570: PetscReal xmin=lg->xmin, xmax=lg->xmax, ymin=lg->ymin, ymax=lg->ymax;
571: PetscInt i, j, dim = lg->dim, nopts = lg->nopts;
577: if (nopts < 1) return(0);
578: if (xmin > xmax || ymin > ymax) return(0);
580: if (!viewer){
581: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)lg),&viewer);
582: }
583: PetscObjectPrintClassNamePrefixType((PetscObject)lg,viewer);
584: for (i = 0; i < dim; i++) {
585: PetscViewerASCIIPrintf(viewer, "Line %D>\n", i);
586: for (j = 0; j < nopts; j++) {
587: PetscViewerASCIIPrintf(viewer, " X: %g Y: %g\n", (double)lg->x[j*dim+i], (double)lg->y[j*dim+i]);
588: }
589: }
590: return(0);
591: }
593: /*@C
594: PetscDrawLGSetOptionsPrefix - Sets the prefix used for searching for all
595: PetscDrawLG options in the database.
597: Logically Collective on PetscDrawLG
599: Input Parameter:
600: + lg - the line graph context
601: - prefix - the prefix to prepend to all option names
603: Level: advanced
605: .keywords: PetscDrawLG, set, options, prefix, database
607: .seealso: PetscDrawLGSetFromOptions(), PetscDrawLGCreate()
608: @*/
609: PetscErrorCode PetscDrawLGSetOptionsPrefix(PetscDrawLG lg,const char prefix[])
610: {
615: PetscObjectSetOptionsPrefix((PetscObject)lg,prefix);
616: return(0);
617: }
619: /*@
620: PetscDrawLGSetFromOptions - Sets options related to the PetscDrawLG
622: Collective on PetscDrawLG
624: Options Database:
626: Level: intermediate
628: Concepts: line graph^creating
630: .seealso: PetscDrawLGDestroy(), PetscDrawLGCreate()
631: @*/
632: PetscErrorCode PetscDrawLGSetFromOptions(PetscDrawLG lg)
633: {
634: PetscErrorCode ierr;
635: PetscBool usemarkers,set;
636: PetscDrawMarkerType markertype;
641: PetscDrawGetMarkerType(lg->win,&markertype);
642: PetscOptionsGetEnum(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_marker_type",PetscDrawMarkerTypes,(PetscEnum*)&markertype,&set);
643: if (set) {
644: PetscDrawLGSetUseMarkers(lg,PETSC_TRUE);
645: PetscDrawSetMarkerType(lg->win,markertype);
646: }
647: usemarkers = lg->use_markers;
648: PetscOptionsGetBool(((PetscObject)lg)->options,((PetscObject)lg)->prefix,"-lg_use_markers",&usemarkers,&set);
649: if (set) {PetscDrawLGSetUseMarkers(lg,usemarkers);}
650: return(0);
651: }