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