1*d279a5e3SJunchao Zhang static const char help[] = "Test PetscSF with derived data types created with MPI large count\n\n"; 2*d279a5e3SJunchao Zhang 3*d279a5e3SJunchao Zhang #include <petscsys.h> 4*d279a5e3SJunchao Zhang #include <petscsf.h> 5*d279a5e3SJunchao Zhang 6*d279a5e3SJunchao Zhang int main(int argc, char **argv) 7*d279a5e3SJunchao Zhang { 8*d279a5e3SJunchao Zhang PetscSF sf; 9*d279a5e3SJunchao Zhang PetscInt i, nroots, nleaves; 10*d279a5e3SJunchao Zhang const PetscInt m = 4, n = 64; 11*d279a5e3SJunchao Zhang PetscSFNode *iremote = NULL; 12*d279a5e3SJunchao Zhang PetscMPIInt rank, size; 13*d279a5e3SJunchao Zhang int *rootdata = NULL, *leafdata = NULL; 14*d279a5e3SJunchao Zhang MPI_Datatype newtype; 15*d279a5e3SJunchao Zhang 16*d279a5e3SJunchao Zhang PetscFunctionBeginUser; 17*d279a5e3SJunchao Zhang PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 18*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 19*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 20*d279a5e3SJunchao Zhang PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "The test can only run with two MPI ranks"); 21*d279a5e3SJunchao Zhang 22*d279a5e3SJunchao Zhang PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sf)); 23*d279a5e3SJunchao Zhang PetscCall(PetscSFSetFromOptions(sf)); 24*d279a5e3SJunchao Zhang 25*d279a5e3SJunchao Zhang if (rank == 0) { 26*d279a5e3SJunchao Zhang nroots = n; 27*d279a5e3SJunchao Zhang nleaves = 0; 28*d279a5e3SJunchao Zhang } else { 29*d279a5e3SJunchao Zhang nroots = 0; 30*d279a5e3SJunchao Zhang nleaves = n; 31*d279a5e3SJunchao Zhang PetscCall(PetscMalloc1(nleaves, &iremote)); 32*d279a5e3SJunchao Zhang for (i = 0; i < nleaves; i++) { 33*d279a5e3SJunchao Zhang iremote[i].rank = 0; 34*d279a5e3SJunchao Zhang iremote[i].index = i; 35*d279a5e3SJunchao Zhang } 36*d279a5e3SJunchao Zhang } 37*d279a5e3SJunchao Zhang PetscCall(PetscSFSetGraph(sf, nroots, nleaves, NULL, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER)); 38*d279a5e3SJunchao Zhang 39*d279a5e3SJunchao Zhang PetscCall(PetscCalloc2(nroots * m, &rootdata, nleaves * m, &leafdata)); // allocate fat nodes to apply a derived data type of m MPI_INTs 40*d279a5e3SJunchao Zhang 41*d279a5e3SJunchao Zhang if (rank == 0) rootdata[nroots * m - 1] = 123; // set the last integer in rootdata and then check on leaves 42*d279a5e3SJunchao Zhang 43*d279a5e3SJunchao Zhang #if defined(PETSC_HAVE_MPI_LARGE_COUNT) 44*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Type_contiguous_c(m, MPI_INT, &newtype)); 45*d279a5e3SJunchao Zhang #else 46*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Type_contiguous(m, MPI_INT, &newtype)); 47*d279a5e3SJunchao Zhang #endif 48*d279a5e3SJunchao Zhang 49*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Type_commit(&newtype)); 50*d279a5e3SJunchao Zhang 51*d279a5e3SJunchao Zhang PetscCall(PetscSFBcastBegin(sf, newtype, rootdata, leafdata, MPI_REPLACE)); // bcast rootdata to leafdata 52*d279a5e3SJunchao Zhang PetscCall(PetscSFBcastEnd(sf, newtype, rootdata, leafdata, MPI_REPLACE)); 53*d279a5e3SJunchao Zhang 54*d279a5e3SJunchao Zhang if (rank == 1) PetscCheck(leafdata[nleaves * m - 1] == 123, PETSC_COMM_SELF, PETSC_ERR_PLIB, "SF: wrong results"); 55*d279a5e3SJunchao Zhang 56*d279a5e3SJunchao Zhang PetscCallMPI(MPI_Type_free(&newtype)); 57*d279a5e3SJunchao Zhang PetscCall(PetscFree2(rootdata, leafdata)); 58*d279a5e3SJunchao Zhang PetscCall(PetscSFDestroy(&sf)); 59*d279a5e3SJunchao Zhang PetscCall(PetscFinalize()); 60*d279a5e3SJunchao Zhang return 0; 61*d279a5e3SJunchao Zhang } 62*d279a5e3SJunchao Zhang 63*d279a5e3SJunchao Zhang /**TEST 64*d279a5e3SJunchao Zhang test: 65*d279a5e3SJunchao Zhang nsize: 2 66*d279a5e3SJunchao Zhang output_file: output/empty.out 67*d279a5e3SJunchao Zhang args: -sf_type {{basic neighbor}} 68*d279a5e3SJunchao Zhang 69*d279a5e3SJunchao Zhang TEST**/ 70