Actual source code: lg.c
2: #include <petsc/private/drawimpl.h>
4: /*@
5: PetscDrawLGAddCommonPoint - Adds another point to each of the line graphs. All the points share
6: the same new X coordinate. The new point must have an X coordinate larger than the old points.
8: Logically Collective
10: Input Parameters:
11: + lg - the line graph context
12: . x - the common x coordinate point
13: - y - the new y coordinate point for each curve.
15: Level: intermediate
17: Notes:
18: You must call `PetscDrawLGDraw()` to display any added points
20: Call `PetscDrawLGReset()` to remove all points
22: .seealso: `PetscDrawLG`, `PetscDrawLGCreate()`, `PetscDrawLGAddPoints()`, `PetscDrawLGAddPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()`
23: @*/
24: PetscErrorCode PetscDrawLGAddCommonPoint(PetscDrawLG lg, const PetscReal x, const PetscReal *y)
25: {
26: PetscInt i;
28: PetscFunctionBegin;
31: if (lg->loc + lg->dim >= lg->len) { /* allocate more space */
32: PetscReal *tmpx, *tmpy;
33: PetscCall(PetscMalloc2(lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpx, lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpy));
34: PetscCall(PetscArraycpy(tmpx, lg->x, lg->len));
35: PetscCall(PetscArraycpy(tmpy, lg->y, lg->len));
36: PetscCall(PetscFree2(lg->x, lg->y));
37: lg->x = tmpx;
38: lg->y = tmpy;
39: lg->len += lg->dim * PETSC_DRAW_LG_CHUNK_SIZE;
40: }
41: for (i = 0; i < lg->dim; i++) {
42: if (x > lg->xmax) lg->xmax = x;
43: if (x < lg->xmin) lg->xmin = x;
44: if (y[i] > lg->ymax) lg->ymax = y[i];
45: if (y[i] < lg->ymin) lg->ymin = y[i];
47: lg->x[lg->loc] = x;
48: lg->y[lg->loc++] = y[i];
49: }
50: lg->nopts++;
51: PetscFunctionReturn(PETSC_SUCCESS);
52: }
54: /*@
55: PetscDrawLGAddPoint - Adds another point to each of the line graphs.
56: The new point must have an X coordinate larger than the old points.
58: Logically Collective
60: Input Parameters:
61: + lg - the line graph context
62: . x - array containing the x coordinate for the point on each curve
63: - y - array containing the y coordinate for the point on each curve
65: Level: intermediate
67: Notes:
68: You must call `PetscDrawLGDraw()` to display any added points
70: Call `PetscDrawLGReset()` to remove all points
72: .seealso: `PetscDrawLG`, `PetscDrawLGCreate()`, `PetscDrawLGAddPoints()`, `PetscDrawLGAddCommonPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()`
73: @*/
74: PetscErrorCode PetscDrawLGAddPoint(PetscDrawLG lg, const PetscReal *x, const PetscReal *y)
75: {
76: PetscInt i;
77: PetscReal xx;
79: PetscFunctionBegin;
82: if (lg->loc + lg->dim >= lg->len) { /* allocate more space */
83: PetscReal *tmpx, *tmpy;
84: PetscCall(PetscMalloc2(lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpx, lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpy));
85: PetscCall(PetscArraycpy(tmpx, lg->x, lg->len));
86: PetscCall(PetscArraycpy(tmpy, lg->y, lg->len));
87: PetscCall(PetscFree2(lg->x, lg->y));
88: lg->x = tmpx;
89: lg->y = tmpy;
90: lg->len += lg->dim * PETSC_DRAW_LG_CHUNK_SIZE;
91: }
92: for (i = 0; i < lg->dim; i++) {
93: if (!x) {
94: xx = lg->nopts;
95: } else {
96: xx = x[i];
97: }
98: if (xx > lg->xmax) lg->xmax = xx;
99: if (xx < lg->xmin) lg->xmin = xx;
100: if (y[i] > lg->ymax) lg->ymax = y[i];
101: if (y[i] < lg->ymin) lg->ymin = y[i];
103: lg->x[lg->loc] = xx;
104: lg->y[lg->loc++] = y[i];
105: }
106: lg->nopts++;
107: PetscFunctionReturn(PETSC_SUCCESS);
108: }
110: /*@C
111: PetscDrawLGAddPoints - Adds several points to each of the line graphs.
112: The new points must have an X coordinate larger than the old points.
114: Logically Collective
116: Input Parameters:
117: + lg - the line graph context
118: . xx - array of pointers that point to arrays containing the new x coordinates for each curve.
119: . yy - array of pointers that point to arrays containing the new y points for each curve.
120: - n - number of points being added
122: Level: intermediate
124: Notes:
125: You must call `PetscDrawLGDraw()` to display any added points
127: Call `PetscDrawLGReset()` to remove all points
129: .seealso: `PetscDrawLG`, `PetscDrawLGCreate()`, `PetscDrawLGAddPoint()`, `PetscDrawLGAddCommonPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()`
130: @*/
131: PetscErrorCode PetscDrawLGAddPoints(PetscDrawLG lg, PetscInt n, PetscReal **xx, PetscReal **yy)
132: {
133: PetscInt i, j, k;
134: PetscReal *x, *y;
136: PetscFunctionBegin;
139: if (lg->loc + n * lg->dim >= lg->len) { /* allocate more space */
140: PetscReal *tmpx, *tmpy;
141: PetscInt chunk = PETSC_DRAW_LG_CHUNK_SIZE;
143: if (n > chunk) chunk = n;
144: PetscCall(PetscMalloc2(lg->len + lg->dim * chunk, &tmpx, lg->len + lg->dim * chunk, &tmpy));
145: PetscCall(PetscArraycpy(tmpx, lg->x, lg->len));
146: PetscCall(PetscArraycpy(tmpy, lg->y, lg->len));
147: PetscCall(PetscFree2(lg->x, lg->y));
148: lg->x = tmpx;
149: lg->y = tmpy;
150: lg->len += lg->dim * chunk;
151: }
152: for (j = 0; j < lg->dim; j++) {
153: x = xx[j];
154: y = yy[j];
155: k = lg->loc + j;
156: for (i = 0; i < n; i++) {
157: if (x[i] > lg->xmax) lg->xmax = x[i];
158: if (x[i] < lg->xmin) lg->xmin = x[i];
159: if (y[i] > lg->ymax) lg->ymax = y[i];
160: if (y[i] < lg->ymin) lg->ymin = y[i];
162: lg->x[k] = x[i];
163: lg->y[k] = y[i];
164: k += lg->dim;
165: }
166: }
167: lg->loc += n * lg->dim;
168: lg->nopts += n;
169: PetscFunctionReturn(PETSC_SUCCESS);
170: }