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