Actual source code: daltol.c


  2: /*
  3:   Code for manipulating distributed regular arrays in parallel.
  4: */

  6: #include <petsc/private/dmdaimpl.h>

  8: /*
  9:    DMLocalToLocalCreate_DA - Creates the local to local scatter

 11:    Collective on da

 13:    Input Parameter:
 14: .  da - the distributed array

 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;

 22:   PetscFunctionBegin;

 25:   if (dd->ltol) PetscFunctionReturn(PETSC_SUCCESS);
 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:   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);

 61:   PetscCall(VecScatterRemap(dd->ltol, idx, NULL));
 62:   PetscCall(PetscFree(idx));
 63:   PetscFunctionReturn(PETSC_SUCCESS);
 64: }

 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().

 71:    Neighbor-wise Collective on da

 73:    Input Parameters:
 74: +  da - the distributed array context
 75: .  g - the original local vector
 76: -  mode - one of INSERT_VALUES or ADD_VALUES

 78:    Output Parameter:
 79: .  l  - the local vector with correct ghost values

 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.

 87: */
 88: PetscErrorCode DMLocalToLocalBegin_DA(DM da, Vec g, InsertMode mode, Vec l)
 89: {
 90:   DM_DA *dd = (DM_DA *)da->data;

 92:   PetscFunctionBegin;
 94:   if (!dd->ltol) PetscCall(DMLocalToLocalCreate_DA(da));
 95:   PetscCall(VecScatterBegin(dd->ltol, g, l, mode, SCATTER_FORWARD));
 96:   PetscFunctionReturn(PETSC_SUCCESS);
 97: }

 99: /*
100:    DMLocalToLocalEnd_DA - Maps from a local vector (including ghost points
101:    that contain irrelevant values) to another local vector where the ghost
102:    points in the second are set correctly.  Must be preceded by
103:    DMLocalToLocalBegin_DA().

105:    Neighbor-wise Collective on da

107:    Input Parameters:
108: +  da - the distributed array context
109: .  g - the original local vector
110: -  mode - one of INSERT_VALUES or ADD_VALUES

112:    Output Parameter:
113: .  l  - the local vector with correct ghost values

115:    Note:
116:    The local vectors used here need not be the same as those
117:    obtained from DMCreateLocalVector(), BUT they
118:    must have the same parallel data layout; they could, for example, be
119:    obtained with VecDuplicate() from the DMDA originating vectors.

121: */
122: PetscErrorCode DMLocalToLocalEnd_DA(DM da, Vec g, InsertMode mode, Vec l)
123: {
124:   DM_DA *dd = (DM_DA *)da->data;

126:   PetscFunctionBegin;
129:   PetscCall(VecScatterEnd(dd->ltol, g, l, mode, SCATTER_FORWARD));
130:   PetscFunctionReturn(PETSC_SUCCESS);
131: }