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: VecStride*() operations can be useful when using DA with dof > 1
64: .keywords: distributed array, create, local, vector
66: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
67: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
68: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector(),
69: VecStrideMax(), VecStrideMin(), VecStrideNorm()
70: @*/
71: int DAGetLocalVector(DA da,Vec* g)
72: {
73: int ierr,i;
77: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
78: if (da->localin[i]) {
79: *g = da->localin[i];
80: da->localin[i] = PETSC_NULL;
81: goto alldone;
82: }
83: }
84: VecDuplicate(da->local,g);
85: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
87: alldone:
88: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
89: if (!da->localout[i]) {
90: da->localout[i] = *g;
91: break;
92: }
93: }
94: return(0);
95: }
97: /*@C
98: DARestoreLocalVector - Returns a Seq PETSc vector that
99: obtained from DAGetLocalVector(). Do not use with vector obtained via
100: DACreateLocalVector().
102: Not Collective
104: Input Parameter:
105: + da - the distributed array
106: - g - the local vector
108: Level: beginner
110: .keywords: distributed array, create, local, vector
112: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
113: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
114: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DAGetLocalVector()
115: @*/
116: int DARestoreLocalVector(DA da,Vec* g)
117: {
118: int ierr,i,j;
122: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
123: if (*g == da->localout[j]) {
124: da->localout[j] = PETSC_NULL;
125: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
126: if (!da->localin[i]) {
127: da->localin[i] = *g;
128: goto alldone;
129: }
130: }
131: }
132: }
133: VecDestroy(*g);
134: alldone:
135: return(0);
136: }
138: /*@C
139: DAGetGlobalVector - Gets a MPI PETSc vector that
140: may be used with the DAXXX routines.
142: Collective on DA
144: Input Parameter:
145: . da - the distributed array
147: Output Parameter:
148: . g - the global vector
150: Level: beginner
152: Note:
153: The output parameter, g, is a regular PETSc vector that should be returned with
154: DARestoreGlobalVector() DO NOT call VecDestroy() on it.
156: VecStride*() operations can be useful when using DA with dof > 1
158: .keywords: distributed array, create, Global, vector
160: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
161: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
162: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector()
163: VecStrideMax(), VecStrideMin(), VecStrideNorm()
165: @*/
166: int DAGetGlobalVector(DA da,Vec* g)
167: {
168: int ierr,i;
172: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
173: if (da->globalin[i]) {
174: *g = da->globalin[i];
175: da->globalin[i] = PETSC_NULL;
176: goto alldone;
177: }
178: }
179: VecDuplicate(da->global,g);
180: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
182: alldone:
183: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
184: if (!da->globalout[i]) {
185: da->globalout[i] = *g;
186: break;
187: }
188: }
189: return(0);
190: }
192: /*@C
193: DARestoreGlobalVector - Returns a Seq PETSc vector that
194: obtained from DAGetGlobalVector(). Do not use with vector obtained via
195: DACreateGlobalVector().
197: Not Collective
199: Input Parameter:
200: + da - the distributed array
201: - g - the global vector
203: Level: beginner
205: .keywords: distributed array, create, global, vector
207: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
208: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToGlobalBegin(),
209: DAGlobalToGlobalEnd(), DAGlobalToGlobal(), DACreateLocalVector(), DAGetGlobalVector()
210: @*/
211: int DARestoreGlobalVector(DA da,Vec* g)
212: {
213: int ierr,i,j;
217: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
218: if (*g == da->globalout[j]) {
219: da->globalout[j] = PETSC_NULL;
220: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
221: if (!da->globalin[i]) {
222: da->globalin[i] = *g;
223: goto alldone;
224: }
225: }
226: }
227: }
228: VecDestroy(*g);
229: alldone:
230: return(0);
231: }
233: /* ------------------------------------------------------------------- */
234: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX)
236: EXTERN_C_BEGIN
237: #include "adic/ad_utils.h"
238: EXTERN_C_END
240: /*@C
241: DAGetAdicArray - Gets an array of derivative types for a DA
242:
243: Input Parameter:
244: + da - information about my local patch
245: - ghosted - do you want arrays for the ghosted or nonghosted patch
247: Output Parameters:
248: + ptr - array data structured to be passed to ad_FormFunctionLocal()
249: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
250: - tdof - total number of degrees of freedom represented in array_start (may be null)
252: Notes: Returns the same type of object as the DAVecGetArray() except its elements are
253: derivative types instead of PetscScalars
255: Level: advanced
257: .seealso: DARestoreAdicArray()
259: @*/
260: int DAGetAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
261: {
262: int ierr,j,i,deriv_type_size,xs,ys,xm,ym,zs,zm,itdof;
263: char *iarray_start;
267: if (ghosted) {
268: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
269: if (da->adarrayghostedin[i]) {
270: *iptr = da->adarrayghostedin[i];
271: iarray_start = (char*)da->adstartghostedin[i];
272: itdof = da->ghostedtdof;
273: da->adarrayghostedin[i] = PETSC_NULL;
274: da->adstartghostedin[i] = PETSC_NULL;
275:
276: goto done;
277: }
278: }
279: xs = da->Xs;
280: ys = da->Ys;
281: zs = da->Zs;
282: xm = da->Xe-da->Xs;
283: ym = da->Ye-da->Ys;
284: zm = da->Ze-da->Zs;
285: } else {
286: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
287: if (da->adarrayin[i]) {
288: *iptr = da->adarrayin[i];
289: iarray_start = (char*)da->adstartin[i];
290: itdof = da->tdof;
291: da->adarrayin[i] = PETSC_NULL;
292: da->adstartin[i] = PETSC_NULL;
293:
294: goto done;
295: }
296: }
297: xs = da->xs;
298: ys = da->ys;
299: zs = da->zs;
300: xm = da->xe-da->xs;
301: ym = da->ye-da->ys;
302: zm = da->ze-da->zs;
303: }
304: deriv_type_size = PetscADGetDerivTypeSize();
306: switch (da->dim) {
307: case 1: {
308: void *ptr;
309: itdof = xm;
311: ierr = PetscMalloc(xm*deriv_type_size,&iarray_start);
312: ierr = PetscMemzero(iarray_start,xm*deriv_type_size);
314: ptr = (void*)(iarray_start - xs*deriv_type_size);
315: *iptr = (void*)ptr;
316: break;}
317: case 2: {
318: void **ptr;
319: itdof = xm*ym;
321: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*deriv_type_size,&iarray_start);
322: ierr = PetscMemzero(iarray_start,xm*ym*deriv_type_size);
324: ptr = (void**)(iarray_start + xm*ym*deriv_type_size - ys*sizeof(void*));
325: for(j=ys;j<ys+ym;j++) {
326: ptr[j] = iarray_start + deriv_type_size*(xm*(j-ys) - xs);
327: }
328: *iptr = (void*)ptr;
329: break;}
330: case 3: {
331: void ***ptr,**bptr;
332: itdof = xm*ym*zm;
334: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*deriv_type_size,&iarray_start);
335: ierr = PetscMemzero(iarray_start,xm*ym*zm*deriv_type_size);
337: ptr = (void***)(iarray_start + xm*ym*zm*deriv_type_size - zs*sizeof(void*));
338: bptr = (void**)(iarray_start + xm*ym*zm*deriv_type_size + zm*sizeof(void**));
339: for(i=zs;i<zs+zm;i++) {
340: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
341: }
342: for (i=zs; i<zs+zm; i++) {
343: for (j=ys; j<ys+ym; j++) {
344: ptr[i][j] = iarray_start + deriv_type_size*(xm*ym*(i-zs) + xm*(j-ys) - xs);
345: }
346: }
348: *iptr = (void*)ptr;
349: break;}
350: default:
351: SETERRQ1(1,"Dimension %d not supported",da->dim);
352: }
354: done:
355: /* add arrays to the checked out list */
356: if (ghosted) {
357: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
358: if (!da->adarrayghostedout[i]) {
359: da->adarrayghostedout[i] = *iptr ;
360: da->adstartghostedout[i] = iarray_start;
361: da->ghostedtdof = itdof;
362: break;
363: }
364: }
365: } else {
366: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
367: if (!da->adarrayout[i]) {
368: da->adarrayout[i] = *iptr ;
369: da->adstartout[i] = iarray_start;
370: da->tdof = itdof;
371: break;
372: }
373: }
374: }
375: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
376: if (tdof) *tdof = itdof;
377: if (array_start) *array_start = iarray_start;
378: return(0);
379: }
381: /*@C
382: DARestoreAdicArray - Restores an array of derivative types for a DA
383:
384: Input Parameter:
385: + da - information about my local patch
386: - ghosted - do you want arrays for the ghosted or nonghosted patch
388: Output Parameters:
389: + ptr - array data structured to be passed to ad_FormFunctionLocal()
390: . array_start - actual start of 1d array of all values that adiC can access directly
391: - tdof - total number of degrees of freedom represented in array_start
393: Level: advanced
395: .seealso: DAGetAdicArray()
397: @*/
398: int DARestoreAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
399: {
400: int i;
401: void *iarray_start = 0;
402:
405: if (ghosted) {
406: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
407: if (da->adarrayghostedout[i] == *iptr) {
408: iarray_start = da->adstartghostedout[i];
409: da->adarrayghostedout[i] = PETSC_NULL;
410: da->adstartghostedout[i] = PETSC_NULL;
411: break;
412: }
413: }
414: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
415: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
416: if (!da->adarrayghostedin[i]){
417: da->adarrayghostedin[i] = *iptr;
418: da->adstartghostedin[i] = iarray_start;
419: break;
420: }
421: }
422: } else {
423: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
424: if (da->adarrayout[i] == *iptr) {
425: iarray_start = da->adstartout[i];
426: da->adarrayout[i] = PETSC_NULL;
427: da->adstartout[i] = PETSC_NULL;
428: break;
429: }
430: }
431: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
432: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
433: if (!da->adarrayin[i]){
434: da->adarrayin[i] = *iptr;
435: da->adstartin[i] = iarray_start;
436: break;
437: }
438: }
439: }
440: return(0);
441: }
443: int ad_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
444: {
447: DAGetAdicArray(da,ghosted,iptr,0,0);
448: return(0);
449: }
451: int ad_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
452: {
455: DARestoreAdicArray(da,ghosted,iptr,0,0);
456: return(0);
457: }
459: #endif
461: /*@C
462: DAGetArray - Gets a work array for a DA
463:
464: Input Parameter:
465: + da - information about my local patch
466: - ghosted - do you want arrays for the ghosted or nonghosted patch
468: Output Parameters:
469: . ptr - array data structured
471: Level: advanced
473: .seealso: DARestoreArray(), DAGetAdicArray()
475: @*/
476: int DAGetArray(DA da,PetscTruth ghosted,void **iptr)
477: {
478: int ierr,j,i,xs,ys,xm,ym,zs,zm;
479: char *iarray_start;
483: if (ghosted) {
484: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
485: if (da->arrayghostedin[i]) {
486: *iptr = da->arrayghostedin[i];
487: iarray_start = (char*)da->startghostedin[i];
488: da->arrayghostedin[i] = PETSC_NULL;
489: da->startghostedin[i] = PETSC_NULL;
490:
491: goto done;
492: }
493: }
494: xs = da->Xs;
495: ys = da->Ys;
496: zs = da->Zs;
497: xm = da->Xe-da->Xs;
498: ym = da->Ye-da->Ys;
499: zm = da->Ze-da->Zs;
500: } else {
501: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
502: if (da->arrayin[i]) {
503: *iptr = da->arrayin[i];
504: iarray_start = (char*)da->startin[i];
505: da->arrayin[i] = PETSC_NULL;
506: da->startin[i] = PETSC_NULL;
507:
508: goto done;
509: }
510: }
511: xs = da->xs;
512: ys = da->ys;
513: zs = da->zs;
514: xm = da->xe-da->xs;
515: ym = da->ye-da->ys;
516: zm = da->ze-da->zs;
517: }
519: switch (da->dim) {
520: case 1: {
521: void *ptr;
523: ierr = PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);
524: ierr = PetscMemzero(iarray_start,xm*sizeof(PetscScalar));
526: ptr = (void*)(iarray_start - xs*sizeof(PetscScalar));
527: *iptr = (void*)ptr;
528: break;}
529: case 2: {
530: void **ptr;
532: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*sizeof(PetscScalar),&iarray_start);
533: ierr = PetscMemzero(iarray_start,xm*ym*sizeof(PetscScalar));
535: ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*));
536: for(j=ys;j<ys+ym;j++) {
537: ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs);
538: }
539: *iptr = (void*)ptr;
540: break;}
541: case 3: {
542: void ***ptr,**bptr;
544: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);
545: ierr = PetscMemzero(iarray_start,xm*ym*zm*sizeof(PetscScalar));
547: ptr = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*));
548: bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**));
549: for(i=zs;i<zs+zm;i++) {
550: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
551: }
552: for (i=zs; i<zs+zm; i++) {
553: for (j=ys; j<ys+ym; j++) {
554: ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
555: }
556: }
558: *iptr = (void*)ptr;
559: break;}
560: default:
561: SETERRQ1(1,"Dimension %d not supported",da->dim);
562: }
564: done:
565: /* add arrays to the checked out list */
566: if (ghosted) {
567: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
568: if (!da->arrayghostedout[i]) {
569: da->arrayghostedout[i] = *iptr ;
570: da->startghostedout[i] = iarray_start;
571: break;
572: }
573: }
574: } else {
575: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
576: if (!da->arrayout[i]) {
577: da->arrayout[i] = *iptr ;
578: da->startout[i] = iarray_start;
579: break;
580: }
581: }
582: }
583: return(0);
584: }
586: /*@C
587: DARestoreArray - Restores an array of derivative types for a DA
588:
589: Input Parameter:
590: + da - information about my local patch
591: . ghosted - do you want arrays for the ghosted or nonghosted patch
592: - ptr - array data structured to be passed to ad_FormFunctionLocal()
594: Level: advanced
596: .seealso: DAGetArray(), DAGetAdicArray()
598: @*/
599: int DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
600: {
601: int i;
602: void *iarray_start = 0;
603:
606: if (ghosted) {
607: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
608: if (da->arrayghostedout[i] == *iptr) {
609: iarray_start = da->startghostedout[i];
610: da->arrayghostedout[i] = PETSC_NULL;
611: da->startghostedout[i] = PETSC_NULL;
612: break;
613: }
614: }
615: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
616: if (!da->arrayghostedin[i]){
617: da->arrayghostedin[i] = *iptr;
618: da->startghostedin[i] = iarray_start;
619: break;
620: }
621: }
622: } else {
623: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
624: if (da->arrayout[i] == *iptr) {
625: iarray_start = da->startout[i];
626: da->arrayout[i] = PETSC_NULL;
627: da->startout[i] = PETSC_NULL;
628: break;
629: }
630: }
631: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
632: if (!da->arrayin[i]){
633: da->arrayin[i] = *iptr;
634: da->startin[i] = iarray_start;
635: break;
636: }
637: }
638: }
639: return(0);
640: }
642: /*@C
643: DAGetAdicMFArray - Gets an array of derivative types for a DA for matrix-free ADIC.
644:
645: Input Parameter:
646: + da - information about my local patch
647: - ghosted - do you want arrays for the ghosted or nonghosted patch?
649: Output Parameters:
650: + ptr - array data structured to be passed to ad_FormFunctionLocal()
651: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
652: - tdof - total number of degrees of freedom represented in array_start (may be null)
654: Notes:
655: This routine returns the same type of object as the DAVecGetArray(), except its
656: elements are derivative types instead of PetscScalars.
658: Level: advanced
660: .seealso: DARestoreAdicMFArray(), DAGetArray(), DAGetAdicArray()
662: @*/
663: int DAGetAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
664: {
665: int ierr,j,i,xs,ys,xm,ym,zs,zm,itdof;
666: char *iarray_start;
670: if (ghosted) {
671: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
672: if (da->admfarrayghostedin[i]) {
673: *iptr = da->admfarrayghostedin[i];
674: iarray_start = (char*)da->admfstartghostedin[i];
675: itdof = da->ghostedtdof;
676: da->admfarrayghostedin[i] = PETSC_NULL;
677: da->admfstartghostedin[i] = PETSC_NULL;
678:
679: goto done;
680: }
681: }
682: xs = da->Xs;
683: ys = da->Ys;
684: zs = da->Zs;
685: xm = da->Xe-da->Xs;
686: ym = da->Ye-da->Ys;
687: zm = da->Ze-da->Zs;
688: } else {
689: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
690: if (da->admfarrayin[i]) {
691: *iptr = da->admfarrayin[i];
692: iarray_start = (char*)da->admfstartin[i];
693: itdof = da->tdof;
694: da->admfarrayin[i] = PETSC_NULL;
695: da->admfstartin[i] = PETSC_NULL;
696:
697: goto done;
698: }
699: }
700: xs = da->xs;
701: ys = da->ys;
702: zs = da->zs;
703: xm = da->xe-da->xs;
704: ym = da->ye-da->ys;
705: zm = da->ze-da->zs;
706: }
708: switch (da->dim) {
709: case 1: {
710: void *ptr;
711: itdof = xm;
713: ierr = PetscMalloc(xm*2*sizeof(PetscScalar),&iarray_start);
714: ierr = PetscMemzero(iarray_start,xm*2*sizeof(PetscScalar));
716: ptr = (void*)(iarray_start - xs*2*sizeof(PetscScalar));
717: *iptr = (void*)ptr;
718: break;}
719: case 2: {
720: void **ptr;
721: itdof = xm*ym;
723: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*2*sizeof(PetscScalar),&iarray_start);
724: ierr = PetscMemzero(iarray_start,xm*ym*2*sizeof(PetscScalar));
726: ptr = (void**)(iarray_start + xm*ym*2*sizeof(PetscScalar) - ys*sizeof(void*));
727: for(j=ys;j<ys+ym;j++) {
728: ptr[j] = iarray_start + 2*sizeof(PetscScalar)*(xm*(j-ys) - xs);
729: }
730: *iptr = (void*)ptr;
731: break;}
732: case 3: {
733: void ***ptr,**bptr;
734: itdof = xm*ym*zm;
736: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*2*sizeof(PetscScalar),&iarray_start);
737: ierr = PetscMemzero(iarray_start,xm*ym*zm*2*sizeof(PetscScalar));
739: ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*));
740: bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**));
741: for(i=zs;i<zs+zm;i++) {
742: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
743: }
744: for (i=zs; i<zs+zm; i++) {
745: for (j=ys; j<ys+ym; j++) {
746: ptr[i][j] = iarray_start + 2*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
747: }
748: }
750: *iptr = (void*)ptr;
751: break;}
752: default:
753: SETERRQ1(1,"Dimension %d not supported",da->dim);
754: }
756: done:
757: /* add arrays to the checked out list */
758: if (ghosted) {
759: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
760: if (!da->admfarrayghostedout[i]) {
761: da->admfarrayghostedout[i] = *iptr ;
762: da->admfstartghostedout[i] = iarray_start;
763: da->ghostedtdof = itdof;
764: break;
765: }
766: }
767: } else {
768: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
769: if (!da->admfarrayout[i]) {
770: da->admfarrayout[i] = *iptr ;
771: da->admfstartout[i] = iarray_start;
772: da->tdof = itdof;
773: break;
774: }
775: }
776: }
777: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
778: if (tdof) *tdof = itdof;
779: if (array_start) *array_start = iarray_start;
780: return(0);
781: }
783: /*@C
784: DARestoreAdicMFArray - Restores an array of derivative types for a DA.
785:
786: Input Parameter:
787: + da - information about my local patch
788: - ghosted - do you want arrays for the ghosted or nonghosted patch?
790: Output Parameters:
791: + ptr - array data structure to be passed to ad_FormFunctionLocal()
792: . array_start - actual start of 1d array of all values that adiC can access directly
793: - tdof - total number of degrees of freedom represented in array_start
795: Level: advanced
797: .seealso: DAGetAdicArray()
799: @*/
800: int DARestoreAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
801: {
802: int i;
803: void *iarray_start = 0;
804:
807: if (ghosted) {
808: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
809: if (da->admfarrayghostedout[i] == *iptr) {
810: iarray_start = da->admfstartghostedout[i];
811: da->admfarrayghostedout[i] = PETSC_NULL;
812: da->admfstartghostedout[i] = PETSC_NULL;
813: break;
814: }
815: }
816: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
817: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
818: if (!da->admfarrayghostedin[i]){
819: da->admfarrayghostedin[i] = *iptr;
820: da->admfstartghostedin[i] = iarray_start;
821: break;
822: }
823: }
824: } else {
825: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
826: if (da->admfarrayout[i] == *iptr) {
827: iarray_start = da->admfstartout[i];
828: da->admfarrayout[i] = PETSC_NULL;
829: da->admfstartout[i] = PETSC_NULL;
830: break;
831: }
832: }
833: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
834: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
835: if (!da->admfarrayin[i]){
836: da->admfarrayin[i] = *iptr;
837: da->admfstartin[i] = iarray_start;
838: break;
839: }
840: }
841: }
842: return(0);
843: }
845: int admf_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
846: {
849: DAGetAdicMFArray(da,ghosted,iptr,0,0);
850: return(0);
851: }
853: int admf_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
854: {
857: DARestoreAdicMFArray(da,ghosted,iptr,0,0);
858: return(0);
859: }