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