Actual source code: lg.c

petsc-3.7.7 2017-09-25
Report Typos and Errors
  2: #include <../src/sys/classes/draw/utils/lgimpl.h>  /*I   "petscdraw.h"  I*/

  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.

 10:    Logically Collective on PetscDrawLG

 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.

 17:    Level: intermediate

 19:    Concepts: line graph^adding points

 21: .seealso: PetscDrawLGAddPoints(), PetscDrawLGAddPoint()
 22: @*/
 23: PetscErrorCode  PetscDrawLGAddCommonPoint(PetscDrawLG lg,const PetscReal x,const PetscReal *y)
 24: {
 26:   PetscInt       i;


 31:   if (lg->loc+lg->dim >= lg->len) { /* allocate more space */
 32:     PetscReal *tmpx,*tmpy;
 33:     PetscMalloc2(lg->len+lg->dim*CHUNCKSIZE,&tmpx,lg->len+lg->dim*CHUNCKSIZE,&tmpy);
 34:     PetscLogObjectMemory((PetscObject)lg,2*lg->dim*CHUNCKSIZE*sizeof(PetscReal));
 35:     PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));
 36:     PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));
 37:     PetscFree2(lg->x,lg->y);
 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];

 48:     lg->x[lg->loc]   = x;
 49:     lg->y[lg->loc++] = y[i];
 50:   }
 51:   lg->nopts++;
 52:   return(0);
 53: }

 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.

 61:    Logically Collective on PetscDrawLG

 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.

 68:    Level: intermediate

 70:    Concepts: line graph^adding points

 72: .seealso: PetscDrawLGAddPoints(), PetscDrawLGAddCommonPoint()
 73: @*/
 74: PetscErrorCode  PetscDrawLGAddPoint(PetscDrawLG lg,const PetscReal *x,const PetscReal *y)
 75: {
 77:   PetscInt       i;
 78:   PetscReal      xx;


 83:   if (lg->loc+lg->dim >= lg->len) { /* allocate more space */
 84:     PetscReal *tmpx,*tmpy;
 85:     PetscMalloc2(lg->len+lg->dim*CHUNCKSIZE,&tmpx,lg->len+lg->dim*CHUNCKSIZE,&tmpy);
 86:     PetscLogObjectMemory((PetscObject)lg,2*lg->dim*CHUNCKSIZE*sizeof(PetscReal));
 87:     PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));
 88:     PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));
 89:     PetscFree2(lg->x,lg->y);
 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];

105:     lg->x[lg->loc]   = xx;
106:     lg->y[lg->loc++] = y[i];
107:   }
108:   lg->nopts++;
109:   return(0);
110: }

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.

118:    Logically Collective on PetscDrawLG

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

126:    Level: intermediate


129:    Concepts: line graph^adding points

131: .seealso: PetscDrawLGAddPoint()
132: @*/
133: PetscErrorCode  PetscDrawLGAddPoints(PetscDrawLG lg,PetscInt n,PetscReal **xx,PetscReal **yy)
134: {
136:   PetscInt       i,j,k;
137:   PetscReal      *x,*y;


142:   if (lg->loc+n*lg->dim >= lg->len) { /* allocate more space */
143:     PetscReal *tmpx,*tmpy;
144:     PetscInt  chunk = CHUNCKSIZE;

146:     if (n > chunk) chunk = n;
147:     PetscMalloc2(lg->len+lg->dim*chunk,&tmpx,lg->len+lg->dim*chunk,&tmpy);
148:     PetscLogObjectMemory((PetscObject)lg,2*lg->dim*chunk*sizeof(PetscReal));
149:     PetscMemcpy(tmpx,lg->x,lg->len*sizeof(PetscReal));
150:     PetscMemcpy(tmpy,lg->y,lg->len*sizeof(PetscReal));
151:     PetscFree2(lg->x,lg->y);
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];

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:   return(0);
173: }