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