xref: /petsc/src/dm/impls/da/daindex.c (revision bebe2cf65d55febe21a5af8db2bd2e168caaa2e7)
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 #undef __FUNCT__
9 #define __FUNCT__ "DMDAGetNatural_Private"
10 /*
11    Gets the natural number for each global number on the process.
12 
13    Used by DMDAGetAO() and DMDAGlobalToNatural_Create()
14 */
15 PetscErrorCode DMDAGetNatural_Private(DM da,PetscInt *outNlocal,IS *isnatural)
16 {
17   PetscErrorCode ierr;
18   PetscInt       Nlocal,i,j,k,*lidx,lict = 0,dim = da->dim;
19   DM_DA          *dd = (DM_DA*)da->data;
20 
21   PetscFunctionBegin;
22   Nlocal = (dd->xe-dd->xs);
23   if (dim > 1) Nlocal *= (dd->ye-dd->ys);
24   if (dim > 2) Nlocal *= (dd->ze-dd->zs);
25 
26   ierr = PetscMalloc1(Nlocal,&lidx);CHKERRQ(ierr);
27 
28   if (dim == 1) {
29     for (i=dd->xs; i<dd->xe; i++) {
30       /*  global number in natural ordering */
31       lidx[lict++] = i;
32     }
33   } else if (dim == 2) {
34     for (j=dd->ys; j<dd->ye; j++) {
35       for (i=dd->xs; i<dd->xe; i++) {
36         /*  global number in natural ordering */
37         lidx[lict++] = i + j*dd->M*dd->w;
38       }
39     }
40   } else if (dim == 3) {
41     for (k=dd->zs; k<dd->ze; k++) {
42       for (j=dd->ys; j<dd->ye; j++) {
43         for (i=dd->xs; i<dd->xe; i++) {
44           lidx[lict++] = i + j*dd->M*dd->w + k*dd->M*dd->N*dd->w;
45         }
46       }
47     }
48   }
49   *outNlocal = Nlocal;
50   ierr       = ISCreateGeneral(PetscObjectComm((PetscObject)da),Nlocal,lidx,PETSC_OWN_POINTER,isnatural);CHKERRQ(ierr);
51   PetscFunctionReturn(0);
52 }
53 
54 #undef __FUNCT__
55 #define __FUNCT__ "DMDASetAOType"
56 /*@
57    DMDASetAOType - Sets the type of application ordering for a distributed array.
58 
59    Collective on DMDA
60 
61    Input Parameter:
62 .  da - the distributed array
63 .  aotype - type of AO
64 
65    Output Parameters:
66 
67    Level: intermediate
68 
69    Notes:
70    It will generate and error if an AO has already been obtained with a call to DMDAGetAO and the user sets a different AOType
71 
72 .keywords: distributed array, get, global, indices, local-to-global
73 
74 .seealso: DMDACreate2d(), DMDAGetAO(), DMDAGetGhostCorners(), DMDAGetCorners(), DMDALocalToGlocal()
75           DMGlobalToLocalBegin(), DMGlobalToLocalEnd(), DMLocalToLocalBegin(), DMLocalToLocalEnd(), DMDAGetGlobalIndices(), DMDAGetOwnershipRanges(),
76           AO, AOPetscToApplication(), AOApplicationToPetsc()
77 @*/
78 PetscErrorCode  DMDASetAOType(DM da,AOType aotype)
79 {
80   DM_DA          *dd;
81   PetscBool      isdmda;
82   PetscErrorCode ierr;
83 
84   PetscFunctionBegin;
85   PetscValidHeaderSpecific(da,DM_CLASSID,1);
86   ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isdmda);CHKERRQ(ierr);
87   if (!isdmda) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_SUP,"Requires a DMDA as input");
88   /* now we can safely dereference */
89   dd = (DM_DA*)da->data;
90   if (dd->ao) { /* check if the already computed AO has the same type as requested */
91     PetscBool match;
92     ierr = PetscObjectTypeCompare((PetscObject)dd->ao,aotype,&match);CHKERRQ(ierr);
93     if (!match) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_SUP,"Cannot change AO type");
94     PetscFunctionReturn(0);
95   }
96   ierr = PetscFree(dd->aotype);CHKERRQ(ierr);
97   ierr = PetscStrallocpy(aotype,(char**)&dd->aotype);CHKERRQ(ierr);
98   PetscFunctionReturn(0);
99 }
100 
101 #undef __FUNCT__
102 #define __FUNCT__ "DMDAGetAO"
103 /*@
104    DMDAGetAO - Gets the application ordering context for a distributed array.
105 
106    Collective on DMDA
107 
108    Input Parameter:
109 .  da - the distributed array
110 
111    Output Parameters:
112 .  ao - the application ordering context for DMDAs
113 
114    Level: intermediate
115 
116    Notes:
117    In this case, the AO maps to the natural grid ordering that would be used
118    for the DMDA if only 1 processor were employed (ordering most rapidly in the
119    x-direction, then y, then z).  Multiple degrees of freedom are numbered
120    for each node (rather than 1 component for the whole grid, then the next
121    component, etc.)
122 
123 .keywords: distributed array, get, global, indices, local-to-global
124 
125 .seealso: DMDACreate2d(), DMDASetAOType(), DMDAGetGhostCorners(), DMDAGetCorners(), DMDALocalToGlocal()
126           DMGlobalToLocalBegin(), DMGlobalToLocalEnd(), DMLocalToLocalBegin(), DMLocalToLocalEnd(),  DMDAGetOwnershipRanges(),
127           AO, AOPetscToApplication(), AOApplicationToPetsc()
128 @*/
129 PetscErrorCode  DMDAGetAO(DM da,AO *ao)
130 {
131   DM_DA          *dd;
132   PetscBool      isdmda;
133   PetscErrorCode ierr;
134 
135   PetscFunctionBegin;
136   PetscValidHeaderSpecific(da,DM_CLASSID,1);
137   PetscValidPointer(ao,2);
138   ierr = PetscObjectTypeCompare((PetscObject)da,DMDA,&isdmda);CHKERRQ(ierr);
139   if (!isdmda) SETERRQ(PetscObjectComm((PetscObject)da),PETSC_ERR_SUP,"Requires a DMDA as input");
140   /* now we can safely dereference */
141   dd = (DM_DA*)da->data;
142 
143   /*
144      Build the natural ordering to PETSc ordering mappings.
145   */
146   if (!dd->ao) {
147     IS             ispetsc,isnatural;
148     PetscErrorCode ierr;
149     PetscInt       Nlocal;
150 
151     ierr = DMDAGetNatural_Private(da,&Nlocal,&isnatural);CHKERRQ(ierr);
152     ierr = ISCreateStride(PetscObjectComm((PetscObject)da),Nlocal,dd->base,1,&ispetsc);CHKERRQ(ierr);
153     ierr = AOCreate(PetscObjectComm((PetscObject)da),&dd->ao);CHKERRQ(ierr);
154     ierr = AOSetIS(dd->ao,isnatural,ispetsc);CHKERRQ(ierr);
155     ierr = AOSetType(dd->ao,dd->aotype);CHKERRQ(ierr);
156     ierr = PetscLogObjectParent((PetscObject)da,(PetscObject)dd->ao);CHKERRQ(ierr);
157     ierr = ISDestroy(&ispetsc);CHKERRQ(ierr);
158     ierr = ISDestroy(&isnatural);CHKERRQ(ierr);
159   }
160   *ao = dd->ao;
161   PetscFunctionReturn(0);
162 }
163 
164 
165