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