Actual source code: ex1f.F90
2: ! Description: A star forest is a simple tree with one root and zero or more leaves.
3: ! Many common communication patterns can be expressed as updates of rootdata using leafdata and vice-versa.
4: ! This example creates a star forest, communicates values using the graph views the graph, then destroys it.
5: !
6: ! This is a copy of ex1.c but currently only tests the broadcast operation
8: program main
9: #include <petsc/finclude/petscvec.h>
10: use petscmpi ! or mpi or mpi_f08
11: use petscvec
12: implicit none
14: PetscErrorCode ierr
15: PetscInt i,nroots,nrootsalloc,nleaves,nleavesalloc,mine(6),stride
16: type(PetscSFNode) remote(6)
17: PetscMPIInt rank,size
18: PetscSF sf
19: PetscInt rootdata(6),leafdata(6)
21: ! used with PetscSFGetGraph()
22: type(PetscSFNode), pointer :: gremote(:)
23: PetscInt, pointer :: gmine(:)
24: PetscInt gnroots,gnleaves;
26: PetscCallA(PetscInitialize(ierr))
27: stride = 2
28: PetscCallMPIA(MPI_Comm_rank(PETSC_COMM_WORLD,rank,ierr))
29: PetscCallMPIA(MPI_Comm_size(PETSC_COMM_WORLD,size,ierr))
31: if (rank == 0) then
32: nroots = 3
33: else
34: nroots = 2
35: endif
36: nrootsalloc = nroots * stride;
37: if (rank > 0) then
38: nleaves = 3
39: else
40: nleaves = 2
41: endif
42: nleavesalloc = nleaves * stride
43: if (stride > 1) then
44: do i=1,nleaves
45: mine(i) = stride * (i-1)
46: enddo
47: endif
49: ! Left periodic neighbor
50: remote(1)%rank = modulo(rank+size-1,size)
51: remote(1)%index = 1 * stride
52: ! Right periodic neighbor
53: remote(2)%rank = modulo(rank+1,size)
54: remote(2)%index = 0 * stride
55: if (rank > 0) then ! All processes reference rank 0, index
56: remote(3)%rank = 0
57: remote(3)%index = 2 * stride
58: endif
60: ! Create a star forest for communication
61: PetscCallA(PetscSFCreate(PETSC_COMM_WORLD,sf,ierr))
62: PetscCallA(PetscSFSetFromOptions(sf,ierr))
63: PetscCallA(PetscSFSetGraph(sf,nrootsalloc,nleaves,mine,PETSC_COPY_VALUES,remote,PETSC_COPY_VALUES,ierr))
64: PetscCallA(PetscSFSetUp(sf,ierr))
66: ! View graph, mostly useful for debugging purposes.
67: PetscCallA(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL,ierr))
68: PetscCallA(PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD,ierr))
69: PetscCallA(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD,ierr))
71: ! Allocate space for send and receive buffers. This example communicates PetscInt, but other types, including
72: ! * user-defined structures, could also be used.
73: ! Set rootdata buffer to be broadcast
74: do i=1,nrootsalloc
75: rootdata(i) = -1
76: enddo
77: do i=1,nroots
78: rootdata(1 + (i-1)*stride) = 100*(rank+1) + i - 1
79: enddo
81: ! Initialize local buffer, these values are never used.
82: do i=1,nleavesalloc
83: leafdata(i) = -1
84: enddo
86: ! Broadcast entries from rootdata to leafdata. Computation or other communication can be performed between the begin and end calls.
87: PetscCallA(PetscSFBcastBegin(sf,MPIU_INTEGER,rootdata,leafdata,MPI_REPLACE,ierr))
88: PetscCallA(PetscSFBcastEnd(sf,MPIU_INTEGER,rootdata,leafdata,MPI_REPLACE,ierr))
89: PetscCallA(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Bcast Rootdata\n",ierr))
90: PetscCallA(PetscIntView(nrootsalloc,rootdata,PETSC_VIEWER_STDOUT_WORLD,ierr))
91: PetscCallA(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"## Bcast Leafdata\n",ierr))
92: PetscCallA(PetscIntView(nleavesalloc,leafdata,PETSC_VIEWER_STDOUT_WORLD,ierr))
94: PetscCallA(PetscSFGetGraph(sf,gnroots,gnleaves,gmine,gremote,ierr))
95: if (gnleaves .ne. nleaves) then; SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'nleaves returned from PetscSFGetGraph() does not match that set with PetscSFSetGraph()'); endif
96: do i=1,nleaves
97: if (gmine(i) .ne. mine(i)) then; SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'Root from PetscSFGetGraph() does not match that set with PetscSFSetGraph()'); endif
98: enddo
99: do i=1,nleaves
100: if (gremote(i)%index .ne. remote(i)%index) then; SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'Leaf from PetscSFGetGraph() does not match that set with PetscSFSetGraph()'); endif
101: enddo
103: deallocate(gremote)
104: ! Clean storage for star forest.
105: PetscCallA(PetscSFDestroy(sf,ierr))
107: ! Create a star forest with continuous leaves and hence no buffer
108: PetscCallA(PetscSFCreate(PETSC_COMM_WORLD,sf,ierr))
109: PetscCallA(PetscSFSetFromOptions(sf,ierr))
110: PetscCallA(PetscSFSetGraph(sf,nrootsalloc,nleaves,PETSC_NULL_INTEGER,PETSC_COPY_VALUES,remote,PETSC_COPY_VALUES,ierr))
111: PetscCallA(PetscSFSetUp(sf,ierr))
113: ! View graph, mostly useful for debugging purposes.
114: PetscCallA(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL,ierr))
115: PetscCallA(PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD,ierr))
116: PetscCallA(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD,ierr))
118: PetscCallA(PetscSFGetGraph(sf,gnroots,gnleaves,gmine,gremote,ierr))
119: if (loc(gmine) .ne. loc(PETSC_NULL_INTEGER)) then; SETERRA(PETSC_COMM_WORLD,PETSC_ERR_PLIB,'Leaves from PetscSFGetGraph() not null as expected'); endif
120: deallocate(gremote)
121: PetscCallA(PetscSFDestroy(sf,ierr))
122: PetscCallA(PetscFinalize(ierr))
123: end
125: !/*TEST
126: ! build:
127: ! requires: defined(PETSC_HAVE_FORTRAN_TYPE_STAR)
128: !
129: ! test:
130: ! nsize: 3
131: !
132: !TEST*/