1 static char help[] = "Test PetscSFCompose() and PetscSFCreateStridedSF() when the ilocal arrays are not identity nor dense\n\n"; 2 3 #include <petsc.h> 4 #include <petscsf.h> 5 6 static PetscErrorCode TestVector(PetscSF sf, const char *sfname) 7 { 8 PetscInt mr, ml; 9 MPI_Comm comm; 10 11 PetscFunctionBeginUser; 12 comm = PetscObjectComm((PetscObject)sf); 13 PetscCall(PetscSFGetGraph(sf, &mr, NULL, NULL, NULL)); 14 PetscCall(PetscSFGetLeafRange(sf, NULL, &ml)); 15 for (PetscInt bs = 1; bs < 3; bs++) { 16 for (PetscInt r = 0; r < 2; r++) { 17 for (PetscInt l = 0; l < 2; l++) { 18 PetscSF vsf; 19 PetscInt *rdata, *ldata, *rdatav, *ldatav; 20 PetscInt ldr = PETSC_DECIDE; 21 PetscInt ldl = PETSC_DECIDE; 22 PetscBool flg; 23 24 if (r == 1) ldr = mr; 25 if (r == 2) ldr = mr + 7; 26 if (l == 1) ldl = ml + 1; 27 if (l == 2) ldl = ml + 5; 28 29 PetscCall(PetscSFCreateStridedSF(sf, bs, ldr, ldl, &vsf)); 30 if (ldr == PETSC_DECIDE) ldr = mr; 31 if (ldl == PETSC_DECIDE) ldl = ml + 1; 32 33 PetscCall(PetscCalloc4(bs * ldr, &rdata, bs * ldl, &ldata, bs * ldr, &rdatav, bs * ldl, &ldatav)); 34 for (PetscInt i = 0; i < bs * ldr; i++) rdata[i] = i + 1; 35 36 for (PetscInt i = 0; i < bs; i++) { 37 PetscCall(PetscSFBcastBegin(sf, MPIU_INT, PetscSafePointerPlusOffset(rdata, i * ldr), PetscSafePointerPlusOffset(ldata, i * ldl), MPI_REPLACE)); 38 PetscCall(PetscSFBcastEnd(sf, MPIU_INT, PetscSafePointerPlusOffset(rdata, i * ldr), PetscSafePointerPlusOffset(ldata, i * ldl), MPI_REPLACE)); 39 } 40 PetscCall(PetscSFBcastBegin(vsf, MPIU_INT, rdata, ldatav, MPI_REPLACE)); 41 PetscCall(PetscSFBcastEnd(vsf, MPIU_INT, rdata, ldatav, MPI_REPLACE)); 42 PetscCall(PetscArraycmp(ldata, ldatav, bs * ldl, &flg)); 43 44 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &flg, 1, MPI_C_BOOL, MPI_LAND, comm)); 45 if (!flg) { 46 PetscCall(PetscPrintf(comm, "Error with Bcast on %s: block size %" PetscInt_FMT ", ldr %" PetscInt_FMT ", ldl %" PetscInt_FMT "\n", sfname, bs, ldr, ldl)); 47 PetscCall(PetscPrintf(comm, "Single SF\n")); 48 PetscCall(PetscIntView(bs * ldl, ldata, PETSC_VIEWER_STDOUT_(comm))); 49 PetscCall(PetscPrintf(comm, "Vector SF\n")); 50 PetscCall(PetscIntView(bs * ldl, ldatav, PETSC_VIEWER_STDOUT_(comm))); 51 } 52 PetscCall(PetscArrayzero(rdata, bs * ldr)); 53 PetscCall(PetscArrayzero(rdatav, bs * ldr)); 54 55 for (PetscInt i = 0; i < bs; i++) { 56 PetscCall(PetscSFReduceBegin(sf, MPIU_INT, PetscSafePointerPlusOffset(ldata, i * ldl), PetscSafePointerPlusOffset(rdata, i * ldr), MPI_SUM)); 57 PetscCall(PetscSFReduceEnd(sf, MPIU_INT, PetscSafePointerPlusOffset(ldata, i * ldl), PetscSafePointerPlusOffset(rdata, i * ldr), MPI_SUM)); 58 } 59 PetscCall(PetscSFReduceBegin(vsf, MPIU_INT, ldata, rdatav, MPI_SUM)); 60 PetscCall(PetscSFReduceEnd(vsf, MPIU_INT, ldata, rdatav, MPI_SUM)); 61 PetscCall(PetscArraycmp(rdata, rdatav, bs * ldr, &flg)); 62 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &flg, 1, MPI_C_BOOL, MPI_LAND, comm)); 63 if (!flg) { 64 PetscCall(PetscPrintf(comm, "Error with Reduce on %s: block size %" PetscInt_FMT ", ldr %" PetscInt_FMT ", ldl %" PetscInt_FMT "\n", sfname, bs, ldr, ldl)); 65 PetscCall(PetscPrintf(comm, "Single SF\n")); 66 PetscCall(PetscIntView(bs * ldr, rdata, PETSC_VIEWER_STDOUT_(comm))); 67 PetscCall(PetscPrintf(comm, "Vector SF\n")); 68 PetscCall(PetscIntView(bs * ldr, rdatav, PETSC_VIEWER_STDOUT_(comm))); 69 } 70 PetscCall(PetscFree4(rdata, ldata, rdatav, ldatav)); 71 PetscCall(PetscSFDestroy(&vsf)); 72 } 73 } 74 } 75 PetscFunctionReturn(PETSC_SUCCESS); 76 } 77 78 int main(int argc, char **argv) 79 { 80 PetscSF sfA, sfB, sfBA, sfAAm, sfBBm, sfAm, sfBm; 81 PetscInt nrootsA, nleavesA, nrootsB, nleavesB; 82 PetscInt *ilocalA, *ilocalB; 83 PetscSFNode *iremoteA, *iremoteB; 84 PetscMPIInt rank, size; 85 PetscInt i, m, n, k, nl = 2, mA, mB, nldataA, nldataB; 86 PetscInt *rdA, *rdB, *ldA, *ldB; 87 PetscBool inverse = PETSC_FALSE, test_vector = PETSC_TRUE; 88 89 PetscFunctionBeginUser; 90 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 91 PetscCall(PetscOptionsGetInt(NULL, NULL, "-nl", &nl, NULL)); 92 PetscCall(PetscOptionsGetBool(NULL, NULL, "-explicit_inverse", &inverse, NULL)); 93 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 94 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 95 96 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfA)); 97 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sfB)); 98 PetscCall(PetscSFSetFromOptions(sfA)); 99 PetscCall(PetscSFSetFromOptions(sfB)); 100 101 // disable vector tests with linux-misc-32bit and sftype window 102 #if (PETSC_SIZEOF_SIZE_T == 4) 103 { 104 PetscBool iswin; 105 106 PetscCall(PetscObjectTypeCompare((PetscObject)sfA, PETSCSFWINDOW, &iswin)); 107 if (iswin) test_vector = PETSC_FALSE; 108 } 109 #endif 110 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_vector", &test_vector, NULL)); 111 112 n = 4 * nl * size; 113 m = 2 * nl; 114 k = nl; 115 116 nldataA = rank == 0 ? n : 0; 117 nldataB = 3 * nl; 118 119 nrootsA = m; 120 nleavesA = rank == 0 ? size * m : 0; 121 nrootsB = rank == 0 ? n : 0; 122 nleavesB = k; 123 124 PetscCall(PetscMalloc1(nleavesA, &ilocalA)); 125 PetscCall(PetscMalloc1(nleavesA, &iremoteA)); 126 PetscCall(PetscMalloc1(nleavesB, &ilocalB)); 127 PetscCall(PetscMalloc1(nleavesB, &iremoteB)); 128 129 /* sf A bcast is equivalent to a sparse gather on process 0 130 process 0 receives data in the middle [nl,3*nl] of the leaf data array for A */ 131 for (i = 0; i < nleavesA; i++) { 132 iremoteA[i].rank = i / m; 133 iremoteA[i].index = i % m; 134 ilocalA[i] = nl + i / m * 4 * nl + i % m; 135 } 136 137 /* sf B bcast is equivalent to a sparse scatter from process 0 138 process 0 sends data from [nl,2*nl] of the leaf data array for A 139 each process receives, in reverse order, in the middle [nl,2*nl] of the leaf data array for B */ 140 for (i = 0; i < nleavesB; i++) { 141 iremoteB[i].rank = 0; 142 iremoteB[i].index = rank * 4 * nl + nl + i % m; 143 ilocalB[i] = 2 * nl - i - 1; 144 } 145 PetscCall(PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER)); 146 PetscCall(PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER)); 147 PetscCall(PetscSFSetUp(sfA)); 148 PetscCall(PetscSFSetUp(sfB)); 149 PetscCall(PetscObjectSetName((PetscObject)sfA, "sfA")); 150 PetscCall(PetscObjectSetName((PetscObject)sfB, "sfB")); 151 PetscCall(PetscSFViewFromOptions(sfA, NULL, "-view")); 152 PetscCall(PetscSFViewFromOptions(sfB, NULL, "-view")); 153 154 if (test_vector) PetscCall(TestVector(sfA, "sfA")); 155 if (test_vector) PetscCall(TestVector(sfB, "sfB")); 156 157 PetscCall(PetscSFGetLeafRange(sfA, NULL, &mA)); 158 PetscCall(PetscSFGetLeafRange(sfB, NULL, &mB)); 159 PetscCall(PetscMalloc2(nrootsA, &rdA, nldataA, &ldA)); 160 PetscCall(PetscMalloc2(nrootsB, &rdB, nldataB, &ldB)); 161 for (i = 0; i < nrootsA; i++) rdA[i] = m * rank + i; 162 for (i = 0; i < nldataA; i++) ldA[i] = -1; 163 for (i = 0; i < nldataB; i++) ldB[i] = -1; 164 165 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastB(BcastA)\n")); 166 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: root data\n")); 167 PetscCall(PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD)); 168 PetscCall(PetscSFBcastBegin(sfA, MPIU_INT, rdA, ldA, MPI_REPLACE)); 169 PetscCall(PetscSFBcastEnd(sfA, MPIU_INT, rdA, ldA, MPI_REPLACE)); 170 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "A: leaf data (all)\n")); 171 PetscCall(PetscIntView(nldataA, ldA, PETSC_VIEWER_STDOUT_WORLD)); 172 PetscCall(PetscSFBcastBegin(sfB, MPIU_INT, ldA, ldB, MPI_REPLACE)); 173 PetscCall(PetscSFBcastEnd(sfB, MPIU_INT, ldA, ldB, MPI_REPLACE)); 174 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "B: leaf data (all)\n")); 175 PetscCall(PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD)); 176 177 PetscCall(PetscSFCompose(sfA, sfB, &sfBA)); 178 PetscCall(PetscSFSetFromOptions(sfBA)); 179 PetscCall(PetscSFSetUp(sfBA)); 180 PetscCall(PetscObjectSetName((PetscObject)sfBA, "sfBA")); 181 PetscCall(PetscSFViewFromOptions(sfBA, NULL, "-view")); 182 if (test_vector) PetscCall(TestVector(sfBA, "sfBA")); 183 184 for (i = 0; i < nldataB; i++) ldB[i] = -1; 185 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BcastBA\n")); 186 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: root data\n")); 187 PetscCall(PetscIntView(nrootsA, rdA, PETSC_VIEWER_STDOUT_WORLD)); 188 PetscCall(PetscSFBcastBegin(sfBA, MPIU_INT, rdA, ldB, MPI_REPLACE)); 189 PetscCall(PetscSFBcastEnd(sfBA, MPIU_INT, rdA, ldB, MPI_REPLACE)); 190 PetscCall(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "BA: leaf data (all)\n")); 191 PetscCall(PetscIntView(nldataB, ldB, PETSC_VIEWER_STDOUT_WORLD)); 192 193 PetscCall(PetscSFCreateInverseSF(sfA, &sfAm)); 194 PetscCall(PetscSFSetFromOptions(sfAm)); 195 PetscCall(PetscObjectSetName((PetscObject)sfAm, "sfAm")); 196 PetscCall(PetscSFViewFromOptions(sfAm, NULL, "-view")); 197 if (test_vector) PetscCall(TestVector(sfAm, "sfAm")); 198 199 if (!inverse) { 200 PetscCall(PetscSFComposeInverse(sfA, sfA, &sfAAm)); 201 } else { 202 PetscCall(PetscSFCompose(sfA, sfAm, &sfAAm)); 203 } 204 PetscCall(PetscSFSetFromOptions(sfAAm)); 205 PetscCall(PetscSFSetUp(sfAAm)); 206 PetscCall(PetscObjectSetName((PetscObject)sfAAm, "sfAAm")); 207 PetscCall(PetscSFViewFromOptions(sfAAm, NULL, "-view")); 208 if (test_vector) PetscCall(TestVector(sfAAm, "sfAAm")); 209 210 PetscCall(PetscSFCreateInverseSF(sfB, &sfBm)); 211 PetscCall(PetscSFSetFromOptions(sfBm)); 212 PetscCall(PetscObjectSetName((PetscObject)sfBm, "sfBm")); 213 PetscCall(PetscSFViewFromOptions(sfBm, NULL, "-view")); 214 if (test_vector) PetscCall(TestVector(sfBm, "sfBm")); 215 216 if (!inverse) { 217 PetscCall(PetscSFComposeInverse(sfB, sfB, &sfBBm)); 218 } else { 219 PetscCall(PetscSFCompose(sfB, sfBm, &sfBBm)); 220 } 221 PetscCall(PetscSFSetFromOptions(sfBBm)); 222 PetscCall(PetscSFSetUp(sfBBm)); 223 PetscCall(PetscObjectSetName((PetscObject)sfBBm, "sfBBm")); 224 PetscCall(PetscSFViewFromOptions(sfBBm, NULL, "-view")); 225 if (test_vector) PetscCall(TestVector(sfBBm, "sfBBm")); 226 227 PetscCall(PetscFree2(rdA, ldA)); 228 PetscCall(PetscFree2(rdB, ldB)); 229 230 PetscCall(PetscSFDestroy(&sfA)); 231 PetscCall(PetscSFDestroy(&sfB)); 232 PetscCall(PetscSFDestroy(&sfBA)); 233 PetscCall(PetscSFDestroy(&sfAm)); 234 PetscCall(PetscSFDestroy(&sfBm)); 235 PetscCall(PetscSFDestroy(&sfAAm)); 236 PetscCall(PetscSFDestroy(&sfBBm)); 237 238 PetscCall(PetscFinalize()); 239 return 0; 240 } 241 242 /*TEST 243 244 test: 245 suffix: 1 246 args: -view -explicit_inverse {{0 1}} 247 248 test: 249 nsize: 7 250 filter: grep -v "type" | grep -v "sort" 251 suffix: 2 252 args: -view -nl 5 -explicit_inverse {{0 1}} 253 254 # 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 255 test: 256 TODO: frequent timeout with the CI job linux-hip-cmplx 257 nsize: 7 258 suffix: 2_window 259 filter: grep -v "type" | grep -v "sort" 260 output_file: output/ex5_2.out 261 args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor {{create allocate}} 262 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW) 263 264 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes 265 test: 266 TODO: frequent timeout with the CI job linux-hip-cmplx 267 nsize: 7 268 suffix: 2_window_shared 269 filter: grep -v "type" | grep -v "sort" 270 output_file: output/ex5_2.out 271 args: -view -nl 5 -explicit_inverse {{0 1}} -sf_type window -sf_window_sync {{fence lock active}} -sf_window_flavor shared 272 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH) defined(PETSC_HAVE_MPI_ONE_SIDED) 273 274 TEST*/ 275