1*8b5d2d36SJed Brown static char help[] = "Test PetscSFSetGraphFromCoordinates()\n\n"; 2*8b5d2d36SJed Brown 3*8b5d2d36SJed Brown #include <petscsf.h> 4*8b5d2d36SJed Brown 5*8b5d2d36SJed Brown int main(int argc, char **argv) 6*8b5d2d36SJed Brown { 7*8b5d2d36SJed Brown PetscSF sf; 8*8b5d2d36SJed Brown MPI_Comm comm; 9*8b5d2d36SJed Brown PetscMPIInt rank, size; 10*8b5d2d36SJed Brown PetscInt height = 2, width = 3, nroots = height, nleaves, dim = 2; 11*8b5d2d36SJed Brown PetscReal *rootcoords, *leafcoords; 12*8b5d2d36SJed Brown PetscViewer viewer; 13*8b5d2d36SJed Brown 14*8b5d2d36SJed Brown PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 15*8b5d2d36SJed Brown comm = PETSC_COMM_WORLD; 16*8b5d2d36SJed Brown PetscCallMPI(MPI_Comm_rank(comm, &rank)); 17*8b5d2d36SJed Brown PetscCallMPI(MPI_Comm_size(comm, &size)); 18*8b5d2d36SJed Brown 19*8b5d2d36SJed Brown nleaves = (width - (rank == 0) - (rank == size - 1)) * height; 20*8b5d2d36SJed Brown PetscCall(PetscMalloc2(nroots * dim, &rootcoords, nleaves * dim, &leafcoords)); 21*8b5d2d36SJed Brown for (PetscInt i = 0; i < height; i++) { 22*8b5d2d36SJed Brown rootcoords[i * dim + 0] = 0.1 * rank; 23*8b5d2d36SJed Brown rootcoords[i * dim + 1] = 1. * i; 24*8b5d2d36SJed Brown for (PetscInt j = 0, l = 0; j < width; j++) { 25*8b5d2d36SJed Brown if (rank + j - 1 < 0 || rank + j - 1 >= size) continue; 26*8b5d2d36SJed Brown leafcoords[(i * nleaves / height + l) * dim + 0] = 0.1 * (rank + j - 1); 27*8b5d2d36SJed Brown leafcoords[(i * nleaves / height + l) * dim + 1] = 1. * i; 28*8b5d2d36SJed Brown l++; 29*8b5d2d36SJed Brown } 30*8b5d2d36SJed Brown } 31*8b5d2d36SJed Brown viewer = PETSC_VIEWER_STDOUT_WORLD; 32*8b5d2d36SJed Brown PetscCall(PetscPrintf(comm, "Roots by rank\n")); 33*8b5d2d36SJed Brown PetscCall(PetscRealView(nroots * dim, rootcoords, viewer)); 34*8b5d2d36SJed Brown PetscCall(PetscPrintf(comm, "Leaves by rank\n")); 35*8b5d2d36SJed Brown PetscCall(PetscRealView(nleaves * dim, leafcoords, viewer)); 36*8b5d2d36SJed Brown 37*8b5d2d36SJed Brown PetscCall(PetscSFCreate(comm, &sf)); 38*8b5d2d36SJed Brown PetscCall(PetscSFSetGraphFromCoordinates(sf, nroots, nleaves, dim, 1e-10, rootcoords, leafcoords)); 39*8b5d2d36SJed Brown 40*8b5d2d36SJed Brown PetscCall(PetscSFViewFromOptions(sf, NULL, "-sf_view")); 41*8b5d2d36SJed Brown PetscCall(PetscFree2(rootcoords, leafcoords)); 42*8b5d2d36SJed Brown PetscCall(PetscSFDestroy(&sf)); 43*8b5d2d36SJed Brown PetscCall(PetscFinalize()); 44*8b5d2d36SJed Brown return 0; 45*8b5d2d36SJed Brown } 46*8b5d2d36SJed Brown 47*8b5d2d36SJed Brown /*TEST 48*8b5d2d36SJed Brown test: 49*8b5d2d36SJed Brown suffix: 1 50*8b5d2d36SJed Brown nsize: 3 51*8b5d2d36SJed Brown args: -sf_view 52*8b5d2d36SJed Brown TEST*/ 53