1 2 /* 3 The memory scalable AO application ordering routines. These store the 4 orderings on each processor for that processor's range of values 5 */ 6 7 #include <../src/vec/is/ao/aoimpl.h> /*I "petscao.h" I*/ 8 9 typedef struct { 10 PetscInt *app_loc; /* app_loc[i] is the partner for the ith local PETSc slot */ 11 PetscInt *petsc_loc; /* petsc_loc[j] is the partner for the jth local app slot */ 12 PetscLayout map; /* determines the local sizes of ao */ 13 } AO_MemoryScalable; 14 15 /* 16 All processors ship the data to process 0 to be printed; note that this is not scalable because 17 process 0 allocates space for all the orderings entry across all the processes 18 */ 19 PetscErrorCode AOView_MemoryScalable(AO ao, PetscViewer viewer) 20 { 21 PetscMPIInt rank, size; 22 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 23 PetscBool iascii; 24 PetscMPIInt tag_app, tag_petsc; 25 PetscLayout map = aomems->map; 26 PetscInt *app, *app_loc, *petsc, *petsc_loc, len, i, j; 27 MPI_Status status; 28 29 PetscFunctionBegin; 30 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii)); 31 PetscCheck(iascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Viewer type %s not supported for AO MemoryScalable", ((PetscObject)viewer)->type_name); 32 33 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)ao), &rank)); 34 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)ao), &size)); 35 36 PetscCall(PetscObjectGetNewTag((PetscObject)ao, &tag_app)); 37 PetscCall(PetscObjectGetNewTag((PetscObject)ao, &tag_petsc)); 38 39 if (rank == 0) { 40 PetscCall(PetscViewerASCIIPrintf(viewer, "Number of elements in ordering %" PetscInt_FMT "\n", ao->N)); 41 PetscCall(PetscViewerASCIIPrintf(viewer, "PETSc->App App->PETSc\n")); 42 43 PetscCall(PetscMalloc2(map->N, &app, map->N, &petsc)); 44 len = map->n; 45 /* print local AO */ 46 PetscCall(PetscViewerASCIIPrintf(viewer, "Process [%d]\n", rank)); 47 for (i = 0; i < len; i++) PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " %3" PetscInt_FMT " %3" PetscInt_FMT " %3" PetscInt_FMT "\n", i, aomems->app_loc[i], i, aomems->petsc_loc[i])); 48 49 /* recv and print off-processor's AO */ 50 for (i = 1; i < size; i++) { 51 len = map->range[i + 1] - map->range[i]; 52 app_loc = app + map->range[i]; 53 petsc_loc = petsc + map->range[i]; 54 PetscCallMPI(MPI_Recv(app_loc, (PetscMPIInt)len, MPIU_INT, i, tag_app, PetscObjectComm((PetscObject)ao), &status)); 55 PetscCallMPI(MPI_Recv(petsc_loc, (PetscMPIInt)len, MPIU_INT, i, tag_petsc, PetscObjectComm((PetscObject)ao), &status)); 56 PetscCall(PetscViewerASCIIPrintf(viewer, "Process [%" PetscInt_FMT "]\n", i)); 57 for (j = 0; j < len; j++) PetscCall(PetscViewerASCIIPrintf(viewer, "%3" PetscInt_FMT " %3" PetscInt_FMT " %3" PetscInt_FMT " %3" PetscInt_FMT "\n", map->range[i] + j, app_loc[j], map->range[i] + j, petsc_loc[j])); 58 } 59 PetscCall(PetscFree2(app, petsc)); 60 61 } else { 62 /* send values */ 63 PetscCallMPI(MPI_Send((void *)aomems->app_loc, map->n, MPIU_INT, 0, tag_app, PetscObjectComm((PetscObject)ao))); 64 PetscCallMPI(MPI_Send((void *)aomems->petsc_loc, map->n, MPIU_INT, 0, tag_petsc, PetscObjectComm((PetscObject)ao))); 65 } 66 PetscCall(PetscViewerFlush(viewer)); 67 PetscFunctionReturn(PETSC_SUCCESS); 68 } 69 70 PetscErrorCode AODestroy_MemoryScalable(AO ao) 71 { 72 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 73 74 PetscFunctionBegin; 75 PetscCall(PetscFree2(aomems->app_loc, aomems->petsc_loc)); 76 PetscCall(PetscLayoutDestroy(&aomems->map)); 77 PetscCall(PetscFree(aomems)); 78 PetscFunctionReturn(PETSC_SUCCESS); 79 } 80 81 /* 82 Input Parameters: 83 + ao - the application ordering context 84 . n - the number of integers in ia[] 85 . ia - the integers; these are replaced with their mapped value 86 - maploc - app_loc or petsc_loc in struct "AO_MemoryScalable" 87 88 Output Parameter: 89 . ia - the mapped interges 90 */ 91 PetscErrorCode AOMap_MemoryScalable_private(AO ao, PetscInt n, PetscInt *ia, const PetscInt *maploc) 92 { 93 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 94 MPI_Comm comm; 95 PetscMPIInt rank, size, tag1, tag2; 96 PetscInt *owner, *start, *sizes, nsends, nreceives; 97 PetscInt nmax, count, *sindices, *rindices, i, j, idx, lastidx, *sindices2, *rindices2; 98 const PetscInt *owners = aomems->map->range; 99 MPI_Request *send_waits, *recv_waits, *send_waits2, *recv_waits2; 100 MPI_Status recv_status; 101 PetscMPIInt nindices, source, widx; 102 PetscInt *rbuf, *sbuf; 103 MPI_Status *send_status, *send_status2; 104 105 PetscFunctionBegin; 106 PetscCall(PetscObjectGetComm((PetscObject)ao, &comm)); 107 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 108 PetscCallMPI(MPI_Comm_size(comm, &size)); 109 110 /* first count number of contributors to each processor */ 111 PetscCall(PetscMalloc1(size, &start)); 112 PetscCall(PetscCalloc2(2 * size, &sizes, n, &owner)); 113 114 j = 0; 115 lastidx = -1; 116 for (i = 0; i < n; i++) { 117 if (ia[i] < 0) owner[i] = -1; /* mark negative entries (which are not to be mapped) with a special negative value */ 118 if (ia[i] >= ao->N) owner[i] = -2; /* mark out of range entries with special negative value */ 119 else { 120 /* if indices are NOT locally sorted, need to start search at the beginning */ 121 if (lastidx > (idx = ia[i])) j = 0; 122 lastidx = idx; 123 for (; j < size; j++) { 124 if (idx >= owners[j] && idx < owners[j + 1]) { 125 sizes[2 * j]++; /* num of indices to be sent */ 126 sizes[2 * j + 1] = 1; /* send to proc[j] */ 127 owner[i] = j; 128 break; 129 } 130 } 131 } 132 } 133 sizes[2 * rank] = sizes[2 * rank + 1] = 0; /* do not receive from self! */ 134 nsends = 0; 135 for (i = 0; i < size; i++) nsends += sizes[2 * i + 1]; 136 137 /* inform other processors of number of messages and max length*/ 138 PetscCall(PetscMaxSum(comm, sizes, &nmax, &nreceives)); 139 140 /* allocate arrays */ 141 PetscCall(PetscObjectGetNewTag((PetscObject)ao, &tag1)); 142 PetscCall(PetscObjectGetNewTag((PetscObject)ao, &tag2)); 143 144 PetscCall(PetscMalloc2(nreceives * nmax, &rindices, nreceives, &recv_waits)); 145 PetscCall(PetscMalloc2(nsends * nmax, &rindices2, nsends, &recv_waits2)); 146 147 PetscCall(PetscMalloc3(n, &sindices, nsends, &send_waits, nsends, &send_status)); 148 PetscCall(PetscMalloc3(n, &sindices2, nreceives, &send_waits2, nreceives, &send_status2)); 149 150 /* post 1st receives: receive others requests 151 since we don't know how long each individual message is we 152 allocate the largest needed buffer for each receive. Potentially 153 this is a lot of wasted space. 154 */ 155 for (i = 0, count = 0; i < nreceives; i++) PetscCallMPI(MPI_Irecv(rindices + nmax * i, nmax, MPIU_INT, MPI_ANY_SOURCE, tag1, comm, recv_waits + count++)); 156 157 /* do 1st sends: 158 1) starts[i] gives the starting index in svalues for stuff going to 159 the ith processor 160 */ 161 start[0] = 0; 162 for (i = 1; i < size; i++) start[i] = start[i - 1] + sizes[2 * i - 2]; 163 for (i = 0; i < n; i++) { 164 j = owner[i]; 165 if (j == -1) continue; /* do not remap negative entries in ia[] */ 166 else if (j == -2) { /* out of range entries get mapped to -1 */ ia[i] = -1; 167 continue; 168 } else if (j != rank) { 169 sindices[start[j]++] = ia[i]; 170 } else { /* compute my own map */ 171 ia[i] = maploc[ia[i] - owners[rank]]; 172 } 173 } 174 175 start[0] = 0; 176 for (i = 1; i < size; i++) start[i] = start[i - 1] + sizes[2 * i - 2]; 177 for (i = 0, count = 0; i < size; i++) { 178 if (sizes[2 * i + 1]) { 179 /* send my request to others */ 180 PetscCallMPI(MPI_Isend(sindices + start[i], sizes[2 * i], MPIU_INT, i, tag1, comm, send_waits + count)); 181 /* post receive for the answer of my request */ 182 PetscCallMPI(MPI_Irecv(sindices2 + start[i], sizes[2 * i], MPIU_INT, i, tag2, comm, recv_waits2 + count)); 183 count++; 184 } 185 } 186 PetscCheck(nsends == count, comm, PETSC_ERR_SUP, "nsends %" PetscInt_FMT " != count %" PetscInt_FMT, nsends, count); 187 188 /* wait on 1st sends */ 189 if (nsends) PetscCallMPI(MPI_Waitall(nsends, send_waits, send_status)); 190 191 /* 1st recvs: other's requests */ 192 for (j = 0; j < nreceives; j++) { 193 PetscCallMPI(MPI_Waitany(nreceives, recv_waits, &widx, &recv_status)); /* idx: index of handle for operation that completed */ 194 PetscCallMPI(MPI_Get_count(&recv_status, MPIU_INT, &nindices)); 195 rbuf = rindices + nmax * widx; /* global index */ 196 source = recv_status.MPI_SOURCE; 197 198 /* compute mapping */ 199 sbuf = rbuf; 200 for (i = 0; i < nindices; i++) sbuf[i] = maploc[rbuf[i] - owners[rank]]; 201 202 /* send mapping back to the sender */ 203 PetscCallMPI(MPI_Isend(sbuf, nindices, MPIU_INT, source, tag2, comm, send_waits2 + widx)); 204 } 205 206 /* wait on 2nd sends */ 207 if (nreceives) PetscCallMPI(MPI_Waitall(nreceives, send_waits2, send_status2)); 208 209 /* 2nd recvs: for the answer of my request */ 210 for (j = 0; j < nsends; j++) { 211 PetscCallMPI(MPI_Waitany(nsends, recv_waits2, &widx, &recv_status)); 212 PetscCallMPI(MPI_Get_count(&recv_status, MPIU_INT, &nindices)); 213 source = recv_status.MPI_SOURCE; 214 /* pack output ia[] */ 215 rbuf = sindices2 + start[source]; 216 count = 0; 217 for (i = 0; i < n; i++) { 218 if (source == owner[i]) ia[i] = rbuf[count++]; 219 } 220 } 221 222 /* free arrays */ 223 PetscCall(PetscFree(start)); 224 PetscCall(PetscFree2(sizes, owner)); 225 PetscCall(PetscFree2(rindices, recv_waits)); 226 PetscCall(PetscFree2(rindices2, recv_waits2)); 227 PetscCall(PetscFree3(sindices, send_waits, send_status)); 228 PetscCall(PetscFree3(sindices2, send_waits2, send_status2)); 229 PetscFunctionReturn(PETSC_SUCCESS); 230 } 231 232 PetscErrorCode AOPetscToApplication_MemoryScalable(AO ao, PetscInt n, PetscInt *ia) 233 { 234 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 235 PetscInt *app_loc = aomems->app_loc; 236 237 PetscFunctionBegin; 238 PetscCall(AOMap_MemoryScalable_private(ao, n, ia, app_loc)); 239 PetscFunctionReturn(PETSC_SUCCESS); 240 } 241 242 PetscErrorCode AOApplicationToPetsc_MemoryScalable(AO ao, PetscInt n, PetscInt *ia) 243 { 244 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 245 PetscInt *petsc_loc = aomems->petsc_loc; 246 247 PetscFunctionBegin; 248 PetscCall(AOMap_MemoryScalable_private(ao, n, ia, petsc_loc)); 249 PetscFunctionReturn(PETSC_SUCCESS); 250 } 251 252 static struct _AOOps AOOps_MemoryScalable = { 253 PetscDesignatedInitializer(view, AOView_MemoryScalable), 254 PetscDesignatedInitializer(destroy, AODestroy_MemoryScalable), 255 PetscDesignatedInitializer(petsctoapplication, AOPetscToApplication_MemoryScalable), 256 PetscDesignatedInitializer(applicationtopetsc, AOApplicationToPetsc_MemoryScalable), 257 }; 258 259 PetscErrorCode AOCreateMemoryScalable_private(MPI_Comm comm, PetscInt napp, const PetscInt from_array[], const PetscInt to_array[], AO ao, PetscInt *aomap_loc) 260 { 261 AO_MemoryScalable *aomems = (AO_MemoryScalable *)ao->data; 262 PetscLayout map = aomems->map; 263 PetscInt n_local = map->n, i, j; 264 PetscMPIInt rank, size, tag; 265 PetscInt *owner, *start, *sizes, nsends, nreceives; 266 PetscInt nmax, count, *sindices, *rindices, idx, lastidx; 267 PetscInt *owners = aomems->map->range; 268 MPI_Request *send_waits, *recv_waits; 269 MPI_Status recv_status; 270 PetscMPIInt nindices, widx; 271 PetscInt *rbuf; 272 PetscInt n = napp, ip, ia; 273 MPI_Status *send_status; 274 275 PetscFunctionBegin; 276 PetscCall(PetscArrayzero(aomap_loc, n_local)); 277 278 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 279 PetscCallMPI(MPI_Comm_size(comm, &size)); 280 281 /* first count number of contributors (of from_array[]) to each processor */ 282 PetscCall(PetscCalloc1(2 * size, &sizes)); 283 PetscCall(PetscMalloc1(n, &owner)); 284 285 j = 0; 286 lastidx = -1; 287 for (i = 0; i < n; i++) { 288 /* if indices are NOT locally sorted, need to start search at the beginning */ 289 if (lastidx > (idx = from_array[i])) j = 0; 290 lastidx = idx; 291 for (; j < size; j++) { 292 if (idx >= owners[j] && idx < owners[j + 1]) { 293 sizes[2 * j] += 2; /* num of indices to be sent - in pairs (ip,ia) */ 294 sizes[2 * j + 1] = 1; /* send to proc[j] */ 295 owner[i] = j; 296 break; 297 } 298 } 299 } 300 sizes[2 * rank] = sizes[2 * rank + 1] = 0; /* do not receive from self! */ 301 nsends = 0; 302 for (i = 0; i < size; i++) nsends += sizes[2 * i + 1]; 303 304 /* inform other processors of number of messages and max length*/ 305 PetscCall(PetscMaxSum(comm, sizes, &nmax, &nreceives)); 306 307 /* allocate arrays */ 308 PetscCall(PetscObjectGetNewTag((PetscObject)ao, &tag)); 309 PetscCall(PetscMalloc2(nreceives * nmax, &rindices, nreceives, &recv_waits)); 310 PetscCall(PetscMalloc3(2 * n, &sindices, nsends, &send_waits, nsends, &send_status)); 311 PetscCall(PetscMalloc1(size, &start)); 312 313 /* post receives: */ 314 for (i = 0; i < nreceives; i++) PetscCallMPI(MPI_Irecv(rindices + nmax * i, nmax, MPIU_INT, MPI_ANY_SOURCE, tag, comm, recv_waits + i)); 315 316 /* do sends: 317 1) starts[i] gives the starting index in svalues for stuff going to 318 the ith processor 319 */ 320 start[0] = 0; 321 for (i = 1; i < size; i++) start[i] = start[i - 1] + sizes[2 * i - 2]; 322 for (i = 0; i < n; i++) { 323 j = owner[i]; 324 if (j != rank) { 325 ip = from_array[i]; 326 ia = to_array[i]; 327 sindices[start[j]++] = ip; 328 sindices[start[j]++] = ia; 329 } else { /* compute my own map */ 330 ip = from_array[i] - owners[rank]; 331 ia = to_array[i]; 332 aomap_loc[ip] = ia; 333 } 334 } 335 336 start[0] = 0; 337 for (i = 1; i < size; i++) start[i] = start[i - 1] + sizes[2 * i - 2]; 338 for (i = 0, count = 0; i < size; i++) { 339 if (sizes[2 * i + 1]) { 340 PetscCallMPI(MPI_Isend(sindices + start[i], sizes[2 * i], MPIU_INT, i, tag, comm, send_waits + count)); 341 count++; 342 } 343 } 344 PetscCheck(nsends == count, comm, PETSC_ERR_SUP, "nsends %" PetscInt_FMT " != count %" PetscInt_FMT, nsends, count); 345 346 /* wait on sends */ 347 if (nsends) PetscCallMPI(MPI_Waitall(nsends, send_waits, send_status)); 348 349 /* recvs */ 350 count = 0; 351 for (j = nreceives; j > 0; j--) { 352 PetscCallMPI(MPI_Waitany(nreceives, recv_waits, &widx, &recv_status)); 353 PetscCallMPI(MPI_Get_count(&recv_status, MPIU_INT, &nindices)); 354 rbuf = rindices + nmax * widx; /* global index */ 355 356 /* compute local mapping */ 357 for (i = 0; i < nindices; i += 2) { /* pack aomap_loc */ 358 ip = rbuf[i] - owners[rank]; /* local index */ 359 ia = rbuf[i + 1]; 360 aomap_loc[ip] = ia; 361 } 362 count++; 363 } 364 365 PetscCall(PetscFree(start)); 366 PetscCall(PetscFree3(sindices, send_waits, send_status)); 367 PetscCall(PetscFree2(rindices, recv_waits)); 368 PetscCall(PetscFree(owner)); 369 PetscCall(PetscFree(sizes)); 370 PetscFunctionReturn(PETSC_SUCCESS); 371 } 372 373 PETSC_EXTERN PetscErrorCode AOCreate_MemoryScalable(AO ao) 374 { 375 IS isapp = ao->isapp, ispetsc = ao->ispetsc; 376 const PetscInt *mypetsc, *myapp; 377 PetscInt napp, n_local, N, i, start, *petsc, *lens, *disp; 378 MPI_Comm comm; 379 AO_MemoryScalable *aomems; 380 PetscLayout map; 381 PetscMPIInt size, rank; 382 383 PetscFunctionBegin; 384 PetscCheck(isapp, PetscObjectComm((PetscObject)ao), PETSC_ERR_ARG_WRONGSTATE, "AOSetIS() must be called before AOSetType()"); 385 /* create special struct aomems */ 386 PetscCall(PetscNew(&aomems)); 387 ao->data = (void *)aomems; 388 ao->ops[0] = AOOps_MemoryScalable; 389 PetscCall(PetscObjectChangeTypeName((PetscObject)ao, AOMEMORYSCALABLE)); 390 391 /* transmit all local lengths of isapp to all processors */ 392 PetscCall(PetscObjectGetComm((PetscObject)isapp, &comm)); 393 PetscCallMPI(MPI_Comm_size(comm, &size)); 394 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 395 PetscCall(PetscMalloc2(size, &lens, size, &disp)); 396 PetscCall(ISGetLocalSize(isapp, &napp)); 397 PetscCallMPI(MPI_Allgather(&napp, 1, MPIU_INT, lens, 1, MPIU_INT, comm)); 398 399 N = 0; 400 for (i = 0; i < size; i++) { 401 disp[i] = N; 402 N += lens[i]; 403 } 404 405 /* If ispetsc is 0 then use "natural" numbering */ 406 if (napp) { 407 if (!ispetsc) { 408 start = disp[rank]; 409 PetscCall(PetscMalloc1(napp + 1, &petsc)); 410 for (i = 0; i < napp; i++) petsc[i] = start + i; 411 } else { 412 PetscCall(ISGetIndices(ispetsc, &mypetsc)); 413 petsc = (PetscInt *)mypetsc; 414 } 415 } else { 416 petsc = NULL; 417 } 418 419 /* create a map with global size N - used to determine the local sizes of ao - shall we use local napp instead of N? */ 420 PetscCall(PetscLayoutCreate(comm, &map)); 421 map->bs = 1; 422 map->N = N; 423 PetscCall(PetscLayoutSetUp(map)); 424 425 ao->N = N; 426 ao->n = map->n; 427 aomems->map = map; 428 429 /* create distributed indices app_loc: petsc->app and petsc_loc: app->petsc */ 430 n_local = map->n; 431 PetscCall(PetscCalloc2(n_local, &aomems->app_loc, n_local, &aomems->petsc_loc)); 432 PetscCall(ISGetIndices(isapp, &myapp)); 433 434 PetscCall(AOCreateMemoryScalable_private(comm, napp, petsc, myapp, ao, aomems->app_loc)); 435 PetscCall(AOCreateMemoryScalable_private(comm, napp, myapp, petsc, ao, aomems->petsc_loc)); 436 437 PetscCall(ISRestoreIndices(isapp, &myapp)); 438 if (napp) { 439 if (ispetsc) { 440 PetscCall(ISRestoreIndices(ispetsc, &mypetsc)); 441 } else { 442 PetscCall(PetscFree(petsc)); 443 } 444 } 445 PetscCall(PetscFree2(lens, disp)); 446 PetscFunctionReturn(PETSC_SUCCESS); 447 } 448 449 /*@C 450 AOCreateMemoryScalable - Creates a memory scalable application ordering using two integer arrays. 451 452 Collective 453 454 Input Parameters: 455 + comm - MPI communicator that is to share the `AO` 456 . napp - size of integer arrays 457 . myapp - integer array that defines an ordering 458 - mypetsc - integer array that defines another ordering (may be NULL to 459 indicate the natural ordering, that is 0,1,2,3,...) 460 461 Output Parameter: 462 . aoout - the new application ordering 463 464 Level: beginner 465 466 Note: 467 The arrays myapp and mypetsc must contain the all the integers 0 to napp-1 with no duplicates; that is there cannot be any "holes" 468 in the indices. Use `AOCreateMapping()` or `AOCreateMappingIS()` if you wish to have "holes" in the indices. 469 Comparing with `AOCreateBasic()`, this routine trades memory with message communication. 470 471 .seealso: [](sec_ao), [](sec_scatter), `AO`, `AOCreateMemoryScalableIS()`, `AODestroy()`, `AOPetscToApplication()`, `AOApplicationToPetsc()` 472 @*/ 473 PetscErrorCode AOCreateMemoryScalable(MPI_Comm comm, PetscInt napp, const PetscInt myapp[], const PetscInt mypetsc[], AO *aoout) 474 { 475 IS isapp, ispetsc; 476 const PetscInt *app = myapp, *petsc = mypetsc; 477 478 PetscFunctionBegin; 479 PetscCall(ISCreateGeneral(comm, napp, app, PETSC_USE_POINTER, &isapp)); 480 if (mypetsc) { 481 PetscCall(ISCreateGeneral(comm, napp, petsc, PETSC_USE_POINTER, &ispetsc)); 482 } else { 483 ispetsc = NULL; 484 } 485 PetscCall(AOCreateMemoryScalableIS(isapp, ispetsc, aoout)); 486 PetscCall(ISDestroy(&isapp)); 487 if (mypetsc) PetscCall(ISDestroy(&ispetsc)); 488 PetscFunctionReturn(PETSC_SUCCESS); 489 } 490 491 /*@C 492 AOCreateMemoryScalableIS - Creates a memory scalable application ordering using two index sets. 493 494 Collective 495 496 Input Parameters: 497 + isapp - index set that defines an ordering 498 - ispetsc - index set that defines another ordering (may be NULL to use the 499 natural ordering) 500 501 Output Parameter: 502 . aoout - the new application ordering 503 504 Level: beginner 505 506 Notes: 507 The index sets isapp and ispetsc must contain the all the integers 0 to napp-1 (where napp is the length of the index sets) with no duplicates; 508 that is there cannot be any "holes". 509 510 Comparing with `AOCreateBasicIS()`, this routine trades memory with message communication. 511 512 .seealso: [](sec_ao), [](sec_scatter), `AO`, `AOCreateBasicIS()`, `AOCreateMemoryScalable()`, `AODestroy()` 513 @*/ 514 PetscErrorCode AOCreateMemoryScalableIS(IS isapp, IS ispetsc, AO *aoout) 515 { 516 MPI_Comm comm; 517 AO ao; 518 519 PetscFunctionBegin; 520 PetscCall(PetscObjectGetComm((PetscObject)isapp, &comm)); 521 PetscCall(AOCreate(comm, &ao)); 522 PetscCall(AOSetIS(ao, isapp, ispetsc)); 523 PetscCall(AOSetType(ao, AOMEMORYSCALABLE)); 524 PetscCall(AOViewFromOptions(ao, NULL, "-ao_view")); 525 *aoout = ao; 526 PetscFunctionReturn(PETSC_SUCCESS); 527 } 528