1 2 static char help[]= "Test leaf sorting in PetscSFSetGraph()\n\n"; 3 4 #include <petscsf.h> 5 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; 15 16 static PetscErrorCode GetOptions(MPI_Comm comm, AppCtx *ctx) 17 { 18 PetscErrorCode ierr; 19 20 PetscFunctionBegin; 21 ctx->comm = comm; 22 ctx->nLeavesPerRank = 4; 23 ctx->leaveStep = 1; 24 ctx->contiguousLeaves = PETSC_FALSE; 25 ctx->localmode = PETSC_OWN_POINTER; 26 ctx->remotemode = PETSC_OWN_POINTER; 27 ctx->ilocal = NULL; 28 ctx->iremote = NULL; 29 ierr = PetscOptionsGetInt(NULL, NULL, "-n_leaves_per_rank", &ctx->nLeavesPerRank, NULL);CHKERRQ(ierr); 30 ierr = PetscOptionsGetInt(NULL, NULL, "-leave_step", &ctx->leaveStep, NULL);CHKERRQ(ierr); 31 ierr = PetscOptionsGetEnum(NULL, NULL, "-localmode", PetscCopyModes, (PetscEnum*) &ctx->localmode, NULL);CHKERRQ(ierr); 32 ierr = PetscOptionsGetEnum(NULL, NULL, "-remotemode", PetscCopyModes, (PetscEnum*) &ctx->remotemode, NULL);CHKERRQ(ierr); 33 ctx->contiguousLeaves = (PetscBool) (ctx->leaveStep == 1); 34 ierr = MPI_Comm_size(comm, &ctx->size);CHKERRMPI(ierr); 35 ierr = MPI_Comm_rank(comm, &ctx->rank);CHKERRMPI(ierr); 36 PetscFunctionReturn(0); 37 } 38 39 static PetscErrorCode PetscSFCheckEqual_Private(PetscSF sf0, PetscSF sf1) 40 { 41 PetscInt nRoot, nLeave; 42 Vec vecRoot0, vecLeave0, vecRoot1, vecLeave1; 43 MPI_Comm comm; 44 PetscBool flg; 45 PetscErrorCode ierr; 46 47 PetscFunctionBegin; 48 ierr = PetscObjectGetComm((PetscObject)sf0, &comm);CHKERRQ(ierr); 49 ierr = PetscSFGetGraph(sf0, &nRoot, NULL, NULL, NULL);CHKERRQ(ierr); 50 ierr = PetscSFGetLeafRange(sf0, NULL, &nLeave);CHKERRQ(ierr); 51 nLeave++; 52 ierr = VecCreateMPI(comm, nRoot, PETSC_DECIDE, &vecRoot0);CHKERRQ(ierr); 53 ierr = VecCreateMPI(comm, nLeave, PETSC_DECIDE, &vecLeave0);CHKERRQ(ierr); 54 ierr = VecDuplicate(vecRoot0, &vecRoot1);CHKERRQ(ierr); 55 ierr = VecDuplicate(vecLeave0, &vecLeave1);CHKERRQ(ierr); 56 { 57 PetscRandom rand; 58 59 ierr = PetscRandomCreate(comm, &rand);CHKERRQ(ierr); 60 ierr = PetscRandomSetFromOptions(rand);CHKERRQ(ierr); 61 ierr = VecSetRandom(vecRoot0, rand);CHKERRQ(ierr); 62 ierr = VecSetRandom(vecLeave0, rand);CHKERRQ(ierr); 63 ierr = VecCopy(vecRoot0, vecRoot1);CHKERRQ(ierr); 64 ierr = VecCopy(vecLeave0, vecLeave1);CHKERRQ(ierr); 65 ierr = PetscRandomDestroy(&rand);CHKERRQ(ierr); 66 } 67 68 ierr = VecScatterBegin(sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 69 ierr = VecScatterEnd( sf0, vecRoot0, vecLeave0, ADD_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 70 ierr = VecScatterBegin(sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 71 ierr = VecScatterEnd( sf1, vecRoot1, vecLeave1, ADD_VALUES, SCATTER_FORWARD);CHKERRQ(ierr); 72 ierr = VecEqual(vecLeave0, vecLeave1, &flg);CHKERRQ(ierr); 73 PetscCheck(flg, comm, PETSC_ERR_PLIB, "leave vectors differ"); 74 75 ierr = VecScatterBegin(sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);CHKERRQ(ierr); 76 ierr = VecScatterEnd( sf0, vecLeave0, vecRoot0, ADD_VALUES, SCATTER_REVERSE);CHKERRQ(ierr); 77 ierr = VecScatterBegin(sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);CHKERRQ(ierr); 78 ierr = VecScatterEnd( sf1, vecLeave1, vecRoot1, ADD_VALUES, SCATTER_REVERSE);CHKERRQ(ierr); 79 ierr = VecEqual(vecRoot0, vecRoot1, &flg);CHKERRQ(ierr); 80 PetscCheck(flg, comm, PETSC_ERR_PLIB, "root vectors differ"); 81 82 ierr = VecDestroy(&vecRoot0);CHKERRQ(ierr); 83 ierr = VecDestroy(&vecRoot1);CHKERRQ(ierr); 84 ierr = VecDestroy(&vecLeave0);CHKERRQ(ierr); 85 ierr = VecDestroy(&vecLeave1);CHKERRQ(ierr); 86 PetscFunctionReturn(0); 87 } 88 89 PetscErrorCode CreateSF0(AppCtx *ctx, PetscSF *sf0) 90 { 91 PetscInt j, k, r; 92 PetscInt nLeaves = ctx->nLeavesPerRank * ctx->size; 93 PetscInt nroots = ctx->nLeavesPerRank; 94 PetscSF sf; 95 PetscInt *ilocal; 96 PetscSFNode *iremote; 97 PetscErrorCode ierr; 98 99 PetscFunctionBegin; 100 ierr = PetscMalloc1(nLeaves+1, &ctx->ilocal);CHKERRQ(ierr); 101 ierr = PetscMalloc1(nLeaves, &ctx->iremote);CHKERRQ(ierr); 102 ilocal = ctx->ilocal; 103 iremote = ctx->iremote; 104 ilocal[nLeaves] = -ctx->leaveStep; 105 ierr = PetscSFCreate(ctx->comm, &sf);CHKERRQ(ierr); 106 for (r=0, j=nLeaves-1; r<ctx->size; r++) { 107 for (k=0; k<ctx->nLeavesPerRank; k++, j--) { 108 ilocal[j] = ilocal[j+1] + ctx->leaveStep; 109 iremote[j].rank = r; 110 iremote[j].index = k; 111 } 112 } 113 ierr = PetscSFSetGraph(sf, nroots, nLeaves, ilocal, ctx->localmode, iremote, ctx->remotemode);CHKERRQ(ierr); 114 { 115 const PetscInt *tlocal; 116 PetscBool sorted; 117 118 ierr = PetscSFGetGraph(sf, NULL, NULL, &tlocal, NULL);CHKERRQ(ierr); 119 PetscCheckFalse(ctx->contiguousLeaves && tlocal,PETSC_COMM_SELF,PETSC_ERR_PLIB,"ilocal=NULL expected for contiguous case"); 120 if (tlocal) { 121 ierr = PetscSortedInt(nLeaves, tlocal, &sorted);CHKERRQ(ierr); 122 PetscCheck(sorted,PETSC_COMM_SELF,PETSC_ERR_PLIB,"ilocal expected to be sorted"); 123 } 124 } 125 *sf0 = sf; 126 PetscFunctionReturn(0); 127 } 128 129 PetscErrorCode CreateSF1(AppCtx *ctx, PetscSF *sf1) 130 { 131 PetscInt j, k, r; 132 PetscInt *ilocal = NULL; 133 PetscSFNode *iremote; 134 PetscInt nLeaves = ctx->nLeavesPerRank * ctx->size; 135 PetscInt nroots = ctx->nLeavesPerRank; 136 PetscSF sf; 137 PetscErrorCode ierr; 138 139 PetscFunctionBegin; 140 ilocal = NULL; 141 if (!ctx->contiguousLeaves) { 142 ierr = PetscCalloc1(nLeaves+1, &ilocal);CHKERRQ(ierr); 143 } 144 ierr = PetscMalloc1(nLeaves, &iremote);CHKERRQ(ierr); 145 ierr = PetscSFCreate(ctx->comm, &sf);CHKERRQ(ierr); 146 for (r=0, j=0; r<ctx->size; r++) { 147 for (k=0; k<ctx->nLeavesPerRank; k++, j++) { 148 if (!ctx->contiguousLeaves) { 149 ilocal[j+1] = ilocal[j] + ctx->leaveStep; 150 } 151 iremote[j].rank = r; 152 iremote[j].index = k; 153 } 154 } 155 PetscCheck(j == nLeaves,PETSC_COMM_SELF,PETSC_ERR_PLIB,"j != nLeaves"); 156 ierr = PetscSFSetGraph(sf, nroots, nLeaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER);CHKERRQ(ierr); 157 if (ctx->contiguousLeaves) { 158 const PetscInt *tlocal; 159 160 ierr = PetscSFGetGraph(sf, NULL, NULL, &tlocal, NULL);CHKERRQ(ierr); 161 PetscCheckFalse(tlocal,PETSC_COMM_SELF,PETSC_ERR_PLIB,"ilocal=NULL expected for contiguous case"); 162 } 163 *sf1 = sf; 164 PetscFunctionReturn(0); 165 } 166 167 int main(int argc, char **argv) 168 { 169 AppCtx ctx; 170 PetscSF sf0, sf1; 171 MPI_Comm comm; 172 PetscErrorCode ierr; 173 174 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 175 comm = PETSC_COMM_WORLD; 176 ierr = GetOptions(comm, &ctx);CHKERRQ(ierr); 177 178 ierr = CreateSF0(&ctx, &sf0);CHKERRQ(ierr); 179 ierr = CreateSF1(&ctx, &sf1);CHKERRQ(ierr); 180 ierr = PetscSFViewFromOptions(sf0, NULL, "-sf0_view"); 181 ierr = PetscSFViewFromOptions(sf1, NULL, "-sf1_view"); 182 ierr = PetscSFCheckEqual_Private(sf0, sf1);CHKERRQ(ierr); 183 184 if (ctx.localmode != PETSC_OWN_POINTER) { 185 ierr = PetscFree(ctx.ilocal);CHKERRQ(ierr); 186 } 187 if (ctx.remotemode != PETSC_OWN_POINTER) { 188 ierr = PetscFree(ctx.iremote);CHKERRQ(ierr); 189 } 190 ierr = PetscSFDestroy(&sf0);CHKERRQ(ierr); 191 ierr = PetscSFDestroy(&sf1);CHKERRQ(ierr); 192 ierr = PetscFinalize(); 193 return ierr; 194 } 195 196 /*TEST 197 testset: 198 suffix: 1 199 nsize: {{1 3}} 200 args: -n_leaves_per_rank {{0 5}} -leave_step {{1 3}} 201 test: 202 suffix: a 203 args: -localmode {{COPY_VALUES OWN_POINTER}} -remotemode {{COPY_VALUES OWN_POINTER}} 204 test: 205 suffix: b 206 args: -localmode USE_POINTER -remotemode {{COPY_VALUES OWN_POINTER USE_POINTER}} 207 TEST*/ 208