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