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 coordiante point 15 - y - the new y coordinate point for each curve. 16 17 Level: intermediate 18 19 Concepts: line graph^adding points 20 21 .seealso: PetscDrawLGAddPoints(), PetscDrawLGAddPoint() 22 @*/ 23 PetscErrorCode PetscDrawLGAddCommonPoint(PetscDrawLG lg,const PetscReal x,const PetscReal *y) 24 { 25 PetscErrorCode ierr; 26 PetscInt i; 27 28 PetscFunctionBegin; 29 if (!lg) PetscFunctionReturn(0); 30 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,1); 31 32 if (lg->loc+lg->dim >= lg->len) { /* allocate more space */ 33 PetscReal *tmpx,*tmpy; 34 ierr = PetscMalloc2(lg->len+lg->dim*CHUNCKSIZE,&tmpx,lg->len+lg->dim*CHUNCKSIZE,&tmpy);CHKERRQ(ierr); 35 ierr = PetscLogObjectMemory((PetscObject)lg,2*lg->dim*CHUNCKSIZE*sizeof(PetscReal));CHKERRQ(ierr); 36 ierr = PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 37 ierr = PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 38 ierr = PetscFree2(lg->x,lg->y);CHKERRQ(ierr); 39 lg->x = tmpx; 40 lg->y = tmpy; 41 lg->len += lg->dim*CHUNCKSIZE; 42 } 43 for (i=0; i<lg->dim; i++) { 44 if (x > lg->xmax) lg->xmax = x; 45 if (x < lg->xmin) lg->xmin = x; 46 if (y[i] > lg->ymax) lg->ymax = y[i]; 47 if (y[i] < lg->ymin) lg->ymin = y[i]; 48 49 lg->x[lg->loc] = x; 50 lg->y[lg->loc++] = y[i]; 51 } 52 lg->nopts++; 53 PetscFunctionReturn(0); 54 } 55 56 #undef __FUNCT__ 57 #define __FUNCT__ "PetscDrawLGAddPoint" 58 /*@ 59 PetscDrawLGAddPoint - Adds another point to each of the line graphs. 60 The new point must have an X coordinate larger than the old points. 61 62 Logically Collective on PetscDrawLG 63 64 Input Parameters: 65 + lg - the LineGraph data structure 66 - x, y - the points to two arrays containing the new x and y 67 point for each curve. 68 69 Level: intermediate 70 71 Concepts: line graph^adding points 72 73 .seealso: PetscDrawLGAddPoints(), PetscDrawLGAddCommonPoint() 74 @*/ 75 PetscErrorCode PetscDrawLGAddPoint(PetscDrawLG lg,const PetscReal *x,const PetscReal *y) 76 { 77 PetscErrorCode ierr; 78 PetscInt i; 79 PetscReal xx; 80 81 PetscFunctionBegin; 82 if (!lg) PetscFunctionReturn(0); 83 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,1); 84 85 if (lg->loc+lg->dim >= lg->len) { /* allocate more space */ 86 PetscReal *tmpx,*tmpy; 87 ierr = PetscMalloc2(lg->len+lg->dim*CHUNCKSIZE,&tmpx,lg->len+lg->dim*CHUNCKSIZE,&tmpy);CHKERRQ(ierr); 88 ierr = PetscLogObjectMemory((PetscObject)lg,2*lg->dim*CHUNCKSIZE*sizeof(PetscReal));CHKERRQ(ierr); 89 ierr = PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 90 ierr = PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 91 ierr = PetscFree2(lg->x,lg->y);CHKERRQ(ierr); 92 lg->x = tmpx; 93 lg->y = tmpy; 94 lg->len += lg->dim*CHUNCKSIZE; 95 } 96 for (i=0; i<lg->dim; i++) { 97 if (!x) { 98 xx = lg->nopts; 99 } else { 100 xx = x[i]; 101 } 102 if (xx > lg->xmax) lg->xmax = xx; 103 if (xx < lg->xmin) lg->xmin = xx; 104 if (y[i] > lg->ymax) lg->ymax = y[i]; 105 if (y[i] < lg->ymin) lg->ymin = y[i]; 106 107 lg->x[lg->loc] = xx; 108 lg->y[lg->loc++] = y[i]; 109 } 110 lg->nopts++; 111 PetscFunctionReturn(0); 112 } 113 114 #undef __FUNCT__ 115 #define __FUNCT__ "PetscDrawLGAddPoints" 116 /*@C 117 PetscDrawLGAddPoints - Adds several points to each of the line graphs. 118 The new points must have an X coordinate larger than the old points. 119 120 Logically Collective on PetscDrawLG 121 122 Input Parameters: 123 + lg - the LineGraph data structure 124 . xx,yy - points to two arrays of pointers that point to arrays 125 containing the new x and y points for each curve. 126 - n - number of points being added 127 128 Level: intermediate 129 130 131 Concepts: line graph^adding points 132 133 .seealso: PetscDrawLGAddPoint() 134 @*/ 135 PetscErrorCode PetscDrawLGAddPoints(PetscDrawLG lg,PetscInt n,PetscReal **xx,PetscReal **yy) 136 { 137 PetscErrorCode ierr; 138 PetscInt i,j,k; 139 PetscReal *x,*y; 140 141 PetscFunctionBegin; 142 if (!lg) PetscFunctionReturn(0); 143 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,1); 144 145 if (lg->loc+n*lg->dim >= lg->len) { /* allocate more space */ 146 PetscReal *tmpx,*tmpy; 147 PetscInt chunk = CHUNCKSIZE; 148 149 if (n > chunk) chunk = n; 150 ierr = PetscMalloc2(lg->len+lg->dim*chunk,&tmpx,lg->len+lg->dim*chunk,&tmpy);CHKERRQ(ierr); 151 ierr = PetscLogObjectMemory((PetscObject)lg,2*lg->dim*chunk*sizeof(PetscReal));CHKERRQ(ierr); 152 ierr = PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 153 ierr = PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));CHKERRQ(ierr); 154 ierr = PetscFree2(lg->x,lg->y);CHKERRQ(ierr); 155 lg->x = tmpx; 156 lg->y = tmpy; 157 lg->len += lg->dim*chunk; 158 } 159 for (j=0; j<lg->dim; j++) { 160 x = xx[j]; y = yy[j]; 161 k = lg->loc + j; 162 for (i=0; i<n; i++) { 163 if (x[i] > lg->xmax) lg->xmax = x[i]; 164 if (x[i] < lg->xmin) lg->xmin = x[i]; 165 if (y[i] > lg->ymax) lg->ymax = y[i]; 166 if (y[i] < lg->ymin) lg->ymin = y[i]; 167 168 lg->x[k] = x[i]; 169 lg->y[k] = y[i]; 170 k += lg->dim; 171 } 172 } 173 lg->loc += n*lg->dim; 174 lg->nopts += n; 175 PetscFunctionReturn(0); 176 } 177 178 #undef __FUNCT__ 179 #define __FUNCT__ "PetscDrawLGSetLimits" 180 /*@ 181 PetscDrawLGSetLimits - Sets the axis limits for a line graph. If more 182 points are added after this call, the limits will be adjusted to 183 include those additional points. 184 185 Logically Collective on PetscDrawLG 186 187 Input Parameters: 188 + xlg - the line graph context 189 - x_min,x_max,y_min,y_max - the limits 190 191 Level: intermediate 192 193 Concepts: line graph^setting axis 194 195 @*/ 196 PetscErrorCode PetscDrawLGSetLimits(PetscDrawLG lg,PetscReal x_min,PetscReal x_max,PetscReal y_min,PetscReal y_max) 197 { 198 PetscFunctionBegin; 199 if (!lg) PetscFunctionReturn(0); 200 PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,1); 201 202 (lg)->xmin = x_min; 203 (lg)->xmax = x_max; 204 (lg)->ymin = y_min; 205 (lg)->ymax = y_max; 206 PetscFunctionReturn(0); 207 } 208 209