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