1 static char help[]= "Test PetscSFFCompose when the ilocal arrays are not identity nor dense\n\n"; 2 3 #include <petsc.h> 4 #include <petscsf.h> 5 6 int main(int argc, char **argv) 7 { 8 PetscErrorCode ierr; 9 PetscSF sfA, sfB, sfBA, sfAAm, sfBBm, sfAm, sfBm; 10 PetscInt nrootsA, nleavesA, nrootsB, nleavesB; 11 PetscInt *ilocalA, *ilocalB; 12 PetscSFNode *iremoteA, *iremoteB; 13 PetscMPIInt rank,size; 14 PetscInt i,m,n,k,nl = 2,mA,mB,nldataA,nldataB; 15 PetscInt *rdA,*rdB,*ldA,*ldB; 16 PetscBool inverse = PETSC_FALSE; 17 18 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 19 ierr = PetscOptionsGetInt(NULL,NULL,"-nl",&nl,NULL);CHKERRQ(ierr); 20 ierr = PetscOptionsGetBool(NULL,NULL,"-explicit_inverse",&inverse,NULL);CHKERRQ(ierr); 21 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 22 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 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 n = 4*nl*size; 30 m = 2*nl; 31 k = nl; 32 33 nldataA = !rank ? n : 0; 34 nldataB = 3*nl; 35 36 nrootsA = m; 37 nleavesA = !rank ? size*m : 0; 38 nrootsB = !rank ? n : 0; 39 nleavesB = k; 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 /* sf A bcast is equivalent to a sparse gather on process 0 47 process 0 receives data in the middle [nl,3*nl] of the leaf data array for A */ 48 for (i = 0; i < nleavesA; i++) { 49 iremoteA[i].rank = i/m; 50 iremoteA[i].index = i%m; 51 ilocalA[i] = nl + i/m * 4*nl + i%m; 52 } 53 54 /* sf B bcast is equivalent to a sparse scatter from process 0 55 process 0 sends data from [nl,2*nl] of the leaf data array for A 56 each process receives, in reverse order, in the middle [nl,2*nl] of the leaf data array for B */ 57 for (i = 0; i < nleavesB; i++) { 58 iremoteB[i].rank = 0; 59 iremoteB[i].index = rank * 4*nl + nl + i%m; 60 ilocalB[i] = 2*nl - i - 1; 61 } 62 ierr = PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER);CHKERRQ(ierr); 63 ierr = PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);CHKERRQ(ierr); 64 ierr = PetscSFSetUp(sfA);CHKERRQ(ierr); 65 ierr = PetscSFSetUp(sfB);CHKERRQ(ierr); 66 ierr = PetscObjectSetName((PetscObject)sfA, "sfA");CHKERRQ(ierr); 67 ierr = PetscObjectSetName((PetscObject)sfB, "sfB");CHKERRQ(ierr); 68 ierr = PetscSFViewFromOptions(sfA, NULL, "-view");CHKERRQ(ierr); 69 ierr = PetscSFViewFromOptions(sfB, NULL, "-view");CHKERRQ(ierr); 70 71 ierr = PetscSFGetLeafRange(sfA, NULL, &mA);CHKERRQ(ierr); 72 ierr = PetscSFGetLeafRange(sfB, NULL, &mB);CHKERRQ(ierr); 73 ierr = PetscMalloc2(nrootsA, &rdA, nldataA, &ldA);CHKERRQ(ierr); 74 ierr = PetscMalloc2(nrootsB, &rdB, nldataB, &ldB);CHKERRQ(ierr); 75 for (i = 0; i < nrootsA; i++) rdA[i] = m*rank + i; 76 for (i = 0; i < nldataA; i++) ldA[i] = -1; 77 for (i = 0; i < nldataB; i++) ldB[i] = -1; 78 79 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastB(BcastA)\n");CHKERRQ(ierr); 80 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: root data\n");CHKERRQ(ierr); 81 ierr = PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 82 ierr = PetscSFBcastBegin(sfA, MPIU_INT, rdA, ldA);CHKERRQ(ierr); 83 ierr = PetscSFBcastEnd(sfA, MPIU_INT, rdA, ldA);CHKERRQ(ierr); 84 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: leaf data (all)\n");CHKERRQ(ierr); 85 ierr = PetscIntView(nldataA, ldA, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 86 ierr = PetscSFBcastBegin(sfB, MPIU_INT, ldA, ldB);CHKERRQ(ierr); 87 ierr = PetscSFBcastEnd(sfB, MPIU_INT, ldA, ldB);CHKERRQ(ierr); 88 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "B: leaf data (all)\n");CHKERRQ(ierr); 89 ierr = PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 90 91 ierr = PetscSFCompose(sfA, sfB, &sfBA);CHKERRQ(ierr); 92 ierr = PetscSFSetFromOptions(sfBA);CHKERRQ(ierr); 93 ierr = PetscSFSetUp(sfBA);CHKERRQ(ierr); 94 ierr = PetscObjectSetName((PetscObject)sfBA, "sfBA");CHKERRQ(ierr); 95 ierr = PetscSFViewFromOptions(sfBA, NULL, "-view");CHKERRQ(ierr); 96 97 for (i = 0; i < nldataB; i++) ldB[i] = -1; 98 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastBA\n");CHKERRQ(ierr); 99 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: root data\n");CHKERRQ(ierr); 100 ierr = PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 101 ierr = PetscSFBcastBegin(sfBA, MPIU_INT, rdA, ldB);CHKERRQ(ierr); 102 ierr = PetscSFBcastEnd(sfBA, MPIU_INT, rdA, ldB);CHKERRQ(ierr); 103 ierr = PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: leaf data (all)\n");CHKERRQ(ierr); 104 ierr = PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 105 106 ierr = PetscSFCreateInverseSF(sfA, &sfAm);CHKERRQ(ierr); 107 ierr = PetscSFSetFromOptions(sfAm);CHKERRQ(ierr); 108 ierr = PetscObjectSetName((PetscObject)sfAm, "sfAm");CHKERRQ(ierr); 109 ierr = PetscSFViewFromOptions(sfAm, NULL, "-view");CHKERRQ(ierr); 110 111 if (!inverse) { 112 ierr = PetscSFComposeInverse(sfA, sfA, &sfAAm);CHKERRQ(ierr); 113 } else { 114 ierr = PetscSFCompose(sfA, sfAm, &sfAAm);CHKERRQ(ierr); 115 } 116 ierr = PetscSFSetFromOptions(sfAAm);CHKERRQ(ierr); 117 ierr = PetscSFSetUp(sfAAm);CHKERRQ(ierr); 118 ierr = PetscObjectSetName((PetscObject)sfAAm, "sfAAm");CHKERRQ(ierr); 119 ierr = PetscSFViewFromOptions(sfAAm, NULL, "-view");CHKERRQ(ierr); 120 121 ierr = PetscSFCreateInverseSF(sfB, &sfBm);CHKERRQ(ierr); 122 ierr = PetscSFSetFromOptions(sfBm);CHKERRQ(ierr); 123 ierr = PetscObjectSetName((PetscObject)sfBm, "sfBm");CHKERRQ(ierr); 124 ierr = PetscSFViewFromOptions(sfBm, NULL, "-view");CHKERRQ(ierr); 125 126 if (!inverse) { 127 ierr = PetscSFComposeInverse(sfB, sfB, &sfBBm);CHKERRQ(ierr); 128 } else { 129 ierr = PetscSFCompose(sfB, sfBm, &sfBBm);CHKERRQ(ierr); 130 } 131 ierr = PetscSFSetFromOptions(sfBBm);CHKERRQ(ierr); 132 ierr = PetscSFSetUp(sfBBm);CHKERRQ(ierr); 133 ierr = PetscObjectSetName((PetscObject)sfBBm, "sfBBm");CHKERRQ(ierr); 134 ierr = PetscSFViewFromOptions(sfBBm, NULL, "-view");CHKERRQ(ierr); 135 136 ierr = PetscFree2(rdA, ldA);CHKERRQ(ierr); 137 ierr = PetscFree2(rdB, ldB);CHKERRQ(ierr); 138 139 ierr = PetscSFDestroy(&sfA);CHKERRQ(ierr); 140 ierr = PetscSFDestroy(&sfB);CHKERRQ(ierr); 141 ierr = PetscSFDestroy(&sfBA);CHKERRQ(ierr); 142 ierr = PetscSFDestroy(&sfAm);CHKERRQ(ierr); 143 ierr = PetscSFDestroy(&sfBm);CHKERRQ(ierr); 144 ierr = PetscSFDestroy(&sfAAm);CHKERRQ(ierr); 145 ierr = PetscSFDestroy(&sfBBm);CHKERRQ(ierr); 146 147 ierr = PetscFinalize(); 148 149 return ierr; 150 } 151 152 /*TEST 153 154 test: 155 suffix: 1 156 args: -view -explicit_inverse {{0 1}} 157 158 test: 159 nsize: 7 160 filter: grep -v "type" | grep -v "sort" 161 suffix: 2 162 args: -view -nl 5 -explicit_inverse {{0 1}} 163 164 # we cannot test for -sf_window_flavor dynamic because SFCompose with sparse leaves may change the root data pointer only locally, and this is not supported by the dynamic case 165 test: 166 nsize: 7 167 suffix: 2_window 168 filter: grep -v "type" | grep -v "sort" 169 output_file: output/ex5_2.out 170 args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor {{create allocate}} 171 requires: define(PETSC_HAVE_MPI_ONE_SIDED) 172 173 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes 174 test: 175 nsize: 7 176 suffix: 2_window_shared 177 filter: grep -v "type" | grep -v "sort" 178 output_file: output/ex5_2.out 179 args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor shared 180 requires: define(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !define(PETSC_HAVE_MPICH_NUMVERSION) define(PETSC_HAVE_MPI_ONE_SIDED) 181 182 TEST*/ 183