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