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