Actual source code: dalocal.c
1: /*$Id: dalocal.c,v 1.35 2001/08/07 03:04:39 balay Exp $*/
2:
3: /*
4: Code for manipulating distributed regular arrays in parallel.
5: */
7: #include src/dm/da/daimpl.h
9: /*@C
10: DACreateLocalVector - Creates a Seq PETSc vector that
11: may be used with the DAXXX routines.
13: Not Collective
15: Input Parameter:
16: . da - the distributed array
18: Output Parameter:
19: . g - the local vector
21: Level: beginner
23: Note:
24: The output parameter, g, is a regular PETSc vector that should be destroyed
25: with a call to VecDestroy() when usage is finished.
27: .keywords: distributed array, create, local, vector
29: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
30: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
31: DAGlobalToLocalEnd(), DALocalToGlobal(), DAGetLocalVector(), DARestoreLocalVector()
32: @*/
33: int DACreateLocalVector(DA da,Vec* g)
34: {
39: VecDuplicate(da->local,g);
40: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
41: return(0);
42: }
44: /*@C
45: DAGetLocalVector - Gets a Seq PETSc vector that
46: may be used with the DAXXX routines.
48: Not Collective
50: Input Parameter:
51: . da - the distributed array
53: Output Parameter:
54: . g - the local vector
56: Level: beginner
58: Note:
59: The output parameter, g, is a regular PETSc vector that should be returned with
60: DARestoreLocalVector() DO NOT call VecDestroy() on it.
62: .keywords: distributed array, create, local, vector
64: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
65: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
66: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector()
67: @*/
68: int DAGetLocalVector(DA da,Vec* g)
69: {
70: int ierr,i;
74: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
75: if (da->localin[i]) {
76: *g = da->localin[i];
77: da->localin[i] = PETSC_NULL;
78: goto alldone;
79: }
80: }
81: VecDuplicate(da->local,g);
82: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
84: alldone:
85: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
86: if (!da->localout[i]) {
87: da->localout[i] = *g;
88: break;
89: }
90: }
91: return(0);
92: }
94: /*@C
95: DARestoreLocalVector - Returns a Seq PETSc vector that
96: obtained from DAGetLocalVector(). Do not use with vector obtained via
97: DACreateLocalVector().
99: Not Collective
101: Input Parameter:
102: + da - the distributed array
103: - g - the local vector
105: Level: beginner
107: .keywords: distributed array, create, local, vector
109: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
110: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
111: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DAGetLocalVector()
112: @*/
113: int DARestoreLocalVector(DA da,Vec* g)
114: {
115: int ierr,i,j;
119: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
120: if (*g == da->localout[j]) {
121: da->localout[j] = PETSC_NULL;
122: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
123: if (!da->localin[i]) {
124: da->localin[i] = *g;
125: goto alldone;
126: }
127: }
128: }
129: }
130: VecDestroy(*g);
131: alldone:
132: return(0);
133: }
135: /*@C
136: DAGetGlobalVector - Gets a MPI PETSc vector that
137: may be used with the DAXXX routines.
139: Collective on DA
141: Input Parameter:
142: . da - the distributed array
144: Output Parameter:
145: . g - the global vector
147: Level: beginner
149: Note:
150: The output parameter, g, is a regular PETSc vector that should be returned with
151: DARestoreGlobalVector() DO NOT call VecDestroy() on it.
153: .keywords: distributed array, create, Global, vector
155: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
156: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
157: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector()
158: @*/
159: int DAGetGlobalVector(DA da,Vec* g)
160: {
161: int ierr,i;
165: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
166: if (da->globalin[i]) {
167: *g = da->globalin[i];
168: da->globalin[i] = PETSC_NULL;
169: goto alldone;
170: }
171: }
172: VecDuplicate(da->global,g);
173: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
175: alldone:
176: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
177: if (!da->globalout[i]) {
178: da->globalout[i] = *g;
179: break;
180: }
181: }
182: return(0);
183: }
185: /*@C
186: DARestoreGlobalVector - Returns a Seq PETSc vector that
187: obtained from DAGetGlobalVector(). Do not use with vector obtained via
188: DACreateGlobalVector().
190: Not Collective
192: Input Parameter:
193: + da - the distributed array
194: - g - the global vector
196: Level: beginner
198: .keywords: distributed array, create, global, vector
200: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
201: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToGlobalBegin(),
202: DAGlobalToGlobalEnd(), DAGlobalToGlobal(), DACreateLocalVector(), DAGetGlobalVector()
203: @*/
204: int DARestoreGlobalVector(DA da,Vec* g)
205: {
206: int ierr,i,j;
210: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
211: if (*g == da->globalout[j]) {
212: da->globalout[j] = PETSC_NULL;
213: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
214: if (!da->globalin[i]) {
215: da->globalin[i] = *g;
216: goto alldone;
217: }
218: }
219: }
220: }
221: VecDestroy(*g);
222: alldone:
223: return(0);
224: }
226: /* ------------------------------------------------------------------- */
227: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX)
229: EXTERN_C_BEGIN
230: #include "adic/ad_utils.h"
231: EXTERN_C_END
233: /*@C
234: DAGetAdicArray - Gets an array of derivative types for a DA
235:
236: Input Parameter:
237: + da - information about my local patch
238: - ghosted - do you want arrays for the ghosted or nonghosted patch
240: Output Parameters:
241: + ptr - array data structured to be passed to ad_FormFunctionLocal()
242: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
243: - tdof - total number of degrees of freedom represented in array_start (may be null)
245: Notes: Returns the same type of object as the DAVecGetArray() except its elements are
246: derivative types instead of PetscScalars
248: Level: advanced
250: .seealso: DARestoreAdicArray()
252: @*/
253: int DAGetAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
254: {
255: int ierr,j,i,deriv_type_size,xs,ys,xm,ym,zs,zm,itdof;
256: char *iarray_start;
260: if (ghosted) {
261: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
262: if (da->adarrayghostedin[i]) {
263: *iptr = da->adarrayghostedin[i];
264: iarray_start = (char*)da->adstartghostedin[i];
265: itdof = da->ghostedtdof;
266: da->adarrayghostedin[i] = PETSC_NULL;
267: da->adstartghostedin[i] = PETSC_NULL;
268:
269: goto done;
270: }
271: }
272: xs = da->Xs;
273: ys = da->Ys;
274: zs = da->Zs;
275: xm = da->Xe-da->Xs;
276: ym = da->Ye-da->Ys;
277: zm = da->Ze-da->Zs;
278: } else {
279: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
280: if (da->adarrayin[i]) {
281: *iptr = da->adarrayin[i];
282: iarray_start = (char*)da->adstartin[i];
283: itdof = da->tdof;
284: da->adarrayin[i] = PETSC_NULL;
285: da->adstartin[i] = PETSC_NULL;
286:
287: goto done;
288: }
289: }
290: xs = da->xs;
291: ys = da->ys;
292: zs = da->zs;
293: xm = da->xe-da->xs;
294: ym = da->ye-da->ys;
295: zm = da->ze-da->zs;
296: }
297: deriv_type_size = PetscADGetDerivTypeSize();
299: switch (da->dim) {
300: case 1: {
301: void *ptr;
302: itdof = xm;
304: ierr = PetscMalloc(xm*deriv_type_size,&iarray_start);
305: ierr = PetscMemzero(iarray_start,xm*deriv_type_size);
307: ptr = (void*)(iarray_start - xs*deriv_type_size);
308: *iptr = (void*)ptr;
309: break;}
310: case 2: {
311: void **ptr;
312: itdof = xm*ym;
314: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*deriv_type_size,&iarray_start);
315: ierr = PetscMemzero(iarray_start,xm*ym*deriv_type_size);
317: ptr = (void**)(iarray_start + xm*ym*deriv_type_size - ys*sizeof(void*));
318: for(j=ys;j<ys+ym;j++) {
319: ptr[j] = iarray_start + deriv_type_size*(xm*(j-ys) - xs);
320: }
321: *iptr = (void*)ptr;
322: break;}
323: case 3: {
324: void ***ptr,**bptr;
325: itdof = xm*ym*zm;
327: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*deriv_type_size,&iarray_start);
328: ierr = PetscMemzero(iarray_start,xm*ym*zm*deriv_type_size);
330: ptr = (void***)(iarray_start + xm*ym*zm*deriv_type_size - zs*sizeof(void*));
331: bptr = (void**)(iarray_start + xm*ym*zm*deriv_type_size + zm*sizeof(void**));
332: for(i=zs;i<zs+zm;i++) {
333: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
334: }
335: for (i=zs; i<zs+zm; i++) {
336: for (j=ys; j<ys+ym; j++) {
337: ptr[i][j] = iarray_start + deriv_type_size*(xm*ym*(i-zs) + xm*(j-ys) - xs);
338: }
339: }
341: *iptr = (void*)ptr;
342: break;}
343: default:
344: SETERRQ1(1,"Dimension %d not supported",da->dim);
345: }
347: done:
348: /* add arrays to the checked out list */
349: if (ghosted) {
350: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
351: if (!da->adarrayghostedout[i]) {
352: da->adarrayghostedout[i] = *iptr ;
353: da->adstartghostedout[i] = iarray_start;
354: da->ghostedtdof = itdof;
355: break;
356: }
357: }
358: } else {
359: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
360: if (!da->adarrayout[i]) {
361: da->adarrayout[i] = *iptr ;
362: da->adstartout[i] = iarray_start;
363: da->tdof = itdof;
364: break;
365: }
366: }
367: }
368: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
369: if (tdof) *tdof = itdof;
370: if (array_start) *array_start = iarray_start;
371: return(0);
372: }
374: /*@C
375: DARestoreAdicArray - Restores an array of derivative types for a DA
376:
377: Input Parameter:
378: + da - information about my local patch
379: - ghosted - do you want arrays for the ghosted or nonghosted patch
381: Output Parameters:
382: + ptr - array data structured to be passed to ad_FormFunctionLocal()
383: . array_start - actual start of 1d array of all values that adiC can access directly
384: - tdof - total number of degrees of freedom represented in array_start
386: Level: advanced
388: .seealso: DAGetAdicArray()
390: @*/
391: int DARestoreAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
392: {
393: int i;
394: void *iarray_start = 0;
395:
398: if (ghosted) {
399: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
400: if (da->adarrayghostedout[i] == *iptr) {
401: iarray_start = da->adstartghostedout[i];
402: da->adarrayghostedout[i] = PETSC_NULL;
403: da->adstartghostedout[i] = PETSC_NULL;
404: break;
405: }
406: }
407: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
408: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
409: if (!da->adarrayghostedin[i]){
410: da->adarrayghostedin[i] = *iptr;
411: da->adstartghostedin[i] = iarray_start;
412: break;
413: }
414: }
415: } else {
416: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
417: if (da->adarrayout[i] == *iptr) {
418: iarray_start = da->adstartout[i];
419: da->adarrayout[i] = PETSC_NULL;
420: da->adstartout[i] = PETSC_NULL;
421: break;
422: }
423: }
424: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
425: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
426: if (!da->adarrayin[i]){
427: da->adarrayin[i] = *iptr;
428: da->adstartin[i] = iarray_start;
429: break;
430: }
431: }
432: }
433: return(0);
434: }
436: int ad_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
437: {
440: DAGetAdicArray(da,ghosted,iptr,0,0);
441: return(0);
442: }
444: int ad_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
445: {
448: DARestoreAdicArray(da,ghosted,iptr,0,0);
449: return(0);
450: }
452: #endif
454: /*@C
455: DAGetArray - Gets a work array for a DA
456:
457: Input Parameter:
458: + da - information about my local patch
459: - ghosted - do you want arrays for the ghosted or nonghosted patch
461: Output Parameters:
462: . ptr - array data structured
464: Level: advanced
466: .seealso: DARestoreArray(), DAGetAdicArray()
468: @*/
469: int DAGetArray(DA da,PetscTruth ghosted,void **iptr)
470: {
471: int ierr,j,i,xs,ys,xm,ym,zs,zm;
472: char *iarray_start;
476: if (ghosted) {
477: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
478: if (da->arrayghostedin[i]) {
479: *iptr = da->arrayghostedin[i];
480: iarray_start = (char*)da->startghostedin[i];
481: da->arrayghostedin[i] = PETSC_NULL;
482: da->startghostedin[i] = PETSC_NULL;
483:
484: goto done;
485: }
486: }
487: xs = da->Xs;
488: ys = da->Ys;
489: zs = da->Zs;
490: xm = da->Xe-da->Xs;
491: ym = da->Ye-da->Ys;
492: zm = da->Ze-da->Zs;
493: } else {
494: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
495: if (da->arrayin[i]) {
496: *iptr = da->arrayin[i];
497: iarray_start = (char*)da->startin[i];
498: da->arrayin[i] = PETSC_NULL;
499: da->startin[i] = PETSC_NULL;
500:
501: goto done;
502: }
503: }
504: xs = da->xs;
505: ys = da->ys;
506: zs = da->zs;
507: xm = da->xe-da->xs;
508: ym = da->ye-da->ys;
509: zm = da->ze-da->zs;
510: }
512: switch (da->dim) {
513: case 1: {
514: void *ptr;
516: ierr = PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);
517: ierr = PetscMemzero(iarray_start,xm*sizeof(PetscScalar));
519: ptr = (void*)(iarray_start - xs*sizeof(PetscScalar));
520: *iptr = (void*)ptr;
521: break;}
522: case 2: {
523: void **ptr;
525: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*sizeof(PetscScalar),&iarray_start);
526: ierr = PetscMemzero(iarray_start,xm*ym*sizeof(PetscScalar));
528: ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*));
529: for(j=ys;j<ys+ym;j++) {
530: ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs);
531: }
532: *iptr = (void*)ptr;
533: break;}
534: case 3: {
535: void ***ptr,**bptr;
537: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);
538: ierr = PetscMemzero(iarray_start,xm*ym*zm*sizeof(PetscScalar));
540: ptr = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*));
541: bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**));
542: for(i=zs;i<zs+zm;i++) {
543: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
544: }
545: for (i=zs; i<zs+zm; i++) {
546: for (j=ys; j<ys+ym; j++) {
547: ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
548: }
549: }
551: *iptr = (void*)ptr;
552: break;}
553: default:
554: SETERRQ1(1,"Dimension %d not supported",da->dim);
555: }
557: done:
558: /* add arrays to the checked out list */
559: if (ghosted) {
560: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
561: if (!da->arrayghostedout[i]) {
562: da->arrayghostedout[i] = *iptr ;
563: da->startghostedout[i] = iarray_start;
564: break;
565: }
566: }
567: } else {
568: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
569: if (!da->arrayout[i]) {
570: da->arrayout[i] = *iptr ;
571: da->startout[i] = iarray_start;
572: break;
573: }
574: }
575: }
576: return(0);
577: }
579: /*@C
580: DARestoreArray - Restores an array of derivative types for a DA
581:
582: Input Parameter:
583: + da - information about my local patch
584: . ghosted - do you want arrays for the ghosted or nonghosted patch
585: - ptr - array data structured to be passed to ad_FormFunctionLocal()
587: Level: advanced
589: .seealso: DAGetArray(), DAGetAdicArray()
591: @*/
592: int DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
593: {
594: int i;
595: void *iarray_start = 0;
596:
599: if (ghosted) {
600: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
601: if (da->arrayghostedout[i] == *iptr) {
602: iarray_start = da->startghostedout[i];
603: da->arrayghostedout[i] = PETSC_NULL;
604: da->startghostedout[i] = PETSC_NULL;
605: break;
606: }
607: }
608: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
609: if (!da->arrayghostedin[i]){
610: da->arrayghostedin[i] = *iptr;
611: da->startghostedin[i] = iarray_start;
612: break;
613: }
614: }
615: } else {
616: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
617: if (da->arrayout[i] == *iptr) {
618: iarray_start = da->startout[i];
619: da->arrayout[i] = PETSC_NULL;
620: da->startout[i] = PETSC_NULL;
621: break;
622: }
623: }
624: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
625: if (!da->arrayin[i]){
626: da->arrayin[i] = *iptr;
627: da->startin[i] = iarray_start;
628: break;
629: }
630: }
631: }
632: return(0);
633: }
635: /*@C
636: DAGetAdicMFArray - Gets an array of derivative types for a DA for matrix-free ADIC.
637:
638: Input Parameter:
639: + da - information about my local patch
640: - ghosted - do you want arrays for the ghosted or nonghosted patch?
642: Output Parameters:
643: + ptr - array data structured to be passed to ad_FormFunctionLocal()
644: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
645: - tdof - total number of degrees of freedom represented in array_start (may be null)
647: Notes:
648: This routine returns the same type of object as the DAVecGetArray(), except its
649: elements are derivative types instead of PetscScalars.
651: Level: advanced
653: .seealso: DARestoreAdicMFArray(), DAGetArray(), DAGetAdicArray()
655: @*/
656: int DAGetAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
657: {
658: int ierr,j,i,xs,ys,xm,ym,zs,zm,itdof;
659: char *iarray_start;
663: if (ghosted) {
664: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
665: if (da->admfarrayghostedin[i]) {
666: *iptr = da->admfarrayghostedin[i];
667: iarray_start = (char*)da->admfstartghostedin[i];
668: itdof = da->ghostedtdof;
669: da->admfarrayghostedin[i] = PETSC_NULL;
670: da->admfstartghostedin[i] = PETSC_NULL;
671:
672: goto done;
673: }
674: }
675: xs = da->Xs;
676: ys = da->Ys;
677: zs = da->Zs;
678: xm = da->Xe-da->Xs;
679: ym = da->Ye-da->Ys;
680: zm = da->Ze-da->Zs;
681: } else {
682: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
683: if (da->admfarrayin[i]) {
684: *iptr = da->admfarrayin[i];
685: iarray_start = (char*)da->admfstartin[i];
686: itdof = da->tdof;
687: da->admfarrayin[i] = PETSC_NULL;
688: da->admfstartin[i] = PETSC_NULL;
689:
690: goto done;
691: }
692: }
693: xs = da->xs;
694: ys = da->ys;
695: zs = da->zs;
696: xm = da->xe-da->xs;
697: ym = da->ye-da->ys;
698: zm = da->ze-da->zs;
699: }
701: switch (da->dim) {
702: case 1: {
703: void *ptr;
704: itdof = xm;
706: ierr = PetscMalloc(xm*2*sizeof(PetscScalar),&iarray_start);
707: ierr = PetscMemzero(iarray_start,xm*2*sizeof(PetscScalar));
709: ptr = (void*)(iarray_start - xs*2*sizeof(PetscScalar));
710: *iptr = (void*)ptr;
711: break;}
712: case 2: {
713: void **ptr;
714: itdof = xm*ym;
716: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*2*sizeof(PetscScalar),&iarray_start);
717: ierr = PetscMemzero(iarray_start,xm*ym*2*sizeof(PetscScalar));
719: ptr = (void**)(iarray_start + xm*ym*2*sizeof(PetscScalar) - ys*sizeof(void*));
720: for(j=ys;j<ys+ym;j++) {
721: ptr[j] = iarray_start + 2*sizeof(PetscScalar)*(xm*(j-ys) - xs);
722: }
723: *iptr = (void*)ptr;
724: break;}
725: case 3: {
726: void ***ptr,**bptr;
727: itdof = xm*ym*zm;
729: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*2*sizeof(PetscScalar),&iarray_start);
730: ierr = PetscMemzero(iarray_start,xm*ym*zm*2*sizeof(PetscScalar));
732: ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*));
733: bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**));
734: for(i=zs;i<zs+zm;i++) {
735: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
736: }
737: for (i=zs; i<zs+zm; i++) {
738: for (j=ys; j<ys+ym; j++) {
739: ptr[i][j] = iarray_start + 2*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
740: }
741: }
743: *iptr = (void*)ptr;
744: break;}
745: default:
746: SETERRQ1(1,"Dimension %d not supported",da->dim);
747: }
749: done:
750: /* add arrays to the checked out list */
751: if (ghosted) {
752: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
753: if (!da->admfarrayghostedout[i]) {
754: da->admfarrayghostedout[i] = *iptr ;
755: da->admfstartghostedout[i] = iarray_start;
756: da->ghostedtdof = itdof;
757: break;
758: }
759: }
760: } else {
761: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
762: if (!da->admfarrayout[i]) {
763: da->admfarrayout[i] = *iptr ;
764: da->admfstartout[i] = iarray_start;
765: da->tdof = itdof;
766: break;
767: }
768: }
769: }
770: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
771: if (tdof) *tdof = itdof;
772: if (array_start) *array_start = iarray_start;
773: return(0);
774: }
776: /*@C
777: DARestoreAdicMFArray - Restores an array of derivative types for a DA.
778:
779: Input Parameter:
780: + da - information about my local patch
781: - ghosted - do you want arrays for the ghosted or nonghosted patch?
783: Output Parameters:
784: + ptr - array data structure to be passed to ad_FormFunctionLocal()
785: . array_start - actual start of 1d array of all values that adiC can access directly
786: - tdof - total number of degrees of freedom represented in array_start
788: Level: advanced
790: .seealso: DAGetAdicArray()
792: @*/
793: int DARestoreAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
794: {
795: int i;
796: void *iarray_start = 0;
797:
800: if (ghosted) {
801: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
802: if (da->admfarrayghostedout[i] == *iptr) {
803: iarray_start = da->admfstartghostedout[i];
804: da->admfarrayghostedout[i] = PETSC_NULL;
805: da->admfstartghostedout[i] = PETSC_NULL;
806: break;
807: }
808: }
809: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
810: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
811: if (!da->admfarrayghostedin[i]){
812: da->admfarrayghostedin[i] = *iptr;
813: da->admfstartghostedin[i] = iarray_start;
814: break;
815: }
816: }
817: } else {
818: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
819: if (da->admfarrayout[i] == *iptr) {
820: iarray_start = da->admfstartout[i];
821: da->admfarrayout[i] = PETSC_NULL;
822: da->admfstartout[i] = PETSC_NULL;
823: break;
824: }
825: }
826: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
827: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
828: if (!da->admfarrayin[i]){
829: da->admfarrayin[i] = *iptr;
830: da->admfstartin[i] = iarray_start;
831: break;
832: }
833: }
834: }
835: return(0);
836: }
838: int admf_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
839: {
842: DAGetAdicMFArray(da,ghosted,iptr,0,0);
843: return(0);
844: }
846: int admf_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
847: {
850: DARestoreAdicMFArray(da,ghosted,iptr,0,0);
851: return(0);
852: }