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