1 2 #include <petsc/private/drawimpl.h> /*I "petscdraw.h" I*/ 3 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. 7 8 Logically Collective on PetscDrawLG 9 10 Input Parameters: 11 + lg - the LineGraph data structure 12 . x - the common x coordinate point 13 - y - the new y coordinate point for each curve. 14 15 Level: intermediate 16 17 Note: You must call PetscDrawLGDraw() to display any added points 18 Call PetscDrawLGReset() to remove all points 19 20 .seealso: `PetscDrawLGCreate()`, `PetscDrawLGAddPoints()`, `PetscDrawLGAddPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()` 21 @*/ 22 PetscErrorCode PetscDrawLGAddCommonPoint(PetscDrawLG lg, const PetscReal x, const PetscReal *y) { 23 PetscInt i; 24 25 PetscFunctionBegin; 26 PetscValidHeaderSpecific(lg, PETSC_DRAWLG_CLASSID, 1); 27 28 if (lg->loc + lg->dim >= lg->len) { /* allocate more space */ 29 PetscReal *tmpx, *tmpy; 30 PetscCall(PetscMalloc2(lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpx, lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpy)); 31 PetscCall(PetscLogObjectMemory((PetscObject)lg, 2 * lg->dim * PETSC_DRAW_LG_CHUNK_SIZE * sizeof(PetscReal))); 32 PetscCall(PetscArraycpy(tmpx, lg->x, lg->len)); 33 PetscCall(PetscArraycpy(tmpy, lg->y, lg->len)); 34 PetscCall(PetscFree2(lg->x, lg->y)); 35 lg->x = tmpx; 36 lg->y = tmpy; 37 lg->len += lg->dim * PETSC_DRAW_LG_CHUNK_SIZE; 38 } 39 for (i = 0; i < lg->dim; i++) { 40 if (x > lg->xmax) lg->xmax = x; 41 if (x < lg->xmin) lg->xmin = x; 42 if (y[i] > lg->ymax) lg->ymax = y[i]; 43 if (y[i] < lg->ymin) lg->ymin = y[i]; 44 45 lg->x[lg->loc] = x; 46 lg->y[lg->loc++] = y[i]; 47 } 48 lg->nopts++; 49 PetscFunctionReturn(0); 50 } 51 52 /*@ 53 PetscDrawLGAddPoint - Adds another point to each of the line graphs. 54 The new point must have an X coordinate larger than the old points. 55 56 Logically Collective on PetscDrawLG 57 58 Input Parameters: 59 + lg - the LineGraph data structure 60 - x, y - the points to two arrays containing the new x and y 61 point for each curve. 62 63 Note: You must call PetscDrawLGDraw() to display any added points 64 Call PetscDrawLGReset() to remove all points 65 66 Level: intermediate 67 68 .seealso: `PetscDrawLGCreate()`, `PetscDrawLGAddPoints()`, `PetscDrawLGAddCommonPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()` 69 @*/ 70 PetscErrorCode PetscDrawLGAddPoint(PetscDrawLG lg, const PetscReal *x, const PetscReal *y) { 71 PetscInt i; 72 PetscReal xx; 73 74 PetscFunctionBegin; 75 PetscValidHeaderSpecific(lg, PETSC_DRAWLG_CLASSID, 1); 76 77 if (lg->loc + lg->dim >= lg->len) { /* allocate more space */ 78 PetscReal *tmpx, *tmpy; 79 PetscCall(PetscMalloc2(lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpx, lg->len + lg->dim * PETSC_DRAW_LG_CHUNK_SIZE, &tmpy)); 80 PetscCall(PetscLogObjectMemory((PetscObject)lg, 2 * lg->dim * PETSC_DRAW_LG_CHUNK_SIZE * sizeof(PetscReal))); 81 PetscCall(PetscArraycpy(tmpx, lg->x, lg->len)); 82 PetscCall(PetscArraycpy(tmpy, lg->y, lg->len)); 83 PetscCall(PetscFree2(lg->x, lg->y)); 84 lg->x = tmpx; 85 lg->y = tmpy; 86 lg->len += lg->dim * PETSC_DRAW_LG_CHUNK_SIZE; 87 } 88 for (i = 0; i < lg->dim; i++) { 89 if (!x) { 90 xx = lg->nopts; 91 } else { 92 xx = x[i]; 93 } 94 if (xx > lg->xmax) lg->xmax = xx; 95 if (xx < lg->xmin) lg->xmin = xx; 96 if (y[i] > lg->ymax) lg->ymax = y[i]; 97 if (y[i] < lg->ymin) lg->ymin = y[i]; 98 99 lg->x[lg->loc] = xx; 100 lg->y[lg->loc++] = y[i]; 101 } 102 lg->nopts++; 103 PetscFunctionReturn(0); 104 } 105 106 /*@C 107 PetscDrawLGAddPoints - Adds several points to each of the line graphs. 108 The new points must have an X coordinate larger than the old points. 109 110 Logically Collective on PetscDrawLG 111 112 Input Parameters: 113 + lg - the LineGraph data structure 114 . xx,yy - points to two arrays of pointers that point to arrays 115 containing the new x and y points for each curve. 116 - n - number of points being added 117 118 Level: intermediate 119 120 Note: You must call PetscDrawLGDraw() to display any added points 121 Call PetscDrawLGReset() to remove all points 122 123 .seealso: `PetscDrawLGCreate()`, `PetscDrawLGAddPoint()`, `PetscDrawLGAddCommonPoint()`, `PetscDrawLGReset()`, `PetscDrawLGDraw()` 124 @*/ 125 PetscErrorCode PetscDrawLGAddPoints(PetscDrawLG lg, PetscInt n, PetscReal **xx, PetscReal **yy) { 126 PetscInt i, j, k; 127 PetscReal *x, *y; 128 129 PetscFunctionBegin; 130 PetscValidHeaderSpecific(lg, PETSC_DRAWLG_CLASSID, 1); 131 132 if (lg->loc + n * lg->dim >= lg->len) { /* allocate more space */ 133 PetscReal *tmpx, *tmpy; 134 PetscInt chunk = PETSC_DRAW_LG_CHUNK_SIZE; 135 136 if (n > chunk) chunk = n; 137 PetscCall(PetscMalloc2(lg->len + lg->dim * chunk, &tmpx, lg->len + lg->dim * chunk, &tmpy)); 138 PetscCall(PetscLogObjectMemory((PetscObject)lg, 2 * lg->dim * chunk * sizeof(PetscReal))); 139 PetscCall(PetscArraycpy(tmpx, lg->x, lg->len)); 140 PetscCall(PetscArraycpy(tmpy, lg->y, lg->len)); 141 PetscCall(PetscFree2(lg->x, lg->y)); 142 lg->x = tmpx; 143 lg->y = tmpy; 144 lg->len += lg->dim * chunk; 145 } 146 for (j = 0; j < lg->dim; j++) { 147 x = xx[j]; 148 y = yy[j]; 149 k = lg->loc + j; 150 for (i = 0; i < n; i++) { 151 if (x[i] > lg->xmax) lg->xmax = x[i]; 152 if (x[i] < lg->xmin) lg->xmin = x[i]; 153 if (y[i] > lg->ymax) lg->ymax = y[i]; 154 if (y[i] < lg->ymin) lg->ymin = y[i]; 155 156 lg->x[k] = x[i]; 157 lg->y[k] = y[i]; 158 k += lg->dim; 159 } 160 } 161 lg->loc += n * lg->dim; 162 lg->nopts += n; 163 PetscFunctionReturn(0); 164 } 165