1 static const char help[] = "Test overlapped communication on a single star forest (PetscSF)\n\n"; 2 3 #include <petscvec.h> 4 #include <petscsf.h> 5 #include <petscviewer.h> 6 7 int main(int argc, char **argv) 8 { 9 PetscInt ierr; 10 PetscSF sf; 11 Vec A,Aout; 12 Vec B,Bout; 13 PetscScalar *bufA; 14 PetscScalar *bufAout; 15 PetscScalar *bufB; 16 PetscScalar *bufBout; 17 PetscMPIInt rank, size; 18 PetscInt nroots, nleaves; 19 PetscInt i; 20 PetscInt *ilocal; 21 PetscSFNode *iremote; 22 23 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 24 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 25 CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 26 27 PetscCheckFalse(size != 2,PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for two MPI processes"); 28 29 CHKERRQ(PetscSFCreate(PETSC_COMM_WORLD,&sf)); 30 CHKERRQ(PetscSFSetFromOptions(sf)); 31 32 nleaves = 2; 33 nroots = 1; 34 CHKERRQ(PetscMalloc1(nleaves,&ilocal)); 35 36 for (i = 0; i<nleaves; i++) { 37 ilocal[i] = i; 38 } 39 40 CHKERRQ(PetscMalloc1(nleaves,&iremote)); 41 if (rank == 0) { 42 iremote[0].rank = 0; 43 iremote[0].index = 0; 44 iremote[1].rank = 1; 45 iremote[1].index = 0; 46 } else { 47 iremote[0].rank = 1; 48 iremote[0].index = 0; 49 iremote[1].rank = 0; 50 iremote[1].index = 0; 51 } 52 CHKERRQ(PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_OWN_POINTER,iremote,PETSC_OWN_POINTER)); 53 CHKERRQ(PetscSFSetUp(sf)); 54 CHKERRQ(PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD)); 55 CHKERRQ(VecCreate(PETSC_COMM_WORLD,&A)); 56 CHKERRQ(VecSetSizes(A,2,PETSC_DETERMINE)); 57 CHKERRQ(VecSetFromOptions(A)); 58 CHKERRQ(VecSetUp(A)); 59 60 CHKERRQ(VecDuplicate(A,&B)); 61 CHKERRQ(VecDuplicate(A,&Aout)); 62 CHKERRQ(VecDuplicate(A,&Bout)); 63 CHKERRQ(VecGetArray(A,&bufA)); 64 CHKERRQ(VecGetArray(B,&bufB)); 65 for (i=0; i<2; i++) { 66 bufA[i] = (PetscScalar)rank; 67 bufB[i] = (PetscScalar)(rank) + 10.0; 68 } 69 CHKERRQ(VecRestoreArray(A,&bufA)); 70 CHKERRQ(VecRestoreArray(B,&bufB)); 71 72 CHKERRQ(VecGetArrayRead(A,(const PetscScalar**)&bufA)); 73 CHKERRQ(VecGetArrayRead(B,(const PetscScalar**)&bufB)); 74 CHKERRQ(VecGetArray(Aout,&bufAout)); 75 CHKERRQ(VecGetArray(Bout,&bufBout)); 76 CHKERRQ(PetscSFBcastBegin(sf,MPIU_SCALAR,(const void*)bufA,(void *)bufAout,MPI_REPLACE)); 77 CHKERRQ(PetscSFBcastBegin(sf,MPIU_SCALAR,(const void*)bufB,(void *)bufBout,MPI_REPLACE)); 78 CHKERRQ(PetscSFBcastEnd(sf,MPIU_SCALAR,(const void*)bufA,(void *)bufAout,MPI_REPLACE)); 79 CHKERRQ(PetscSFBcastEnd(sf,MPIU_SCALAR,(const void*)bufB,(void *)bufBout,MPI_REPLACE)); 80 CHKERRQ(VecRestoreArrayRead(A,(const PetscScalar**)&bufA)); 81 CHKERRQ(VecRestoreArrayRead(B,(const PetscScalar**)&bufB)); 82 CHKERRQ(VecRestoreArray(Aout,&bufAout)); 83 CHKERRQ(VecRestoreArray(Bout,&bufBout)); 84 85 CHKERRQ(VecView(Aout,PETSC_VIEWER_STDOUT_WORLD)); 86 CHKERRQ(VecView(Bout,PETSC_VIEWER_STDOUT_WORLD)); 87 CHKERRQ(VecDestroy(&A)); 88 CHKERRQ(VecDestroy(&B)); 89 CHKERRQ(VecDestroy(&Aout)); 90 CHKERRQ(VecDestroy(&Bout)); 91 CHKERRQ(PetscSFDestroy(&sf)); 92 93 ierr = PetscFinalize(); 94 return ierr; 95 } 96 97 /*TEST 98 99 test: 100 suffix: basic 101 nsize: 2 102 filter: grep -v "type" | grep -v "sort" 103 args: -sf_type basic 104 105 test: 106 suffix: window 107 nsize: 2 108 filter: grep -v "type" | grep -v "sort" 109 output_file: output/ex2_basic.out 110 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create dynamic allocate}} 111 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW) 112 113 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes 114 test: 115 suffix: window_shared 116 nsize: 2 117 filter: grep -v "type" | grep -v "sort" 118 output_file: output/ex2_basic.out 119 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared 120 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH_NUMVERSION) defined(PETSC_HAVE_MPI_ONE_SIDED) 121 122 TEST*/ 123