1 2 /* 3 Code for manipulating distributed regular arrays in parallel. 4 */ 5 6 #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/ 7 8 /* 9 DMLocalToLocalCreate_DA - Creates the local to local scatter 10 11 Collective on da 12 13 Input Parameter: 14 . da - the distributed array 15 16 */ 17 PetscErrorCode DMLocalToLocalCreate_DA(DM da) 18 { 19 PetscInt *idx,left,j,count,up,down,i,bottom,top,k,dim=da->dim; 20 DM_DA *dd = (DM_DA*)da->data; 21 22 PetscFunctionBegin; 23 PetscValidHeaderSpecific(da,DM_CLASSID,1); 24 25 if (dd->ltol) PetscFunctionReturn(0); 26 /* 27 We simply remap the values in the from part of 28 global to local to read from an array with the ghost values 29 rather then from the plain array. 30 */ 31 PetscCall(VecScatterCopy(dd->gtol,&dd->ltol)); 32 PetscCall(PetscLogObjectParent((PetscObject)da,(PetscObject)dd->ltol)); 33 if (dim == 1) { 34 left = dd->xs - dd->Xs; 35 PetscCall(PetscMalloc1(dd->xe-dd->xs,&idx)); 36 for (j=0; j<dd->xe-dd->xs; j++) idx[j] = left + j; 37 } else if (dim == 2) { 38 left = dd->xs - dd->Xs; down = dd->ys - dd->Ys; up = down + dd->ye-dd->ys; 39 PetscCall(PetscMalloc1((dd->xe-dd->xs)*(up - down),&idx)); 40 count = 0; 41 for (i=down; i<up; i++) { 42 for (j=0; j<dd->xe-dd->xs; j++) { 43 idx[count++] = left + i*(dd->Xe-dd->Xs) + j; 44 } 45 } 46 } else if (dim == 3) { 47 left = dd->xs - dd->Xs; 48 bottom = dd->ys - dd->Ys; top = bottom + dd->ye-dd->ys; 49 down = dd->zs - dd->Zs; up = down + dd->ze-dd->zs; 50 count = (dd->xe-dd->xs)*(top-bottom)*(up-down); 51 PetscCall(PetscMalloc1(count,&idx)); 52 count = 0; 53 for (i=down; i<up; i++) { 54 for (j=bottom; j<top; j++) { 55 for (k=0; k<dd->xe-dd->xs; k++) { 56 idx[count++] = (left+j*(dd->Xe-dd->Xs))+i*(dd->Xe-dd->Xs)*(dd->Ye-dd->Ys) + k; 57 } 58 } 59 } 60 } else SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_ARG_CORRUPT,"DMDA has invalid dimension %" PetscInt_FMT,dim); 61 62 PetscCall(VecScatterRemap(dd->ltol,idx,NULL)); 63 PetscCall(PetscFree(idx)); 64 PetscFunctionReturn(0); 65 } 66 67 /* 68 DMLocalToLocalBegin_DA - Maps from a local vector (including ghost points 69 that contain irrelevant values) to another local vector where the ghost 70 points in the second are set correctly. Must be followed by DMLocalToLocalEnd_DA(). 71 72 Neighbor-wise Collective on da 73 74 Input Parameters: 75 + da - the distributed array context 76 . g - the original local vector 77 - mode - one of INSERT_VALUES or ADD_VALUES 78 79 Output Parameter: 80 . l - the local vector with correct ghost values 81 82 Notes: 83 The local vectors used here need not be the same as those 84 obtained from DMCreateLocalVector(), BUT they 85 must have the same parallel data layout; they could, for example, be 86 obtained with VecDuplicate() from the DMDA originating vectors. 87 88 */ 89 PetscErrorCode DMLocalToLocalBegin_DA(DM da,Vec g,InsertMode mode,Vec l) 90 { 91 DM_DA *dd = (DM_DA*)da->data; 92 93 PetscFunctionBegin; 94 PetscValidHeaderSpecific(da,DM_CLASSID,1); 95 if (!dd->ltol) { 96 PetscCall(DMLocalToLocalCreate_DA(da)); 97 } 98 PetscCall(VecScatterBegin(dd->ltol,g,l,mode,SCATTER_FORWARD)); 99 PetscFunctionReturn(0); 100 } 101 102 /* 103 DMLocalToLocalEnd_DA - Maps from a local vector (including ghost points 104 that contain irrelevant values) to another local vector where the ghost 105 points in the second are set correctly. Must be preceded by 106 DMLocalToLocalBegin_DA(). 107 108 Neighbor-wise Collective on da 109 110 Input Parameters: 111 + da - the distributed array context 112 . g - the original local vector 113 - mode - one of INSERT_VALUES or ADD_VALUES 114 115 Output Parameter: 116 . l - the local vector with correct ghost values 117 118 Note: 119 The local vectors used here need not be the same as those 120 obtained from DMCreateLocalVector(), BUT they 121 must have the same parallel data layout; they could, for example, be 122 obtained with VecDuplicate() from the DMDA originating vectors. 123 124 */ 125 PetscErrorCode DMLocalToLocalEnd_DA(DM da,Vec g,InsertMode mode,Vec l) 126 { 127 DM_DA *dd = (DM_DA*)da->data; 128 129 PetscFunctionBegin; 130 PetscValidHeaderSpecific(da,DM_CLASSID,1); 131 PetscValidHeaderSpecific(g,VEC_CLASSID,2); 132 PetscCall(VecScatterEnd(dd->ltol,g,l,mode,SCATTER_FORWARD)); 133 PetscFunctionReturn(0); 134 } 135