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