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