1 static const char help[] = "Test embedded sf with one leaf data item connnected to multiple roots\n\n"; 2 3 #include <petscsf.h> 4 5 int main(int argc,char **argv) 6 { 7 PetscSF sf,newsf; 8 PetscInt i,nroots,nleaves,ilocal[2],leafdata,rootdata[2],nselected,selected,errors=0; 9 PetscSFNode iremote[2]; 10 PetscMPIInt myrank,next,nproc; 11 PetscErrorCode ierr; 12 MPI_Info info; 13 14 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 15 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&nproc);CHKERRQ(ierr); 16 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&myrank);CHKERRQ(ierr); 17 18 /* Create an SF that each process has two roots and two leaves. The two leaves are connected to the two 19 roots on next process. The special thing is each process only has one data item in its leaf data buffer. 20 */ 21 nroots = 2; 22 nleaves = 2; 23 ilocal[0] = 0; /* One leaf data item serves two leaves */ 24 ilocal[1] = 0; 25 next = (myrank+1)%nproc; 26 for (i=0; i<nleaves; i++) { 27 iremote[i].rank = next; 28 iremote[i].index = i; 29 } 30 31 ierr = PetscSFCreate(PETSC_COMM_WORLD,&sf);CHKERRQ(ierr); 32 ierr = PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_COPY_VALUES,iremote,PETSC_COPY_VALUES);CHKERRQ(ierr); 33 ierr = PetscSFSetFromOptions(sf);CHKERRQ(ierr); 34 35 /* Test PetscSFWindowSetInfo */ 36 ierr = MPI_Info_create(&info);CHKERRQ(ierr); 37 ierr = MPI_Info_set(info,"petsc_test_info1","info1");CHKERRQ(ierr); 38 ierr = MPI_Info_set(info,"petsc_test_info2","info2");CHKERRQ(ierr); 39 ierr = PetscSFWindowSetInfo(sf,info);CHKERRQ(ierr); 40 ierr = MPI_Info_free(&info);CHKERRQ(ierr); 41 42 ierr = PetscSFSetUp(sf);CHKERRQ(ierr); 43 ierr = PetscSFViewFromOptions(sf,NULL,"-sf_view");CHKERRQ(ierr); 44 45 /* Create an embedded sf by only selecting the first root on each process */ 46 nselected = 1; 47 selected = 0; 48 ierr = PetscSFCreateEmbeddedSF(sf,nselected,&selected,&newsf);CHKERRQ(ierr); 49 ierr = PetscSFViewFromOptions(newsf,NULL,"-esf_view");CHKERRQ(ierr); 50 51 /* Do reduce */ 52 leafdata = 1; 53 rootdata[0] = rootdata[1] = 0; 54 55 ierr = PetscSFReduceBegin(newsf,MPIU_INT,&leafdata,rootdata,MPIU_REPLACE);CHKERRQ(ierr); 56 ierr = PetscSFReduceEnd(newsf,MPIU_INT,&leafdata,rootdata,MPIU_REPLACE);CHKERRQ(ierr); 57 58 /* Check rootdata */ 59 if (rootdata[0] != 1 || rootdata[1] != 0) errors = 1; 60 ierr = MPI_Allreduce(MPI_IN_PLACE,&errors,1,MPIU_INT,MPI_SUM,PETSC_COMM_WORLD);CHKERRQ(ierr); 61 if (errors) {ierr = PetscPrintf(PETSC_COMM_WORLD,"Error: Unexpected rootdata on processors\n");CHKERRQ(ierr);} 62 63 ierr = PetscSFDestroy(&sf);CHKERRQ(ierr); 64 ierr = PetscSFDestroy(&newsf);CHKERRQ(ierr); 65 ierr = PetscFinalize(); 66 return ierr; 67 } 68 69 /*TEST 70 71 build: 72 requires: mpi 73 74 test: 75 nsize: 2 76 77 test: 78 suffix: 1_window 79 nsize: 2 80 filter: grep -v "type" | grep -v "sort" 81 output_file: output/ex2_1_window.out 82 args: -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor {{create dynamic allocate}} -sf_view ::ascii_info_detail -esf_view ::ascii_info_detail 83 requires: define(PETSC_HAVE_MPI_ONE_SIDED) 84 85 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes 86 test: 87 suffix: 1_window_shared 88 nsize: 2 89 output_file: output/ex2_1.out 90 args: -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor shared 91 requires: define(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !define(PETSC_HAVE_MPICH_NUMVERSION) define(PETSC_HAVE_MPI_ONE_SIDED) 92 93 TEST*/ 94