Actual source code: dagetarray.c
1: #define PETSCDM_DLL
2:
3: #include petscda.h
7: /*@C
8: DAVecGetArray - Returns a multiple dimension array that shares data with
9: the underlying vector and is indexed using the global dimensions.
11: Not Collective
13: Input Parameter:
14: + da - the distributed array
15: - vec - the vector, either a vector the same size as one obtained with
16: DACreateGlobalVector() or DACreateLocalVector()
17:
18: Output Parameter:
19: . array - the array
21: Notes:
22: Call DAVecRestoreArray() once you have finished accessing the vector entries.
24: Level: intermediate
26: .keywords: distributed array, get, corners, nodes, local indices, coordinates
28: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecRestoreArray(), DAVecRestoreArrayDOF()
29: DAVecGetarrayDOF()
30: @*/
31: PetscErrorCode DAVecGetArray(DA da,Vec vec,void *array)
32: {
34: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
37: DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
38: DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
39: DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);
41: /* Handle case where user passes in global vector as opposed to local */
42: VecGetLocalSize(vec,&N);
43: if (N == xm*ym*zm*dof) {
44: gxm = xm;
45: gym = ym;
46: gzm = zm;
47: gxs = xs;
48: gys = ys;
49: gzs = zs;
50: } else if (N != gxm*gym*gzm*dof) {
51: 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);
52: }
54: if (dim == 1) {
55: VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);
56: } else if (dim == 2) {
57: VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
58: } else if (dim == 3) {
59: VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
60: } else {
61: SETERRQ1(PETSC_ERR_ARG_CORRUPT,"DA dimension not 1, 2, or 3, it is %D\n",dim);
62: }
64: return(0);
65: }
69: /*@
70: DAVecRestoreArray - Restores a multiple dimension array obtained with DAVecGetArray()
72: Not Collective
74: Input Parameter:
75: + da - the distributed array
76: . vec - the vector, either a vector the same size as one obtained with
77: DACreateGlobalVector() or DACreateLocalVector()
78: - array - the array
80: Level: intermediate
82: .keywords: distributed array, get, corners, nodes, local indices, coordinates
84: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecGetArray()
85: @*/
86: PetscErrorCode DAVecRestoreArray(DA da,Vec vec,void *array)
87: {
89: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
92: DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
93: DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
94: DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);
96: /* Handle case where user passes in global vector as opposed to local */
97: VecGetLocalSize(vec,&N);
98: if (N == xm*ym*zm*dof) {
99: gxm = xm;
100: gym = ym;
101: gzm = zm;
102: gxs = xs;
103: gys = ys;
104: gzs = zs;
105: } else if (N != gxm*gym*gzm*dof) {
106: 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);
107: }
109: if (dim == 1) {
110: VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);
111: } else if (dim == 2) {
112: VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
113: } else if (dim == 3) {
114: VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
115: } else {
116: SETERRQ1(PETSC_ERR_ARG_CORRUPT,"DA dimension not 1, 2, or 3, it is %D\n",dim);
117: }
118: return(0);
119: }
123: /*@C
124: DAVecGetArrayDOF - Returns a multiple dimension array that shares data with
125: the underlying vector and is indexed using the global dimensions.
127: Not Collective
129: Input Parameter:
130: + da - the distributed array
131: - vec - the vector, either a vector the same size as one obtained with
132: DACreateGlobalVector() or DACreateLocalVector()
133:
134: Output Parameter:
135: . array - the array
137: Notes:
138: Call DAVecRestoreArray() once you have finished accessing the vector entries.
140: Level: intermediate
142: .keywords: distributed array, get, corners, nodes, local indices, coordinates
144: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecRestoreArray(), DAVecGetArray(), DAVecRestoreArrayDOF()
145: @*/
146: PetscErrorCode DAVecGetArrayDOF(DA da,Vec vec,void *array)
147: {
149: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
152: DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
153: DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
154: DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);
156: /* Handle case where user passes in global vector as opposed to local */
157: VecGetLocalSize(vec,&N);
158: if (N == xm*ym*zm*dof) {
159: gxm = xm;
160: gym = ym;
161: gzm = zm;
162: gxs = xs;
163: gys = ys;
164: gzs = zs;
165: } else if (N != gxm*gym*gzm*dof) {
166: 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);
167: }
169: if (dim == 1) {
170: VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar ***)array);
171: } else if (dim == 2) {
172: VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
173: } else if (dim == 3) {
174: VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
175: } else {
176: SETERRQ1(PETSC_ERR_ARG_CORRUPT,"DA dimension not 1, 2, or 3, it is %D\n",dim);
177: }
178: return(0);
179: }
183: /*@
184: DAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DAVecGetArrayDOF()
186: Not Collective
188: Input Parameter:
189: + da - the distributed array
190: . vec - the vector, either a vector the same size as one obtained with
191: DACreateGlobalVector() or DACreateLocalVector()
192: - array - the array
194: Level: intermediate
196: .keywords: distributed array, get, corners, nodes, local indices, coordinates
198: .seealso: DAGetGhostCorners(), DAGetCorners(), VecGetArray(), VecRestoreArray(), DAVecGetArray(), DAVecGetArrayDOF(), DAVecRestoreArrayDOF()
199: @*/
200: PetscErrorCode DAVecRestoreArrayDOF(DA da,Vec vec,void *array)
201: {
203: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
206: DAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
207: DAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
208: DAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0);
210: /* Handle case where user passes in global vector as opposed to local */
211: VecGetLocalSize(vec,&N);
212: if (N == xm*ym*zm*dof) {
213: gxm = xm;
214: gym = ym;
215: gzm = zm;
216: gxs = xs;
217: gys = ys;
218: gzs = zs;
219: }
221: if (dim == 1) {
222: VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
223: } else if (dim == 2) {
224: VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
225: } else if (dim == 3) {
226: VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
227: } else {
228: SETERRQ1(PETSC_ERR_ARG_CORRUPT,"DA dimension not 1, 2, or 3, it is %D\n",dim);
229: }
230: return(0);
231: }