Actual source code: dagetarray.c

  1: /*$Id: dagetarray.c,v 1.13 2001/08/06 21:18:33 bsmith Exp $*/
  2: 
 3:  #include petscda.h

  5: /*@
  6:    DAVecGetArray - Returns a multiple dimension array that shares data with 
  7:       the underlying vector and is indexed using the global dimensions.

  9:    Not Collective

 11:    Input Parameter:
 12: +  da - the distributed array
 13: -  vec - the vector, either a vector the same size as one obtained with 
 14:          DACreateGlobalVector() or DACreateLocalVector()
 15:    
 16:    Output Parameter:
 17: .  array - the array

 19:    Notes:
 20:     Call DAVecRestoreArray() once you have finished accessing the vector entries.

 22:   Level: intermediate

 24: .keywords: distributed array, get, corners, nodes, local indices, coordinates

 26: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecRestoreArray()
 27: @*/
 28: int DAVecGetArray(DA da,Vec vec,void **array)
 29: {
 30:   int ierr,xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

 33:   DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
 34:   DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
 35:   DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);

 37:   /* Handle case where user passes in global vector as opposed to local */
 38:   VecGetLocalSize(vec,&N);
 39:   if (N == xm*ym*zm*dof) {
 40:     gxm = xm;
 41:     gym = ym;
 42:     gzm = zm;
 43:     gxs = xs;
 44:     gys = ys;
 45:     gzs = zs;
 46:   } else if (N != gxm*gym*gzm*dof) {
 47:     SETERRQ3(1,"Vector local size %d is not compatible with DA local sizes %d %dn",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
 48:   }

 50:   if (dim == 1) {
 51:     VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);
 52:   } else if (dim == 2) {
 53:     VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
 54:   } else if (dim == 3) {
 55:     VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
 56:   } else {
 57:     SETERRQ1(1,"DA dimension not 1, 2, or 3, it is %dn",dim);
 58:   }

 60:   return(0);
 61: }

 63: /*@
 64:    DAVecRestoreArray - Restores a multiple dimension array obtained with DAVecGetArray()

 66:    Not Collective

 68:    Input Parameter:
 69: +  da - the distributed array
 70: .  vec - the vector, either a vector the same size as one obtained with 
 71:          DACreateGlobalVector() or DACreateLocalVector()
 72: -  array - the array

 74:   Level: intermediate

 76: .keywords: distributed array, get, corners, nodes, local indices, coordinates

 78: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecGetArray()
 79: @*/
 80: int DAVecRestoreArray(DA da,Vec vec,void **array)
 81: {
 82:   int ierr,xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

 85:   DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
 86:   DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
 87:   DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);

 89:   /* Handle case where user passes in global vector as opposed to local */
 90:   VecGetLocalSize(vec,&N);
 91:   if (N == xm*ym*zm*dof) {
 92:     gxm = xm;
 93:     gym = ym;
 94:     gzm = zm;
 95:     gxs = xs;
 96:     gys = ys;
 97:     gzs = zs;
 98:   }

100:   if (dim == 1) {
101:     VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);
102:   } else if (dim == 2) {
103:     VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
104:   } else if (dim == 3) {
105:     VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
106:   } else {
107:     SETERRQ1(1,"DA dimension not 1, 2, or 3, it is %dn",dim);
108:   }
109:   return(0);
110: }