1 static char help[]= "Test PetscSFFCompose when the ilocal array is not the identity\n\n"; 2 3 #include <petscsf.h> 4 5 int main(int argc, char **argv) 6 { 7 PetscErrorCode ierr; 8 PetscSF sfA, sfB, sfBA; 9 PetscInt nrootsA, nleavesA, nrootsB, nleavesB; 10 PetscInt *ilocalA, *ilocalB; 11 PetscSFNode *iremoteA, *iremoteB; 12 Vec a, b, ba; 13 const PetscScalar *arrayR; 14 PetscScalar *arrayW; 15 PetscMPIInt size; 16 PetscInt i; 17 PetscInt maxleafB; 18 PetscBool flag = PETSC_FALSE; 19 20 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 21 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRMPI(ierr); 22 PetscCheck(size == 1,PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for one MPI process"); 23 24 ierr = PetscSFCreate(PETSC_COMM_WORLD, &sfA);CHKERRQ(ierr); 25 ierr = PetscSFCreate(PETSC_COMM_WORLD, &sfB);CHKERRQ(ierr); 26 ierr = PetscSFSetFromOptions(sfA);CHKERRQ(ierr); 27 ierr = PetscSFSetFromOptions(sfB);CHKERRQ(ierr); 28 29 ierr = PetscOptionsGetBool(NULL,NULL,"-sparse_sfB",&flag,NULL);CHKERRQ(ierr); 30 31 if (flag) { 32 /* sfA permutes indices, sfB has sparse leaf space. */ 33 nrootsA = 3; 34 nleavesA = 3; 35 nrootsB = 3; 36 nleavesB = 2; 37 } else { 38 /* sfA reverses indices, sfB is identity */ 39 nrootsA = nrootsB = nleavesA = nleavesB = 4; 40 } 41 ierr = PetscMalloc1(nleavesA, &ilocalA);CHKERRQ(ierr); 42 ierr = PetscMalloc1(nleavesA, &iremoteA);CHKERRQ(ierr); 43 ierr = PetscMalloc1(nleavesB, &ilocalB);CHKERRQ(ierr); 44 ierr = PetscMalloc1(nleavesB, &iremoteB);CHKERRQ(ierr); 45 46 for (i = 0; i < nleavesA; i++) { 47 iremoteA[i].rank = 0; 48 iremoteA[i].index = i; 49 if (flag) { 50 ilocalA[i] = (i + 1) % nleavesA; 51 } else { 52 ilocalA[i] = nleavesA - i - 1; 53 } 54 } 55 56 for (i = 0; i < nleavesB; i++) { 57 iremoteB[i].rank = 0; 58 if (flag) { 59 ilocalB[i] = nleavesB - i; 60 iremoteB[i].index = nleavesB - i - 1; 61 } else { 62 ilocalB[i] = i; 63 iremoteB[i].index = i; 64 } 65 } 66 67 ierr = PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER);CHKERRQ(ierr); 68 ierr = PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);CHKERRQ(ierr); 69 ierr = PetscSFSetUp(sfA);CHKERRQ(ierr); 70 ierr = PetscSFSetUp(sfB);CHKERRQ(ierr); 71 ierr = PetscObjectSetName((PetscObject)sfA, "sfA");CHKERRQ(ierr); 72 ierr = PetscObjectSetName((PetscObject)sfB, "sfB");CHKERRQ(ierr); 73 74 ierr = VecCreateSeq(PETSC_COMM_WORLD, nrootsA, &a);CHKERRQ(ierr); 75 ierr = VecCreateSeq(PETSC_COMM_WORLD, nleavesA, &b);CHKERRQ(ierr); 76 ierr = PetscSFGetLeafRange(sfB, NULL, &maxleafB);CHKERRQ(ierr); 77 ierr = VecCreateSeq(PETSC_COMM_WORLD, maxleafB+1, &ba);CHKERRQ(ierr); 78 ierr = VecGetArray(a, &arrayW);CHKERRQ(ierr); 79 for (i = 0; i < nrootsA; i++) { 80 arrayW[i] = (PetscScalar)i; 81 } 82 ierr = VecRestoreArray(a, &arrayW);CHKERRQ(ierr); 83 84 ierr = PetscPrintf(PETSC_COMM_WORLD, "Initial Vec A\n");CHKERRQ(ierr); 85 ierr = VecView(a, NULL);CHKERRQ(ierr); 86 ierr = VecGetArrayRead(a, &arrayR);CHKERRQ(ierr); 87 ierr = VecGetArray(b, &arrayW);CHKERRQ(ierr); 88 89 ierr = PetscSFBcastBegin(sfA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 90 ierr = PetscSFBcastEnd(sfA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 91 ierr = VecRestoreArray(b, &arrayW);CHKERRQ(ierr); 92 ierr = VecRestoreArrayRead(a, &arrayR);CHKERRQ(ierr); 93 ierr = PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->B over sfA\n");CHKERRQ(ierr); 94 ierr = VecView(b, NULL);CHKERRQ(ierr); 95 96 ierr = VecGetArrayRead(b, &arrayR);CHKERRQ(ierr); 97 ierr = VecGetArray(ba, &arrayW);CHKERRQ(ierr); 98 arrayW[0] = 10.0; /* Not touched by bcast */ 99 ierr = PetscSFBcastBegin(sfB, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 100 ierr = PetscSFBcastEnd(sfB, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 101 ierr = VecRestoreArrayRead(b, &arrayR);CHKERRQ(ierr); 102 ierr = VecRestoreArray(ba, &arrayW);CHKERRQ(ierr); 103 104 ierr = PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast B->BA over sfB\n");CHKERRQ(ierr); 105 ierr = VecView(ba, NULL);CHKERRQ(ierr); 106 107 ierr = PetscSFCompose(sfA, sfB, &sfBA);CHKERRQ(ierr); 108 ierr = PetscSFSetFromOptions(sfBA);CHKERRQ(ierr); 109 ierr = PetscObjectSetName((PetscObject)sfBA, "(sfB o sfA)");CHKERRQ(ierr); 110 ierr = VecGetArrayRead(a, &arrayR);CHKERRQ(ierr); 111 ierr = VecGetArray(ba, &arrayW);CHKERRQ(ierr); 112 arrayW[0] = 11.0; /* Not touched by bcast */ 113 ierr = PetscSFBcastBegin(sfBA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 114 ierr = PetscSFBcastEnd(sfBA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);CHKERRQ(ierr); 115 ierr = VecRestoreArray(ba, &arrayW);CHKERRQ(ierr); 116 ierr = VecRestoreArrayRead(a, &arrayR);CHKERRQ(ierr); 117 ierr = PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->BA over sfBA (sfB o sfA)\n");CHKERRQ(ierr); 118 ierr = VecView(ba, NULL);CHKERRQ(ierr); 119 120 ierr = VecDestroy(&ba);CHKERRQ(ierr); 121 ierr = VecDestroy(&b);CHKERRQ(ierr); 122 ierr = VecDestroy(&a);CHKERRQ(ierr); 123 124 ierr = PetscSFView(sfA, NULL);CHKERRQ(ierr); 125 ierr = PetscSFView(sfB, NULL);CHKERRQ(ierr); 126 ierr = PetscSFView(sfBA, NULL);CHKERRQ(ierr); 127 ierr = PetscSFDestroy(&sfA);CHKERRQ(ierr); 128 ierr = PetscSFDestroy(&sfB);CHKERRQ(ierr); 129 ierr = PetscSFDestroy(&sfBA);CHKERRQ(ierr); 130 131 ierr = PetscFinalize(); 132 return ierr; 133 } 134 135 /*TEST 136 137 test: 138 suffix: 1 139 140 test: 141 suffix: 2 142 filter: grep -v "type" | grep -v "sort" 143 args: -sparse_sfB 144 145 test: 146 suffix: 2_window 147 filter: grep -v "type" | grep -v "sort" 148 output_file: output/ex4_2.out 149 args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create dynamic allocate}} 150 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW) 151 152 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes 153 test: 154 suffix: 2_window_shared 155 filter: grep -v "type" | grep -v "sort" 156 output_file: output/ex4_2.out 157 args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared 158 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH_NUMVERSION) defined(PETSC_HAVE_MPI_ONE_SIDED) 159 160 TEST*/ 161