Actual source code: dalocal.c

petsc-master 2019-08-22
Report Typos and Errors

  2: /*
  3:   Code for manipulating distributed regular arrays in parallel.
  4: */

  6:  #include <petsc/private/dmdaimpl.h>
  7:  #include <petscbt.h>
  8:  #include <petscsf.h>
  9:  #include <petscds.h>
 10:  #include <petscfe.h>

 12: /*
 13:    This allows the DMDA vectors to properly tell MATLAB their dimensions
 14: */
 15: #if defined(PETSC_HAVE_MATLAB_ENGINE)
 16: #include <engine.h>   /* MATLAB include file */
 17: #include <mex.h>      /* MATLAB include file */
 18: static PetscErrorCode  VecMatlabEnginePut_DA2d(PetscObject obj,void *mengine)
 19: {
 21:   PetscInt       n,m;
 22:   Vec            vec = (Vec)obj;
 23:   PetscScalar    *array;
 24:   mxArray        *mat;
 25:   DM             da;

 28:   VecGetDM(vec, &da);
 29:   if (!da) SETERRQ(PetscObjectComm((PetscObject)vec),PETSC_ERR_ARG_WRONGSTATE,"Vector not associated with a DMDA");
 30:   DMDAGetGhostCorners(da,0,0,0,&m,&n,0);

 32:   VecGetArray(vec,&array);
 33: #if !defined(PETSC_USE_COMPLEX)
 34:   mat = mxCreateDoubleMatrix(m,n,mxREAL);
 35: #else
 36:   mat = mxCreateDoubleMatrix(m,n,mxCOMPLEX);
 37: #endif
 38:   PetscArraycpy(mxGetPr(mat),array,n*m);
 39:   PetscObjectName(obj);
 40:   engPutVariable((Engine*)mengine,obj->name,mat);

 42:   VecRestoreArray(vec,&array);
 43:   return(0);
 44: }
 45: #endif

 47: PetscErrorCode  DMCreateLocalVector_DA(DM da,Vec *g)
 48: {
 50:   DM_DA          *dd = (DM_DA*)da->data;

 55:   VecCreate(PETSC_COMM_SELF,g);
 56:   VecSetSizes(*g,dd->nlocal,PETSC_DETERMINE);
 57:   VecSetBlockSize(*g,dd->w);
 58:   VecSetType(*g,da->vectype);
 59:   VecSetDM(*g, da);
 60: #if defined(PETSC_HAVE_MATLAB_ENGINE)
 61:   if (dd->w == 1  && da->dim == 2) {
 62:     PetscObjectComposeFunction((PetscObject)*g,"PetscMatlabEnginePut_C",VecMatlabEnginePut_DA2d);
 63:   }
 64: #endif
 65:   return(0);
 66: }

 68: /*@
 69:   DMDAGetNumCells - Get the number of cells in the local piece of the DMDA. This includes ghost cells.

 71:   Input Parameter:
 72: . dm - The DM object

 74:   Output Parameters:
 75: + numCellsX - The number of local cells in the x-direction
 76: . numCellsY - The number of local cells in the y-direction
 77: . numCellsZ - The number of local cells in the z-direction
 78: - numCells - The number of local cells

 80:   Level: developer

 82: .seealso: DMDAGetCellPoint()
 83: @*/
 84: PetscErrorCode DMDAGetNumCells(DM dm, PetscInt *numCellsX, PetscInt *numCellsY, PetscInt *numCellsZ, PetscInt *numCells)
 85: {
 86:   DM_DA         *da  = (DM_DA*) dm->data;
 87:   const PetscInt dim = dm->dim;
 88:   const PetscInt mx  = (da->Xe - da->Xs)/da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
 89:   const PetscInt nC  = (mx)*(dim > 1 ? (my)*(dim > 2 ? (mz) : 1) : 1);

 93:   if (numCellsX) {
 95:     *numCellsX = mx;
 96:   }
 97:   if (numCellsY) {
 99:     *numCellsY = my;
100:   }
101:   if (numCellsZ) {
103:     *numCellsZ = mz;
104:   }
105:   if (numCells) {
107:     *numCells = nC;
108:   }
109:   return(0);
110: }

112: /*@
113:   DMDAGetCellPoint - Get the DM point corresponding to the tuple (i, j, k) in the DMDA

115:   Input Parameters:
116: + dm - The DM object
117: - i,j,k - The global indices for the cell

119:   Output Parameters:
120: . point - The local DM point

122:   Level: developer

124: .seealso: DMDAGetNumCells()
125: @*/
126: PetscErrorCode DMDAGetCellPoint(DM dm, PetscInt i, PetscInt j, PetscInt k, PetscInt *point)
127: {
128:   const PetscInt dim = dm->dim;
129:   DMDALocalInfo  info;

135:   DMDAGetLocalInfo(dm, &info);
136:   if (dim > 0) {if ((i < info.gxs) || (i >= info.gxs+info.gxm)) SETERRQ3(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "X index %d not in [%d, %d)", i, info.gxs, info.gxs+info.gxm);}
137:   if (dim > 1) {if ((j < info.gys) || (j >= info.gys+info.gym)) SETERRQ3(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Y index %d not in [%d, %d)", j, info.gys, info.gys+info.gym);}
138:   if (dim > 2) {if ((k < info.gzs) || (k >= info.gzs+info.gzm)) SETERRQ3(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Z index %d not in [%d, %d)", k, info.gzs, info.gzs+info.gzm);}
139:   *point = i + (dim > 1 ? (j + (dim > 2 ? k*info.gym : 0))*info.gxm : 0);
140:   return(0);
141: }

143: PetscErrorCode DMDAGetNumVertices(DM dm, PetscInt *numVerticesX, PetscInt *numVerticesY, PetscInt *numVerticesZ, PetscInt *numVertices)
144: {
145:   DM_DA          *da = (DM_DA*) dm->data;
146:   const PetscInt dim = dm->dim;
147:   const PetscInt mx  = (da->Xe - da->Xs)/da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
148:   const PetscInt nVx = mx+1;
149:   const PetscInt nVy = dim > 1 ? (my+1) : 1;
150:   const PetscInt nVz = dim > 2 ? (mz+1) : 1;
151:   const PetscInt nV  = nVx*nVy*nVz;

154:   if (numVerticesX) {
156:     *numVerticesX = nVx;
157:   }
158:   if (numVerticesY) {
160:     *numVerticesY = nVy;
161:   }
162:   if (numVerticesZ) {
164:     *numVerticesZ = nVz;
165:   }
166:   if (numVertices) {
168:     *numVertices = nV;
169:   }
170:   return(0);
171: }

173: PetscErrorCode DMDAGetNumFaces(DM dm, PetscInt *numXFacesX, PetscInt *numXFaces, PetscInt *numYFacesY, PetscInt *numYFaces, PetscInt *numZFacesZ, PetscInt *numZFaces)
174: {
175:   DM_DA          *da = (DM_DA*) dm->data;
176:   const PetscInt dim = dm->dim;
177:   const PetscInt mx  = (da->Xe - da->Xs)/da->w, my = da->Ye - da->Ys, mz = da->Ze - da->Zs;
178:   const PetscInt nxF = (dim > 1 ? (my)*(dim > 2 ? (mz) : 1) : 1);
179:   const PetscInt nXF = (mx+1)*nxF;
180:   const PetscInt nyF = mx*(dim > 2 ? mz : 1);
181:   const PetscInt nYF = dim > 1 ? (my+1)*nyF : 0;
182:   const PetscInt nzF = mx*(dim > 1 ? my : 0);
183:   const PetscInt nZF = dim > 2 ? (mz+1)*nzF : 0;

186:   if (numXFacesX) {
188:     *numXFacesX = nxF;
189:   }
190:   if (numXFaces) {
192:     *numXFaces = nXF;
193:   }
194:   if (numYFacesY) {
196:     *numYFacesY = nyF;
197:   }
198:   if (numYFaces) {
200:     *numYFaces = nYF;
201:   }
202:   if (numZFacesZ) {
204:     *numZFacesZ = nzF;
205:   }
206:   if (numZFaces) {
208:     *numZFaces = nZF;
209:   }
210:   return(0);
211: }

213: PetscErrorCode DMDAGetHeightStratum(DM dm, PetscInt height, PetscInt *pStart, PetscInt *pEnd)
214: {
215:   const PetscInt dim = dm->dim;
216:   PetscInt       nC, nV, nXF, nYF, nZF;

222:   DMDAGetNumCells(dm, NULL, NULL, NULL, &nC);
223:   DMDAGetNumVertices(dm, NULL, NULL, NULL, &nV);
224:   DMDAGetNumFaces(dm, NULL, &nXF, NULL, &nYF, NULL, &nZF);
225:   if (height == 0) {
226:     /* Cells */
227:     if (pStart) *pStart = 0;
228:     if (pEnd)   *pEnd   = nC;
229:   } else if (height == 1) {
230:     /* Faces */
231:     if (pStart) *pStart = nC+nV;
232:     if (pEnd)   *pEnd   = nC+nV+nXF+nYF+nZF;
233:   } else if (height == dim) {
234:     /* Vertices */
235:     if (pStart) *pStart = nC;
236:     if (pEnd)   *pEnd   = nC+nV;
237:   } else if (height < 0) {
238:     /* All points */
239:     if (pStart) *pStart = 0;
240:     if (pEnd)   *pEnd   = nC+nV+nXF+nYF+nZF;
241:   } else SETERRQ1(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "No points of height %d in the DA", height);
242:   return(0);
243: }

245: PetscErrorCode DMDAGetDepthStratum(DM dm, PetscInt depth, PetscInt *pStart, PetscInt *pEnd)
246: {
247:   const PetscInt dim = dm->dim;
248:   PetscInt       nC, nV, nXF, nYF, nZF;

254:   DMDAGetNumCells(dm, NULL, NULL, NULL, &nC);
255:   DMDAGetNumVertices(dm, NULL, NULL, NULL, &nV);
256:   DMDAGetNumFaces(dm, NULL, &nXF, NULL, &nYF, NULL, &nZF);
257:   if (depth == dim) {
258:     /* Cells */
259:     if (pStart) *pStart = 0;
260:     if (pEnd)   *pEnd   = nC;
261:   } else if (depth == dim-1) {
262:     /* Faces */
263:     if (pStart) *pStart = nC+nV;
264:     if (pEnd)   *pEnd   = nC+nV+nXF+nYF+nZF;
265:   } else if (depth == 0) {
266:     /* Vertices */
267:     if (pStart) *pStart = nC;
268:     if (pEnd)   *pEnd   = nC+nV;
269:   } else if (depth < 0) {
270:     /* All points */
271:     if (pStart) *pStart = 0;
272:     if (pEnd)   *pEnd   = nC+nV+nXF+nYF+nZF;
273:   } else SETERRQ1(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "No points of depth %d in the DA", depth);
274:   return(0);
275: }

277: PetscErrorCode DMDAGetConeSize(DM dm, PetscInt p, PetscInt *coneSize)
278: {
279:   const PetscInt dim = dm->dim;
280:   PetscInt       nC, nV, nXF, nYF, nZF;

284:   *coneSize = 0;
285:   DMDAGetNumCells(dm, NULL, NULL, NULL, &nC);
286:   DMDAGetNumVertices(dm, NULL, NULL, NULL, &nV);
287:   DMDAGetNumFaces(dm, NULL, &nXF, NULL, &nYF, NULL, &nZF);
288:   switch (dim) {
289:   case 2:
290:     if (p >= 0) {
291:       if (p < nC) {
292:         *coneSize = 4;
293:       } else if (p < nC+nV) {
294:         *coneSize = 0;
295:       } else if (p < nC+nV+nXF+nYF+nZF) {
296:         *coneSize = 2;
297:       } else SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Point %d should be in [0, %d)", p, nC+nV+nXF+nYF+nZF);
298:     } else SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Negative point %d is invalid", p);
299:     break;
300:   case 3:
301:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Too lazy to do 3D");
302:     break;
303:   }
304:   return(0);
305: }

307: PetscErrorCode DMDAGetCone(DM dm, PetscInt p, PetscInt *cone[])
308: {
309:   const PetscInt dim = dm->dim;
310:   PetscInt       nCx, nCy, nCz, nC, nVx, nVy, nVz, nV, nxF, nyF, nzF, nXF, nYF, nZF;

314:   if (!cone) {DMGetWorkArray(dm, 6, MPIU_INT, cone);}
315:   DMDAGetNumCells(dm, &nCx, &nCy, &nCz, &nC);
316:   DMDAGetNumVertices(dm, &nVx, &nVy, &nVz, &nV);
317:   DMDAGetNumFaces(dm, &nxF, &nXF, &nyF, &nYF, &nzF, &nZF);
318:   switch (dim) {
319:   case 2:
320:     if (p >= 0) {
321:       if (p < nC) {
322:         const PetscInt cy = p / nCx;
323:         const PetscInt cx = p % nCx;

325:         (*cone)[0] = cy*nxF + cx + nC+nV;
326:         (*cone)[1] = cx*nyF + cy + nyF + nC+nV+nXF;
327:         (*cone)[2] = cy*nxF + cx + nxF + nC+nV;
328:         (*cone)[3] = cx*nyF + cy + nC+nV+nXF;
329:         SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Too lazy to do cell cones");
330:       } else if (p < nC+nV) {
331:       } else if (p < nC+nV+nXF) {
332:         const PetscInt fy = (p - nC+nV) / nxF;
333:         const PetscInt fx = (p - nC+nV) % nxF;

335:         (*cone)[0] = fy*nVx + fx + nC;
336:         (*cone)[1] = fy*nVx + fx + 1 + nC;
337:       } else if (p < nC+nV+nXF+nYF) {
338:         const PetscInt fx = (p - nC+nV+nXF) / nyF;
339:         const PetscInt fy = (p - nC+nV+nXF) % nyF;

341:         (*cone)[0] = fy*nVx + fx + nC;
342:         (*cone)[1] = fy*nVx + fx + nVx + nC;
343:       } else SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Point %d should be in [0, %d)", p, nC+nV+nXF+nYF+nZF);
344:     } else SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Negative point %d is invalid", p);
345:     break;
346:   case 3:
347:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "Too lazy to do 3D");
348:     break;
349:   }
350:   return(0);
351: }

353: PetscErrorCode DMDARestoreCone(DM dm, PetscInt p, PetscInt *cone[])
354: {

358:   DMGetWorkArray(dm, 6, MPIU_INT, cone);
359:   return(0);
360: }

362: PetscErrorCode DMDASetVertexCoordinates(DM dm, PetscReal xl, PetscReal xu, PetscReal yl, PetscReal yu, PetscReal zl, PetscReal zu)
363: {
364:   DM_DA         *da = (DM_DA *) dm->data;
365:   Vec            coordinates;
366:   PetscSection   section;
367:   PetscScalar   *coords;
368:   PetscReal      h[3];
369:   PetscInt       dim, size, M, N, P, nVx, nVy, nVz, nV, vStart, vEnd, v, i, j, k;

374:   DMDAGetInfo(dm, &dim, &M, &N, &P, 0,0,0,0,0,0,0,0,0);
375:   if (dim > 3) SETERRQ(PetscObjectComm((PetscObject)dm),PETSC_ERR_PLIB,"The following code only works for dim <= 3");
376:   h[0] = (xu - xl)/M;
377:   h[1] = (yu - yl)/N;
378:   h[2] = (zu - zl)/P;
379:   DMDAGetDepthStratum(dm, 0, &vStart, &vEnd);
380:   DMDAGetNumVertices(dm, &nVx, &nVy, &nVz, &nV);
381:   PetscSectionCreate(PetscObjectComm((PetscObject) dm), &section);
382:   PetscSectionSetNumFields(section, 1);
383:   PetscSectionSetFieldComponents(section, 0, dim);
384:   PetscSectionSetChart(section, vStart, vEnd);
385:   for (v = vStart; v < vEnd; ++v) {
386:     PetscSectionSetDof(section, v, dim);
387:   }
388:   PetscSectionSetUp(section);
389:   PetscSectionGetStorageSize(section, &size);
390:   VecCreateSeq(PETSC_COMM_SELF, size, &coordinates);
391:   PetscObjectSetName((PetscObject)coordinates,"coordinates");
392:   VecGetArray(coordinates, &coords);
393:   for (k = 0; k < nVz; ++k) {
394:     PetscInt ind[3], d, off;

396:     ind[0] = 0;
397:     ind[1] = 0;
398:     ind[2] = k + da->zs;
399:     for (j = 0; j < nVy; ++j) {
400:       ind[1] = j + da->ys;
401:       for (i = 0; i < nVx; ++i) {
402:         const PetscInt vertex = (k*nVy + j)*nVx + i + vStart;

404:         PetscSectionGetOffset(section, vertex, &off);
405:         ind[0] = i + da->xs;
406:         for (d = 0; d < dim; ++d) {
407:           coords[off+d] = h[d]*ind[d];
408:         }
409:       }
410:     }
411:   }
412:   VecRestoreArray(coordinates, &coords);
413:   DMSetCoordinateSection(dm, PETSC_DETERMINE, section);
414:   DMSetCoordinatesLocal(dm, coordinates);
415:   PetscSectionDestroy(&section);
416:   VecDestroy(&coordinates);
417:   return(0);
418: }

420: /* ------------------------------------------------------------------- */

422: /*@C
423:      DMDAGetArray - Gets a work array for a DMDA

425:     Input Parameter:
426: +    da - information about my local patch
427: -    ghosted - do you want arrays for the ghosted or nonghosted patch

429:     Output Parameters:
430: .    vptr - array data structured

432:     Note:  The vector values are NOT initialized and may have garbage in them, so you may need
433:            to zero them.

435:   Level: advanced

437: .seealso: DMDARestoreArray()

439: @*/
440: PetscErrorCode  DMDAGetArray(DM da,PetscBool ghosted,void *vptr)
441: {
443:   PetscInt       j,i,xs,ys,xm,ym,zs,zm;
444:   char           *iarray_start;
445:   void           **iptr = (void**)vptr;
446:   DM_DA          *dd    = (DM_DA*)da->data;

450:   if (ghosted) {
451:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
452:       if (dd->arrayghostedin[i]) {
453:         *iptr                 = dd->arrayghostedin[i];
454:         iarray_start          = (char*)dd->startghostedin[i];
455:         dd->arrayghostedin[i] = NULL;
456:         dd->startghostedin[i] = NULL;

458:         goto done;
459:       }
460:     }
461:     xs = dd->Xs;
462:     ys = dd->Ys;
463:     zs = dd->Zs;
464:     xm = dd->Xe-dd->Xs;
465:     ym = dd->Ye-dd->Ys;
466:     zm = dd->Ze-dd->Zs;
467:   } else {
468:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
469:       if (dd->arrayin[i]) {
470:         *iptr          = dd->arrayin[i];
471:         iarray_start   = (char*)dd->startin[i];
472:         dd->arrayin[i] = NULL;
473:         dd->startin[i] = NULL;

475:         goto done;
476:       }
477:     }
478:     xs = dd->xs;
479:     ys = dd->ys;
480:     zs = dd->zs;
481:     xm = dd->xe-dd->xs;
482:     ym = dd->ye-dd->ys;
483:     zm = dd->ze-dd->zs;
484:   }

486:   switch (da->dim) {
487:   case 1: {
488:     void *ptr;

490:     PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);

492:     ptr   = (void*)(iarray_start - xs*sizeof(PetscScalar));
493:     *iptr = (void*)ptr;
494:     break;
495:   }
496:   case 2: {
497:     void **ptr;

499:     PetscMalloc((ym+1)*sizeof(void*)+xm*ym*sizeof(PetscScalar),&iarray_start);

501:     ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*));
502:     for (j=ys; j<ys+ym; j++) ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs);
503:     *iptr = (void*)ptr;
504:     break;
505:   }
506:   case 3: {
507:     void ***ptr,**bptr;

509:     PetscMalloc((zm+1)*sizeof(void**)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);

511:     ptr  = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*));
512:     bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**));
513:     for (i=zs; i<zs+zm; i++) ptr[i] = bptr + ((i-zs)*ym - ys);
514:     for (i=zs; i<zs+zm; i++) {
515:       for (j=ys; j<ys+ym; j++) {
516:         ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
517:       }
518:     }
519:     *iptr = (void*)ptr;
520:     break;
521:   }
522:   default:
523:     SETERRQ1(PetscObjectComm((PetscObject)da),PETSC_ERR_SUP,"Dimension %D not supported",da->dim);
524:   }

526: done:
527:   /* add arrays to the checked out list */
528:   if (ghosted) {
529:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
530:       if (!dd->arrayghostedout[i]) {
531:         dd->arrayghostedout[i] = *iptr;
532:         dd->startghostedout[i] = iarray_start;
533:         break;
534:       }
535:     }
536:   } else {
537:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
538:       if (!dd->arrayout[i]) {
539:         dd->arrayout[i] = *iptr;
540:         dd->startout[i] = iarray_start;
541:         break;
542:       }
543:     }
544:   }
545:   return(0);
546: }

548: /*@C
549:      DMDARestoreArray - Restores an array of derivative types for a DMDA

551:     Input Parameter:
552: +    da - information about my local patch
553: .    ghosted - do you want arrays for the ghosted or nonghosted patch
554: -    vptr - array data structured

556:      Level: advanced

558: .seealso: DMDAGetArray()

560: @*/
561: PetscErrorCode  DMDARestoreArray(DM da,PetscBool ghosted,void *vptr)
562: {
563:   PetscInt i;
564:   void     **iptr = (void**)vptr,*iarray_start = 0;
565:   DM_DA    *dd    = (DM_DA*)da->data;

569:   if (ghosted) {
570:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
571:       if (dd->arrayghostedout[i] == *iptr) {
572:         iarray_start           = dd->startghostedout[i];
573:         dd->arrayghostedout[i] = NULL;
574:         dd->startghostedout[i] = NULL;
575:         break;
576:       }
577:     }
578:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
579:       if (!dd->arrayghostedin[i]) {
580:         dd->arrayghostedin[i] = *iptr;
581:         dd->startghostedin[i] = iarray_start;
582:         break;
583:       }
584:     }
585:   } else {
586:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
587:       if (dd->arrayout[i] == *iptr) {
588:         iarray_start    = dd->startout[i];
589:         dd->arrayout[i] = NULL;
590:         dd->startout[i] = NULL;
591:         break;
592:       }
593:     }
594:     for (i=0; i<DMDA_MAX_WORK_ARRAYS; i++) {
595:       if (!dd->arrayin[i]) {
596:         dd->arrayin[i] = *iptr;
597:         dd->startin[i] = iarray_start;
598:         break;
599:       }
600:     }
601:   }
602:   return(0);
603: }