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: /*
10: This allows the DA vectors to properly tell Matlab their dimensions
11: */
12: #if defined(PETSC_HAVE_MATLAB_ENGINE) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
13: #include "engine.h" /* Matlab include file */
14: #include "mex.h" /* Matlab include file */
15: EXTERN_C_BEGIN
16: #undef __FUNCT__
18: int VecMatlabEnginePut_DA2d(PetscObject obj,void *mengine)
19: {
20: int ierr,n,m;
21: Vec vec = (Vec)obj;
22: PetscScalar *array;
23: mxArray *mat;
24: DA da;
27: PetscObjectQuery((PetscObject)vec,"DA",(PetscObject*)&da);
28: if (!da) SETERRQ(1,"Vector not associated with a DA");
29: DAGetGhostCorners(da,0,0,0,&m,&n,0);
31: VecGetArray(vec,&array);
32: #if !defined(PETSC_USE_COMPLEX)
33: mat = mxCreateDoubleMatrix(m,n,mxREAL);
34: #else
35: mat = mxCreateDoubleMatrix(m,n,mxCOMPLEX);
36: #endif
37: PetscMemcpy(mxGetPr(mat),array,n*m*sizeof(PetscScalar));
38: PetscObjectName(obj);
39: mxSetName(mat,obj->name);
40: engPutArray((Engine *)mengine,mat);
41:
42: VecRestoreArray(vec,&array);
43: return(0);
44: }
45: EXTERN_C_END
46: #endif
49: #undef __FUNCT__
51: /*@C
52: DACreateLocalVector - Creates a Seq PETSc vector that
53: may be used with the DAXXX routines.
55: Not Collective
57: Input Parameter:
58: . da - the distributed array
60: Output Parameter:
61: . g - the local vector
63: Level: beginner
65: Note:
66: The output parameter, g, is a regular PETSc vector that should be destroyed
67: with a call to VecDestroy() when usage is finished.
69: .keywords: distributed array, create, local, vector
71: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
72: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
73: DAGlobalToLocalEnd(), DALocalToGlobal(), DAGetLocalVector(), DARestoreLocalVector()
74: @*/
75: int DACreateLocalVector(DA da,Vec* g)
76: {
81: VecCreateSeq(PETSC_COMM_SELF,da->nlocal,g);
82: VecSetBlockSize(*g,da->w);
83: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
84: #if defined(PETSC_HAVE_MATLAB_ENGINE) && !defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_SINGLE)
85: if (da->w == 1 && da->dim == 2) {
86: PetscObjectComposeFunctionDynamic((PetscObject)*g,"PetscMatlabEnginePut_C","VecMatlabEnginePut_DA2d",VecMatlabEnginePut_DA2d);
87: }
88: #endif
89: return(0);
90: }
92: #undef __FUNCT__
94: /*@C
95: DAGetLocalVector - Gets a Seq PETSc vector that
96: may be used with the DAXXX routines.
98: Not Collective
100: Input Parameter:
101: . da - the distributed array
103: Output Parameter:
104: . g - the local vector
106: Level: beginner
108: Note:
109: The output parameter, g, is a regular PETSc vector that should be returned with
110: DARestoreLocalVector() DO NOT call VecDestroy() on it.
112: VecStride*() operations can be useful when using DA with dof > 1
114: .keywords: distributed array, create, local, vector
116: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
117: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
118: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector(),
119: VecStrideMax(), VecStrideMin(), VecStrideNorm()
120: @*/
121: int DAGetLocalVector(DA da,Vec* g)
122: {
123: int ierr,i;
127: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
128: if (da->localin[i]) {
129: *g = da->localin[i];
130: da->localin[i] = PETSC_NULL;
131: goto alldone;
132: }
133: }
134: DACreateLocalVector(da,g);
136: alldone:
137: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
138: if (!da->localout[i]) {
139: da->localout[i] = *g;
140: break;
141: }
142: }
143: return(0);
144: }
146: #undef __FUNCT__
148: /*@C
149: DARestoreLocalVector - Returns a Seq PETSc vector that
150: obtained from DAGetLocalVector(). Do not use with vector obtained via
151: DACreateLocalVector().
153: Not Collective
155: Input Parameter:
156: + da - the distributed array
157: - g - the local vector
159: Level: beginner
161: .keywords: distributed array, create, local, vector
163: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
164: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
165: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DAGetLocalVector()
166: @*/
167: int DARestoreLocalVector(DA da,Vec* g)
168: {
169: int ierr,i,j;
173: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
174: if (*g == da->localout[j]) {
175: da->localout[j] = PETSC_NULL;
176: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
177: if (!da->localin[i]) {
178: da->localin[i] = *g;
179: goto alldone;
180: }
181: }
182: }
183: }
184: VecDestroy(*g);
185: alldone:
186: return(0);
187: }
189: #undef __FUNCT__
191: /*@C
192: DAGetGlobalVector - Gets a MPI PETSc vector that
193: may be used with the DAXXX routines.
195: Collective on DA
197: Input Parameter:
198: . da - the distributed array
200: Output Parameter:
201: . g - the global vector
203: Level: beginner
205: Note:
206: The output parameter, g, is a regular PETSc vector that should be returned with
207: DARestoreGlobalVector() DO NOT call VecDestroy() on it.
209: VecStride*() operations can be useful when using DA with dof > 1
211: .keywords: distributed array, create, Global, vector
213: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
214: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
215: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateLocalVector(), DARestoreLocalVector()
216: VecStrideMax(), VecStrideMin(), VecStrideNorm()
218: @*/
219: int DAGetGlobalVector(DA da,Vec* g)
220: {
221: int ierr,i;
225: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
226: if (da->globalin[i]) {
227: *g = da->globalin[i];
228: da->globalin[i] = PETSC_NULL;
229: goto alldone;
230: }
231: }
232: DACreateGlobalVector(da,g);
234: alldone:
235: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
236: if (!da->globalout[i]) {
237: da->globalout[i] = *g;
238: break;
239: }
240: }
241: return(0);
242: }
244: #undef __FUNCT__
246: /*@C
247: DARestoreGlobalVector - Returns a Seq PETSc vector that
248: obtained from DAGetGlobalVector(). Do not use with vector obtained via
249: DACreateGlobalVector().
251: Not Collective
253: Input Parameter:
254: + da - the distributed array
255: - g - the global vector
257: Level: beginner
259: .keywords: distributed array, create, global, vector
261: .seealso: DACreateGlobalVector(), VecDuplicate(), VecDuplicateVecs(),
262: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToGlobalBegin(),
263: DAGlobalToGlobalEnd(), DAGlobalToGlobal(), DACreateLocalVector(), DAGetGlobalVector()
264: @*/
265: int DARestoreGlobalVector(DA da,Vec* g)
266: {
267: int ierr,i,j;
271: for (j=0; j<DA_MAX_WORK_VECTORS; j++) {
272: if (*g == da->globalout[j]) {
273: da->globalout[j] = PETSC_NULL;
274: for (i=0; i<DA_MAX_WORK_VECTORS; i++) {
275: if (!da->globalin[i]) {
276: da->globalin[i] = *g;
277: goto alldone;
278: }
279: }
280: }
281: }
282: VecDestroy(*g);
283: alldone:
284: return(0);
285: }
287: /* ------------------------------------------------------------------- */
288: #if defined(PETSC_HAVE_ADIC) && !defined(PETSC_USE_COMPLEX)
290: EXTERN_C_BEGIN
291: #include "adic/ad_utils.h"
292: EXTERN_C_END
294: #undef __FUNCT__
296: /*@C
297: DAGetAdicArray - Gets an array of derivative types for a DA
298:
299: Input Parameter:
300: + da - information about my local patch
301: - ghosted - do you want arrays for the ghosted or nonghosted patch
303: Output Parameters:
304: + ptr - array data structured to be passed to ad_FormFunctionLocal()
305: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
306: - tdof - total number of degrees of freedom represented in array_start (may be null)
308: Notes: Returns the same type of object as the DAVecGetArray() except its elements are
309: derivative types instead of PetscScalars
311: Level: advanced
313: .seealso: DARestoreAdicArray()
315: @*/
316: int DAGetAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
317: {
318: int ierr,j,i,deriv_type_size,xs,ys,xm,ym,zs,zm,itdof;
319: char *iarray_start;
323: if (ghosted) {
324: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
325: if (da->adarrayghostedin[i]) {
326: *iptr = da->adarrayghostedin[i];
327: iarray_start = (char*)da->adstartghostedin[i];
328: itdof = da->ghostedtdof;
329: da->adarrayghostedin[i] = PETSC_NULL;
330: da->adstartghostedin[i] = PETSC_NULL;
331:
332: goto done;
333: }
334: }
335: xs = da->Xs;
336: ys = da->Ys;
337: zs = da->Zs;
338: xm = da->Xe-da->Xs;
339: ym = da->Ye-da->Ys;
340: zm = da->Ze-da->Zs;
341: } else {
342: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
343: if (da->adarrayin[i]) {
344: *iptr = da->adarrayin[i];
345: iarray_start = (char*)da->adstartin[i];
346: itdof = da->tdof;
347: da->adarrayin[i] = PETSC_NULL;
348: da->adstartin[i] = PETSC_NULL;
349:
350: goto done;
351: }
352: }
353: xs = da->xs;
354: ys = da->ys;
355: zs = da->zs;
356: xm = da->xe-da->xs;
357: ym = da->ye-da->ys;
358: zm = da->ze-da->zs;
359: }
360: deriv_type_size = PetscADGetDerivTypeSize();
362: switch (da->dim) {
363: case 1: {
364: void *ptr;
365: itdof = xm;
367: ierr = PetscMalloc(xm*deriv_type_size,&iarray_start);
368: ierr = PetscMemzero(iarray_start,xm*deriv_type_size);
370: ptr = (void*)(iarray_start - xs*deriv_type_size);
371: *iptr = (void*)ptr;
372: break;}
373: case 2: {
374: void **ptr;
375: itdof = xm*ym;
377: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*deriv_type_size,&iarray_start);
378: ierr = PetscMemzero(iarray_start,xm*ym*deriv_type_size);
380: ptr = (void**)(iarray_start + xm*ym*deriv_type_size - ys*sizeof(void*));
381: for(j=ys;j<ys+ym;j++) {
382: ptr[j] = iarray_start + deriv_type_size*(xm*(j-ys) - xs);
383: }
384: *iptr = (void*)ptr;
385: break;}
386: case 3: {
387: void ***ptr,**bptr;
388: itdof = xm*ym*zm;
390: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*deriv_type_size,&iarray_start);
391: ierr = PetscMemzero(iarray_start,xm*ym*zm*deriv_type_size);
393: ptr = (void***)(iarray_start + xm*ym*zm*deriv_type_size - zs*sizeof(void*));
394: bptr = (void**)(iarray_start + xm*ym*zm*deriv_type_size + zm*sizeof(void**));
395: for(i=zs;i<zs+zm;i++) {
396: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
397: }
398: for (i=zs; i<zs+zm; i++) {
399: for (j=ys; j<ys+ym; j++) {
400: ptr[i][j] = iarray_start + deriv_type_size*(xm*ym*(i-zs) + xm*(j-ys) - xs);
401: }
402: }
404: *iptr = (void*)ptr;
405: break;}
406: default:
407: SETERRQ1(1,"Dimension %d not supported",da->dim);
408: }
410: done:
411: /* add arrays to the checked out list */
412: if (ghosted) {
413: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
414: if (!da->adarrayghostedout[i]) {
415: da->adarrayghostedout[i] = *iptr ;
416: da->adstartghostedout[i] = iarray_start;
417: da->ghostedtdof = itdof;
418: break;
419: }
420: }
421: } else {
422: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
423: if (!da->adarrayout[i]) {
424: da->adarrayout[i] = *iptr ;
425: da->adstartout[i] = iarray_start;
426: da->tdof = itdof;
427: break;
428: }
429: }
430: }
431: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
432: if (tdof) *tdof = itdof;
433: if (array_start) *array_start = iarray_start;
434: return(0);
435: }
437: #undef __FUNCT__
439: /*@C
440: DARestoreAdicArray - Restores an array of derivative types for a DA
441:
442: Input Parameter:
443: + da - information about my local patch
444: - ghosted - do you want arrays for the ghosted or nonghosted patch
446: Output Parameters:
447: + ptr - array data structured to be passed to ad_FormFunctionLocal()
448: . array_start - actual start of 1d array of all values that adiC can access directly
449: - tdof - total number of degrees of freedom represented in array_start
451: Level: advanced
453: .seealso: DAGetAdicArray()
455: @*/
456: int DARestoreAdicArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
457: {
458: int i;
459: void *iarray_start = 0;
460:
463: if (ghosted) {
464: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
465: if (da->adarrayghostedout[i] == *iptr) {
466: iarray_start = da->adstartghostedout[i];
467: da->adarrayghostedout[i] = PETSC_NULL;
468: da->adstartghostedout[i] = PETSC_NULL;
469: break;
470: }
471: }
472: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
473: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
474: if (!da->adarrayghostedin[i]){
475: da->adarrayghostedin[i] = *iptr;
476: da->adstartghostedin[i] = iarray_start;
477: break;
478: }
479: }
480: } else {
481: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
482: if (da->adarrayout[i] == *iptr) {
483: iarray_start = da->adstartout[i];
484: da->adarrayout[i] = PETSC_NULL;
485: da->adstartout[i] = PETSC_NULL;
486: break;
487: }
488: }
489: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
490: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
491: if (!da->adarrayin[i]){
492: da->adarrayin[i] = *iptr;
493: da->adstartin[i] = iarray_start;
494: break;
495: }
496: }
497: }
498: return(0);
499: }
501: #undef __FUNCT__
503: int ad_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
504: {
507: DAGetAdicArray(da,ghosted,iptr,0,0);
508: return(0);
509: }
511: #undef __FUNCT__
513: int ad_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
514: {
517: DARestoreAdicArray(da,ghosted,iptr,0,0);
518: return(0);
519: }
521: #endif
523: #undef __FUNCT__
525: /*@C
526: DAGetArray - Gets a work array for a DA
527:
528: Input Parameter:
529: + da - information about my local patch
530: - ghosted - do you want arrays for the ghosted or nonghosted patch
532: Output Parameters:
533: . ptr - array data structured
535: Level: advanced
537: .seealso: DARestoreArray(), DAGetAdicArray()
539: @*/
540: int DAGetArray(DA da,PetscTruth ghosted,void **iptr)
541: {
542: int ierr,j,i,xs,ys,xm,ym,zs,zm;
543: char *iarray_start;
547: if (ghosted) {
548: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
549: if (da->arrayghostedin[i]) {
550: *iptr = da->arrayghostedin[i];
551: iarray_start = (char*)da->startghostedin[i];
552: da->arrayghostedin[i] = PETSC_NULL;
553: da->startghostedin[i] = PETSC_NULL;
554:
555: goto done;
556: }
557: }
558: xs = da->Xs;
559: ys = da->Ys;
560: zs = da->Zs;
561: xm = da->Xe-da->Xs;
562: ym = da->Ye-da->Ys;
563: zm = da->Ze-da->Zs;
564: } else {
565: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
566: if (da->arrayin[i]) {
567: *iptr = da->arrayin[i];
568: iarray_start = (char*)da->startin[i];
569: da->arrayin[i] = PETSC_NULL;
570: da->startin[i] = PETSC_NULL;
571:
572: goto done;
573: }
574: }
575: xs = da->xs;
576: ys = da->ys;
577: zs = da->zs;
578: xm = da->xe-da->xs;
579: ym = da->ye-da->ys;
580: zm = da->ze-da->zs;
581: }
583: switch (da->dim) {
584: case 1: {
585: void *ptr;
587: ierr = PetscMalloc(xm*sizeof(PetscScalar),&iarray_start);
588: ierr = PetscMemzero(iarray_start,xm*sizeof(PetscScalar));
590: ptr = (void*)(iarray_start - xs*sizeof(PetscScalar));
591: *iptr = (void*)ptr;
592: break;}
593: case 2: {
594: void **ptr;
596: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*sizeof(PetscScalar),&iarray_start);
597: ierr = PetscMemzero(iarray_start,xm*ym*sizeof(PetscScalar));
599: ptr = (void**)(iarray_start + xm*ym*sizeof(PetscScalar) - ys*sizeof(void*));
600: for(j=ys;j<ys+ym;j++) {
601: ptr[j] = iarray_start + sizeof(PetscScalar)*(xm*(j-ys) - xs);
602: }
603: *iptr = (void*)ptr;
604: break;}
605: case 3: {
606: void ***ptr,**bptr;
608: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*sizeof(PetscScalar),&iarray_start);
609: ierr = PetscMemzero(iarray_start,xm*ym*zm*sizeof(PetscScalar));
611: ptr = (void***)(iarray_start + xm*ym*zm*sizeof(PetscScalar) - zs*sizeof(void*));
612: bptr = (void**)(iarray_start + xm*ym*zm*sizeof(PetscScalar) + zm*sizeof(void**));
613: for(i=zs;i<zs+zm;i++) {
614: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
615: }
616: for (i=zs; i<zs+zm; i++) {
617: for (j=ys; j<ys+ym; j++) {
618: ptr[i][j] = iarray_start + sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
619: }
620: }
622: *iptr = (void*)ptr;
623: break;}
624: default:
625: SETERRQ1(1,"Dimension %d not supported",da->dim);
626: }
628: done:
629: /* add arrays to the checked out list */
630: if (ghosted) {
631: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
632: if (!da->arrayghostedout[i]) {
633: da->arrayghostedout[i] = *iptr ;
634: da->startghostedout[i] = iarray_start;
635: break;
636: }
637: }
638: } else {
639: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
640: if (!da->arrayout[i]) {
641: da->arrayout[i] = *iptr ;
642: da->startout[i] = iarray_start;
643: break;
644: }
645: }
646: }
647: return(0);
648: }
650: #undef __FUNCT__
652: /*@C
653: DARestoreArray - Restores an array of derivative types for a DA
654:
655: Input Parameter:
656: + da - information about my local patch
657: . ghosted - do you want arrays for the ghosted or nonghosted patch
658: - ptr - array data structured to be passed to ad_FormFunctionLocal()
660: Level: advanced
662: .seealso: DAGetArray(), DAGetAdicArray()
664: @*/
665: int DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
666: {
667: int i;
668: void *iarray_start = 0;
669:
672: if (ghosted) {
673: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
674: if (da->arrayghostedout[i] == *iptr) {
675: iarray_start = da->startghostedout[i];
676: da->arrayghostedout[i] = PETSC_NULL;
677: da->startghostedout[i] = PETSC_NULL;
678: break;
679: }
680: }
681: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
682: if (!da->arrayghostedin[i]){
683: da->arrayghostedin[i] = *iptr;
684: da->startghostedin[i] = iarray_start;
685: break;
686: }
687: }
688: } else {
689: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
690: if (da->arrayout[i] == *iptr) {
691: iarray_start = da->startout[i];
692: da->arrayout[i] = PETSC_NULL;
693: da->startout[i] = PETSC_NULL;
694: break;
695: }
696: }
697: for (i=0; i<DA_MAX_WORK_ARRAYS; i++) {
698: if (!da->arrayin[i]){
699: da->arrayin[i] = *iptr;
700: da->startin[i] = iarray_start;
701: break;
702: }
703: }
704: }
705: return(0);
706: }
708: #undef __FUNCT__
710: /*@C
711: DAGetAdicMFArray - Gets an array of derivative types for a DA for matrix-free ADIC.
712:
713: Input Parameter:
714: + da - information about my local patch
715: - ghosted - do you want arrays for the ghosted or nonghosted patch?
717: Output Parameters:
718: + ptr - array data structured to be passed to ad_FormFunctionLocal()
719: . array_start - actual start of 1d array of all values that adiC can access directly (may be null)
720: - tdof - total number of degrees of freedom represented in array_start (may be null)
722: Notes:
723: This routine returns the same type of object as the DAVecGetArray(), except its
724: elements are derivative types instead of PetscScalars.
726: Level: advanced
728: .seealso: DARestoreAdicMFArray(), DAGetArray(), DAGetAdicArray()
730: @*/
731: int DAGetAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
732: {
733: int ierr,j,i,xs,ys,xm,ym,zs,zm,itdof;
734: char *iarray_start;
738: if (ghosted) {
739: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
740: if (da->admfarrayghostedin[i]) {
741: *iptr = da->admfarrayghostedin[i];
742: iarray_start = (char*)da->admfstartghostedin[i];
743: itdof = da->ghostedtdof;
744: da->admfarrayghostedin[i] = PETSC_NULL;
745: da->admfstartghostedin[i] = PETSC_NULL;
746:
747: goto done;
748: }
749: }
750: xs = da->Xs;
751: ys = da->Ys;
752: zs = da->Zs;
753: xm = da->Xe-da->Xs;
754: ym = da->Ye-da->Ys;
755: zm = da->Ze-da->Zs;
756: } else {
757: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
758: if (da->admfarrayin[i]) {
759: *iptr = da->admfarrayin[i];
760: iarray_start = (char*)da->admfstartin[i];
761: itdof = da->tdof;
762: da->admfarrayin[i] = PETSC_NULL;
763: da->admfstartin[i] = PETSC_NULL;
764:
765: goto done;
766: }
767: }
768: xs = da->xs;
769: ys = da->ys;
770: zs = da->zs;
771: xm = da->xe-da->xs;
772: ym = da->ye-da->ys;
773: zm = da->ze-da->zs;
774: }
776: switch (da->dim) {
777: case 1: {
778: void *ptr;
779: itdof = xm;
781: ierr = PetscMalloc(xm*2*sizeof(PetscScalar),&iarray_start);
782: ierr = PetscMemzero(iarray_start,xm*2*sizeof(PetscScalar));
784: ptr = (void*)(iarray_start - xs*2*sizeof(PetscScalar));
785: *iptr = (void*)ptr;
786: break;}
787: case 2: {
788: void **ptr;
789: itdof = xm*ym;
791: ierr = PetscMalloc((ym+1)*sizeof(void *)+xm*ym*2*sizeof(PetscScalar),&iarray_start);
792: ierr = PetscMemzero(iarray_start,xm*ym*2*sizeof(PetscScalar));
794: ptr = (void**)(iarray_start + xm*ym*2*sizeof(PetscScalar) - ys*sizeof(void*));
795: for(j=ys;j<ys+ym;j++) {
796: ptr[j] = iarray_start + 2*sizeof(PetscScalar)*(xm*(j-ys) - xs);
797: }
798: *iptr = (void*)ptr;
799: break;}
800: case 3: {
801: void ***ptr,**bptr;
802: itdof = xm*ym*zm;
804: ierr = PetscMalloc((zm+1)*sizeof(void **)+(ym*zm+1)*sizeof(void*)+xm*ym*zm*2*sizeof(PetscScalar),&iarray_start);
805: ierr = PetscMemzero(iarray_start,xm*ym*zm*2*sizeof(PetscScalar));
807: ptr = (void***)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) - zs*sizeof(void*));
808: bptr = (void**)(iarray_start + xm*ym*zm*2*sizeof(PetscScalar) + zm*sizeof(void**));
809: for(i=zs;i<zs+zm;i++) {
810: ptr[i] = bptr + ((i-zs)*ym* - ys)*sizeof(void*);
811: }
812: for (i=zs; i<zs+zm; i++) {
813: for (j=ys; j<ys+ym; j++) {
814: ptr[i][j] = iarray_start + 2*sizeof(PetscScalar)*(xm*ym*(i-zs) + xm*(j-ys) - xs);
815: }
816: }
818: *iptr = (void*)ptr;
819: break;}
820: default:
821: SETERRQ1(1,"Dimension %d not supported",da->dim);
822: }
824: done:
825: /* add arrays to the checked out list */
826: if (ghosted) {
827: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
828: if (!da->admfarrayghostedout[i]) {
829: da->admfarrayghostedout[i] = *iptr ;
830: da->admfstartghostedout[i] = iarray_start;
831: da->ghostedtdof = itdof;
832: break;
833: }
834: }
835: } else {
836: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
837: if (!da->admfarrayout[i]) {
838: da->admfarrayout[i] = *iptr ;
839: da->admfstartout[i] = iarray_start;
840: da->tdof = itdof;
841: break;
842: }
843: }
844: }
845: if (i == DA_MAX_AD_ARRAYS+1) SETERRQ(1,"Too many DA ADIC arrays obtained");
846: if (tdof) *tdof = itdof;
847: if (array_start) *array_start = iarray_start;
848: return(0);
849: }
851: #undef __FUNCT__
853: /*@C
854: DARestoreAdicMFArray - Restores an array of derivative types for a DA.
855:
856: Input Parameter:
857: + da - information about my local patch
858: - ghosted - do you want arrays for the ghosted or nonghosted patch?
860: Output Parameters:
861: + ptr - array data structure to be passed to ad_FormFunctionLocal()
862: . array_start - actual start of 1d array of all values that adiC can access directly
863: - tdof - total number of degrees of freedom represented in array_start
865: Level: advanced
867: .seealso: DAGetAdicArray()
869: @*/
870: int DARestoreAdicMFArray(DA da,PetscTruth ghosted,void **iptr,void **array_start,int *tdof)
871: {
872: int i;
873: void *iarray_start = 0;
874:
877: if (ghosted) {
878: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
879: if (da->admfarrayghostedout[i] == *iptr) {
880: iarray_start = da->admfstartghostedout[i];
881: da->admfarrayghostedout[i] = PETSC_NULL;
882: da->admfstartghostedout[i] = PETSC_NULL;
883: break;
884: }
885: }
886: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
887: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
888: if (!da->admfarrayghostedin[i]){
889: da->admfarrayghostedin[i] = *iptr;
890: da->admfstartghostedin[i] = iarray_start;
891: break;
892: }
893: }
894: } else {
895: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
896: if (da->admfarrayout[i] == *iptr) {
897: iarray_start = da->admfstartout[i];
898: da->admfarrayout[i] = PETSC_NULL;
899: da->admfstartout[i] = PETSC_NULL;
900: break;
901: }
902: }
903: if (!iarray_start) SETERRQ(1,"Could not find array in checkout list");
904: for (i=0; i<DA_MAX_AD_ARRAYS; i++) {
905: if (!da->admfarrayin[i]){
906: da->admfarrayin[i] = *iptr;
907: da->admfstartin[i] = iarray_start;
908: break;
909: }
910: }
911: }
912: return(0);
913: }
915: #undef __FUNCT__
917: int admf_DAGetArray(DA da,PetscTruth ghosted,void **iptr)
918: {
921: DAGetAdicMFArray(da,ghosted,iptr,0,0);
922: return(0);
923: }
925: #undef __FUNCT__
927: int admf_DARestoreArray(DA da,PetscTruth ghosted,void **iptr)
928: {
931: DARestoreAdicMFArray(da,ghosted,iptr,0,0);
932: return(0);
933: }