Actual source code: dagetarray.c

  1: 
 2:  #include petscda.h

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

 10:    Not Collective

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

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

 23:   Level: intermediate

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

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

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

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

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

 62:   return(0);
 63: }

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

 70:    Not Collective

 72:    Input Parameter:
 73: +  da - the distributed array
 74: .  vec - the vector, either a vector the same size as one obtained with 
 75:          DACreateGlobalVector() or DACreateLocalVector()
 76: -  array - the array

 78:   Level: intermediate

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

 82: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecGetArray()
 83: @*/
 84: PetscErrorCode DAVecRestoreArray(DA da,Vec vec,void *array)
 85: {
 87:   PetscInt  xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

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

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

105:   if (dim == 1) {
106:     VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);
107:   } else if (dim == 2) {
108:     VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
109:   } else if (dim == 3) {
110:     VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
111:   } else {
112:     SETERRQ1(PETSC_ERR_ARG_CORRUPT,"DA dimension not 1, 2, or 3, it is %D\n",dim);
113:   }
114:   return(0);
115: }