Actual source code: ex19.c


  2: static char help[]= "Test leaf sorting in PetscSFSetGraph()\n\n";

  4: #include <petscsf.h>

  6: typedef struct {
  7:   MPI_Comm          comm;
  8:   PetscMPIInt       rank, size;
  9:   PetscInt          leaveStep, nLeavesPerRank;
 10:   PetscBool         contiguousLeaves;
 11:   PetscCopyMode     localmode, remotemode;
 12:   PetscInt         *ilocal;
 13:   PetscSFNode      *iremote;
 14: } AppCtx;

 16: static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx)
 17: {
 18:   ctx->comm = comm;
 19:   ctx->nLeavesPerRank = 4;
 20:   ctx->leaveStep = 1;
 21:   ctx->contiguousLeaves = PETSC_FALSE;
 22:   ctx->localmode = PETSC_OWN_POINTER;
 23:   ctx->remotemode = PETSC_OWN_POINTER;
 24:   ctx->ilocal = NULL;
 25:   ctx->iremote = NULL;
 26:   PetscOptionsGetInt(NULL, NULL, "-n_leaves_per_rank", &ctx->nLeavesPerRank, NULL);
 27:   PetscOptionsGetInt(NULL, NULL, "-leave_step", &ctx->leaveStep, NULL);
 28:   PetscOptionsGetEnum(NULL, NULL, "-localmode", PetscCopyModes, (PetscEnum*) &ctx->localmode, NULL);
 29:   PetscOptionsGetEnum(NULL, NULL, "-remotemode", PetscCopyModes, (PetscEnum*) &ctx->remotemode, NULL);
 30:   ctx->contiguousLeaves = (PetscBool) (ctx->leaveStep == 1);
 31:   MPI_Comm_size(comm, &ctx->size);
 32:   MPI_Comm_rank(comm, &ctx->rank);
 33:   return 0;
 34: }

 36: static PetscErrorCode PetscSFCheckEqual_Private(PetscSF sf0, PetscSF sf1)
 37: {
 38:   PetscInt          nRoot, nLeave;
 39:   Vec               vecRoot0, vecLeave0, vecRoot1, vecLeave1;
 40:   MPI_Comm          comm;
 41:   PetscBool         flg;

 43:   PetscObjectGetComm((PetscObject)sf0, &comm);
 44:   PetscSFGetGraph(sf0, &nRoot, NULL, NULL, NULL);
 45:   PetscSFGetLeafRange(sf0, NULL, &nLeave);
 46:   nLeave++;
 47:   VecCreateMPI(comm, nRoot, PETSC_DECIDE, &vecRoot0);
 48:   VecCreateMPI(comm, nLeave, PETSC_DECIDE, &vecLeave0);
 49:   VecDuplicate(vecRoot0, &vecRoot1);
 50:   VecDuplicate(vecLeave0, &vecLeave1);
 51:   {
 52:     PetscRandom       rand;

 54:     PetscRandomCreate(comm, &rand);
 55:     PetscRandomSetFromOptions(rand);
 56:     VecSetRandom(vecRoot0, rand);
 57:     VecSetRandom(vecLeave0, rand);
 58:     VecCopy(vecRoot0, vecRoot1);
 59:     VecCopy(vecLeave0, vecLeave1);
 60:     PetscRandomDestroy(&rand);
 61:   }

 63:   VecScatterBegin(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);
 64:   VecScatterEnd(  sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);
 65:   VecScatterBegin(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);
 66:   VecScatterEnd(  sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);
 67:   VecEqual(vecLeave0, vecLeave1, &flg);

 70:   VecScatterBegin(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);
 71:   VecScatterEnd(  sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);
 72:   VecScatterBegin(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);
 73:   VecScatterEnd(  sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);
 74:   VecEqual(vecRoot0, vecRoot1, &flg);

 77:   VecDestroy(&vecRoot0);
 78:   VecDestroy(&vecRoot1);
 79:   VecDestroy(&vecLeave0);
 80:   VecDestroy(&vecLeave1);
 81:   return 0;
 82: }

 84: PetscErrorCode CreateSF0(AppCtx *ctx, PetscSF *sf0)
 85: {
 86:   PetscInt          j, k, r;
 87:   PetscInt          nLeaves = ctx->nLeavesPerRank * ctx->size;
 88:   PetscInt          nroots  = ctx->nLeavesPerRank;
 89:   PetscSF           sf;
 90:   PetscInt         *ilocal;
 91:   PetscSFNode      *iremote;

 93:   PetscMalloc1(nLeaves+1, &ctx->ilocal);
 94:   PetscMalloc1(nLeaves, &ctx->iremote);
 95:   ilocal = ctx->ilocal;
 96:   iremote = ctx->iremote;
 97:   ilocal[nLeaves] = -ctx->leaveStep;
 98:   PetscSFCreate(ctx->comm, &sf);
 99:   for (r=0, j=nLeaves-1; r<ctx->size; r++) {
100:     for (k=0; k<ctx->nLeavesPerRank; k++, j--) {
101:       ilocal[j] = ilocal[j+1] + ctx->leaveStep;
102:       iremote[j].rank = r;
103:       iremote[j].index = k;
104:     }
105:   }
106:   PetscSFSetGraph(sf, nroots, nLeaves, ilocal, ctx->localmode, iremote, ctx->remotemode);
107:   {
108:     const PetscInt *tlocal;
109:     PetscBool       sorted;

111:     PetscSFGetGraph(sf, NULL, NULL, &tlocal, NULL);
113:     if (tlocal) {
114:       PetscSortedInt(nLeaves, tlocal, &sorted);
116:     }
117:   }
118:   *sf0 = sf;
119:   return 0;
120: }

122: PetscErrorCode CreateSF1(AppCtx *ctx, PetscSF *sf1)
123: {
124:   PetscInt          j, k, r;
125:   PetscInt         *ilocal = NULL;
126:   PetscSFNode      *iremote;
127:   PetscInt          nLeaves = ctx->nLeavesPerRank * ctx->size;
128:   PetscInt          nroots  = ctx->nLeavesPerRank;
129:   PetscSF           sf;

131:   ilocal = NULL;
132:   if (!ctx->contiguousLeaves) {
133:     PetscCalloc1(nLeaves+1, &ilocal);
134:   }
135:   PetscMalloc1(nLeaves, &iremote);
136:   PetscSFCreate(ctx->comm, &sf);
137:   for (r=0, j=0; r<ctx->size; r++) {
138:     for (k=0; k<ctx->nLeavesPerRank; k++, j++) {
139:       if (!ctx->contiguousLeaves) {
140:         ilocal[j+1] = ilocal[j] + ctx->leaveStep;
141:       }
142:       iremote[j].rank = r;
143:       iremote[j].index = k;
144:     }
145:   }
147:   PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);
148:   if (ctx->contiguousLeaves) {
149:     const PetscInt *tlocal;

151:     PetscSFGetGraph(sf, NULL, NULL, &tlocal, NULL);
153:   }
154:   *sf1 = sf;
155:   return 0;
156: }

158: int main(int argc, char **argv)
159: {
160:   AppCtx   ctx;
161:   PetscSF  sf0, sf1;
162:   MPI_Comm comm;

164:   PetscInitialize(&argc,&argv,NULL,help);
165:   comm = PETSC_COMM_WORLD;
166:   GetOptions(comm, &ctx);

168:   CreateSF0(&ctx, &sf0);
169:   CreateSF1(&ctx, &sf1);
170:   PetscSFViewFromOptions(sf0, NULL, "-sf0_view");
171:   PetscSFViewFromOptions(sf1, NULL, "-sf1_view");
172:   PetscSFCheckEqual_Private(sf0, sf1);

174:   if (ctx.localmode != PETSC_OWN_POINTER)  PetscFree(ctx.ilocal);
175:   if (ctx.remotemode != PETSC_OWN_POINTER) PetscFree(ctx.iremote);
176:   PetscSFDestroy(&sf0);
177:   PetscSFDestroy(&sf1);
178:   PetscFinalize();
179:   return 0;
180: }

182: /*TEST
183:   testset:
184:     suffix: 1
185:     nsize: {{1 3}}
186:     args: -n_leaves_per_rank {{0 5}} -leave_step {{1 3}}
187:     test:
188:       suffix: a
189:       args: -localmode {{COPY_VALUES OWN_POINTER}} -remotemode {{COPY_VALUES OWN_POINTER}}
190:     test:
191:       suffix: b
192:       args: -localmode USE_POINTER -remotemode {{COPY_VALUES OWN_POINTER USE_POINTER}}
193: TEST*/