1 // Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2 // the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3 // reserved. See files LICENSE and NOTICE for details. 4 // 5 // This file is part of CEED, a collection of benchmarks, miniapps, software 6 // libraries and APIs for efficient high-order finite element and spectral 7 // element discretizations for exascale applications. For more information and 8 // source code availability see http://github.com/ceed. 9 // 10 // The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11 // a collaborative effort of two U.S. Department of Energy organizations (Office 12 // of Science and the National Nuclear Security Administration) responsible for 13 // the planning and preparation of a capable exascale ecosystem, including 14 // software, applications, hardware, advanced system engineering and early 15 // testbed platforms, in support of the nation's exascale computing imperative. 16 17 #include <ceed-impl.h> 18 #include <ceed-backend.h> 19 20 /// @file 21 /// Implementation of CeedElemRestriction interfaces 22 23 /// ---------------------------------------------------------------------------- 24 /// CeedElemRestriction Library Internal Functions 25 /// ---------------------------------------------------------------------------- 26 /// @addtogroup CeedElemRestrictionDeveloper 27 /// @{ 28 29 /** 30 @brief Permute and pad offsets for a blocked restriction 31 32 @param offsets Array of shape [@a nelem, @a elemsize]. Row i holds the 33 ordered list of the offsets (into the input CeedVector) 34 for the unknowns corresponding to element i, where 35 0 <= i < @a nelem. All offsets must be in the range 36 [0, @a lsize - 1]. 37 @param blkoffsets Array of permuted and padded offsets of 38 shape [@a nblk, @a elemsize, @a blksize]. 39 @param nblk Number of blocks 40 @param nelem Number of elements 41 @param blksize Number of elements in a block 42 @param elemsize Size of each element 43 44 @return An error code: 0 - success, otherwise - failure 45 46 @ref Utility 47 **/ 48 int CeedPermutePadOffsets(const CeedInt *offsets, CeedInt *blkoffsets, 49 CeedInt nblk, CeedInt nelem, CeedInt blksize, 50 CeedInt elemsize) { 51 for (CeedInt e = 0; e < nblk*blksize; e+=blksize) 52 for (int j = 0; j < blksize; j++) 53 for (int k = 0; k < elemsize; k++) 54 blkoffsets[e*elemsize + k*blksize + j] 55 = offsets[CeedIntMin(e+j,nelem-1)*elemsize + k]; 56 return 0; 57 } 58 59 /// @} 60 61 /// ---------------------------------------------------------------------------- 62 /// CeedElemRestriction Backend API 63 /// ---------------------------------------------------------------------------- 64 /// @addtogroup CeedElemRestrictionBackend 65 /// @{ 66 67 /** 68 @brief Get the Ceed associated with a CeedElemRestriction 69 70 @param rstr CeedElemRestriction 71 @param[out] ceed Variable to store Ceed 72 73 @return An error code: 0 - success, otherwise - failure 74 75 @ref Backend 76 **/ 77 int CeedElemRestrictionGetCeed(CeedElemRestriction rstr, Ceed *ceed) { 78 *ceed = rstr->ceed; 79 return 0; 80 } 81 82 /** 83 84 @brief Get the strides of a strided CeedElemRestriction 85 86 @param rstr CeedElemRestriction 87 @param[out] strides Variable to store strides array 88 89 @return An error code: 0 - success, otherwise - failure 90 91 @ref Backend 92 **/ 93 int CeedElemRestrictionGetStrides(CeedElemRestriction rstr, 94 CeedInt (*strides)[3]) { 95 if (!rstr->strides) 96 // LCOV_EXCL_START 97 return CeedError(rstr->ceed, 1, "ElemRestriction has no stride data"); 98 // LCOV_EXCL_STOP 99 100 for (int i = 0; i<3; i++) 101 (*strides)[i] = rstr->strides[i]; 102 return 0; 103 } 104 105 /** 106 @brief Get read-only access to a CeedElemRestriction offsets array by memtype 107 108 @param rstr CeedElemRestriction to retrieve offsets 109 @param mtype Memory type on which to access the array. If the backend 110 uses a different memory type, this will perform a copy 111 (possibly cached). 112 @param[out] offsets Array on memory type mtype 113 114 @return An error code: 0 - success, otherwise - failure 115 116 @ref User 117 **/ 118 int CeedElemRestrictionGetOffsets(CeedElemRestriction rstr, CeedMemType mtype, 119 const CeedInt **offsets) { 120 int ierr; 121 122 if (!rstr->GetOffsets) 123 // LCOV_EXCL_START 124 return CeedError(rstr->ceed, 1, "Backend does not support GetOffsets"); 125 // LCOV_EXCL_STOP 126 127 ierr = rstr->GetOffsets(rstr, mtype, offsets); CeedChk(ierr); 128 rstr->numreaders++; 129 return 0; 130 } 131 132 /** 133 @brief Restore an offsets array obtained using CeedElemRestrictionGetOffsets() 134 135 @param rstr CeedElemRestriction to restore 136 @param offsets Array of offset data 137 138 @return An error code: 0 - success, otherwise - failure 139 140 @ref User 141 **/ 142 int CeedElemRestrictionRestoreOffsets(CeedElemRestriction rstr, 143 const CeedInt **offsets) { 144 *offsets = NULL; 145 rstr->numreaders--; 146 return 0; 147 } 148 149 /** 150 @brief Get the strided status of a CeedElemRestriction 151 152 @param rstr CeedElemRestriction 153 @param[out] status Variable to store strided status, 1 if strided else 0 154 155 @return An error code: 0 - success, otherwise - failure 156 157 @ref Backend 158 **/ 159 int CeedElemRestrictionGetStridedStatus(CeedElemRestriction rstr, 160 bool *status) { 161 *status = rstr->strides ? 1 : 0; 162 return 0; 163 } 164 165 /** 166 @brief Get the backend stride status of a CeedElemRestriction 167 168 @param rstr CeedElemRestriction 169 @param[out] status Variable to store stride status 170 171 @return An error code: 0 - success, otherwise - failure 172 173 @ref Backend 174 **/ 175 int CeedElemRestrictionGetBackendStridesStatus(CeedElemRestriction rstr, 176 bool *status) { 177 if (!rstr->strides) 178 // LCOV_EXCL_START 179 return CeedError(rstr->ceed, 1, "ElemRestriction has no stride data"); 180 // LCOV_EXCL_STOP 181 182 *status = ((rstr->strides[0] == CEED_STRIDES_BACKEND[0]) && 183 (rstr->strides[1] == CEED_STRIDES_BACKEND[1]) && 184 (rstr->strides[2] == CEED_STRIDES_BACKEND[2])); 185 return 0; 186 } 187 188 /** 189 190 @brief Get the E-vector layout of a CeedElemRestriction 191 192 @param rstr CeedElemRestriction 193 @param[out] layout Variable to store layout array, 194 stored as [nodes, components, elements]. 195 The data for node i, component j, element k in the 196 E-vector is given by 197 i*layout[0] + j*layout[1] + k*layout[2] 198 199 @return An error code: 0 - success, otherwise - failure 200 201 @ref Backend 202 **/ 203 int CeedElemRestrictionGetELayout(CeedElemRestriction rstr, 204 CeedInt (*layout)[3]) { 205 if (!rstr->layout[0]) 206 // LCOV_EXCL_START 207 return CeedError(rstr->ceed, 1, "ElemRestriction has no layout data"); 208 // LCOV_EXCL_STOP 209 210 for (int i = 0; i<3; i++) 211 (*layout)[i] = rstr->layout[i]; 212 return 0; 213 } 214 215 /** 216 217 @brief Set the E-vector layout of a CeedElemRestriction 218 219 @param rstr CeedElemRestriction 220 @param layout Variable to containing layout array, 221 stored as [nodes, components, elements]. 222 The data for node i, component j, element k in the 223 E-vector is given by 224 i*layout[0] + j*layout[1] + k*layout[2] 225 226 @return An error code: 0 - success, otherwise - failure 227 228 @ref Backend 229 **/ 230 int CeedElemRestrictionSetELayout(CeedElemRestriction rstr, 231 CeedInt layout[3]) { 232 for (int i = 0; i<3; i++) 233 rstr->layout[i] = layout[i]; 234 return 0; 235 } 236 237 /** 238 @brief Get the backend data of a CeedElemRestriction 239 240 @param rstr CeedElemRestriction 241 @param[out] data Variable to store data 242 243 @return An error code: 0 - success, otherwise - failure 244 245 @ref Backend 246 **/ 247 int CeedElemRestrictionGetData(CeedElemRestriction rstr, void **data) { 248 *data = rstr->data; 249 return 0; 250 } 251 252 /** 253 @brief Set the backend data of a CeedElemRestriction 254 255 @param[out] rstr CeedElemRestriction 256 @param data Data to set 257 258 @return An error code: 0 - success, otherwise - failure 259 260 @ref Backend 261 **/ 262 int CeedElemRestrictionSetData(CeedElemRestriction rstr, void **data) { 263 rstr->data = *data; 264 return 0; 265 } 266 267 /// @} 268 269 /// @cond DOXYGEN_SKIP 270 static struct CeedElemRestriction_private ceed_elemrestriction_none; 271 /// @endcond 272 273 /// ---------------------------------------------------------------------------- 274 /// CeedElemRestriction Public API 275 /// ---------------------------------------------------------------------------- 276 /// @addtogroup CeedElemRestrictionUser 277 /// @{ 278 279 /// Indicate that the stride is determined by the backend 280 const CeedInt CEED_STRIDES_BACKEND[3] = {}; 281 282 /// Indicate that no CeedElemRestriction is provided by the user 283 const CeedElemRestriction CEED_ELEMRESTRICTION_NONE = 284 &ceed_elemrestriction_none; 285 286 /** 287 @brief Create a CeedElemRestriction 288 289 @param ceed A Ceed object where the CeedElemRestriction will be created 290 @param nelem Number of elements described in the @a offsets array 291 @param elemsize Size (number of "nodes") per element 292 @param ncomp Number of field components per interpolation node 293 (1 for scalar fields) 294 @param compstride Stride between components for the same L-vector "node". 295 Data for node i, component j, element k can be found in 296 the L-vector at index 297 offsets[i + j*elemsize] + k*compstride. 298 @param lsize The size of the L-vector. This vector may be larger than 299 the elements and fields given by this restriction. 300 @param mtype Memory type of the @a offsets array, see CeedMemType 301 @param cmode Copy mode for the @a offsets array, see CeedCopyMode 302 @param offsets Array of shape [@a nelem, @a elemsize]. Row i holds the 303 ordered list of the offsets (into the input CeedVector) 304 for the unknowns corresponding to element i, where 305 0 <= i < @a nelem. All offsets must be in the range 306 [0, @a lsize - 1]. 307 @param[out] rstr Address of the variable where the newly created 308 CeedElemRestriction will be stored 309 310 @return An error code: 0 - success, otherwise - failure 311 312 @ref User 313 **/ 314 int CeedElemRestrictionCreate(Ceed ceed, CeedInt nelem, CeedInt elemsize, 315 CeedInt ncomp, CeedInt compstride, 316 CeedInt lsize, CeedMemType mtype, 317 CeedCopyMode cmode, const CeedInt *offsets, 318 CeedElemRestriction *rstr) { 319 int ierr; 320 321 if (!ceed->ElemRestrictionCreate) { 322 Ceed delegate; 323 ierr = CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction"); 324 CeedChk(ierr); 325 326 if (!delegate) 327 // LCOV_EXCL_START 328 return CeedError(ceed, 1, "Backend does not support ElemRestrictionCreate"); 329 // LCOV_EXCL_STOP 330 331 ierr = CeedElemRestrictionCreate(delegate, nelem, elemsize, ncomp, 332 compstride, lsize, mtype, cmode, 333 offsets, rstr); CeedChk(ierr); 334 return 0; 335 } 336 337 ierr = CeedCalloc(1, rstr); CeedChk(ierr); 338 (*rstr)->ceed = ceed; 339 ceed->refcount++; 340 (*rstr)->refcount = 1; 341 (*rstr)->nelem = nelem; 342 (*rstr)->elemsize = elemsize; 343 (*rstr)->ncomp = ncomp; 344 (*rstr)->compstride = compstride; 345 (*rstr)->lsize = lsize; 346 (*rstr)->nblk = nelem; 347 (*rstr)->blksize = 1; 348 ierr = ceed->ElemRestrictionCreate(mtype, cmode, offsets, *rstr); 349 CeedChk(ierr); 350 return 0; 351 } 352 353 /** 354 @brief Create a strided CeedElemRestriction 355 356 @param ceed A Ceed object where the CeedElemRestriction will be created 357 @param nelem Number of elements described by the restriction 358 @param elemsize Size (number of "nodes") per element 359 @param ncomp Number of field components per interpolation "node" 360 (1 for scalar fields) 361 @param lsize The size of the L-vector. This vector may be larger than 362 the elements and fields given by this restriction. 363 @param strides Array for strides between [nodes, components, elements]. 364 Data for node i, component j, element k can be found in 365 the L-vector at index 366 i*strides[0] + j*strides[1] + k*strides[2]. 367 @a CEED_STRIDES_BACKEND may be used with vectors created 368 by a Ceed backend. 369 @param rstr Address of the variable where the newly created 370 CeedElemRestriction will be stored 371 372 @return An error code: 0 - success, otherwise - failure 373 374 @ref User 375 **/ 376 int CeedElemRestrictionCreateStrided(Ceed ceed, CeedInt nelem, CeedInt elemsize, 377 CeedInt ncomp, CeedInt lsize, 378 const CeedInt strides[3], 379 CeedElemRestriction *rstr) { 380 int ierr; 381 382 if (!ceed->ElemRestrictionCreate) { 383 Ceed delegate; 384 ierr = CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction"); 385 CeedChk(ierr); 386 387 if (!delegate) 388 // LCOV_EXCL_START 389 return CeedError(ceed, 1, "Backend does not support ElemRestrictionCreate"); 390 // LCOV_EXCL_STOP 391 392 ierr = CeedElemRestrictionCreateStrided(delegate, nelem, elemsize, ncomp, 393 lsize, strides, rstr); 394 CeedChk(ierr); 395 return 0; 396 } 397 398 ierr = CeedCalloc(1, rstr); CeedChk(ierr); 399 (*rstr)->ceed = ceed; 400 ceed->refcount++; 401 (*rstr)->refcount = 1; 402 (*rstr)->nelem = nelem; 403 (*rstr)->elemsize = elemsize; 404 (*rstr)->ncomp = ncomp; 405 (*rstr)->lsize = lsize; 406 (*rstr)->nblk = nelem; 407 (*rstr)->blksize = 1; 408 ierr = CeedMalloc(3, &(*rstr)->strides); CeedChk(ierr); 409 for (int i = 0; i<3; i++) 410 (*rstr)->strides[i] = strides[i]; 411 ierr = ceed->ElemRestrictionCreate(CEED_MEM_HOST, CEED_OWN_POINTER, NULL, 412 *rstr); 413 CeedChk(ierr); 414 return 0; 415 } 416 417 /** 418 @brief Create a blocked CeedElemRestriction, typically only called by backends 419 420 @param ceed A Ceed object where the CeedElemRestriction will be created. 421 @param nelem Number of elements described in the @a offsets array. 422 @param elemsize Size (number of unknowns) per element 423 @param blksize Number of elements in a block 424 @param ncomp Number of field components per interpolation node 425 (1 for scalar fields) 426 @param compstride Stride between components for the same L-vector "node". 427 Data for node i, component j, element k can be found in 428 the L-vector at index 429 offsets[i + j*elemsize] + k*compstride. 430 @param lsize The size of the L-vector. This vector may be larger than 431 the elements and fields given by this restriction. 432 @param mtype Memory type of the @a offsets array, see CeedMemType 433 @param cmode Copy mode for the @a offsets array, see CeedCopyMode 434 @param offsets Array of shape [@a nelem, @a elemsize]. Row i holds the 435 ordered list of the offsets (into the input CeedVector) 436 for the unknowns corresponding to element i, where 437 0 <= i < @a nelem. All offsets must be in the range 438 [0, @a lsize - 1]. The backend will permute and pad this 439 array to the desired ordering for the blocksize, which is 440 typically given by the backend. The default reordering is 441 to interlace elements. 442 @param rstr Address of the variable where the newly created 443 CeedElemRestriction will be stored 444 445 @return An error code: 0 - success, otherwise - failure 446 447 @ref Backend 448 **/ 449 int CeedElemRestrictionCreateBlocked(Ceed ceed, CeedInt nelem, CeedInt elemsize, 450 CeedInt blksize, CeedInt ncomp, 451 CeedInt compstride, CeedInt lsize, 452 CeedMemType mtype, CeedCopyMode cmode, 453 const CeedInt *offsets, 454 CeedElemRestriction *rstr) { 455 int ierr; 456 CeedInt *blkoffsets; 457 CeedInt nblk = (nelem / blksize) + !!(nelem % blksize); 458 459 if (!ceed->ElemRestrictionCreateBlocked) { 460 Ceed delegate; 461 ierr = CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction"); 462 CeedChk(ierr); 463 464 if (!delegate) 465 // LCOV_EXCL_START 466 return CeedError(ceed, 1, "Backend does not support " 467 "ElemRestrictionCreateBlocked"); 468 // LCOV_EXCL_STOP 469 470 ierr = CeedElemRestrictionCreateBlocked(delegate, nelem, elemsize, blksize, 471 ncomp, compstride, lsize, mtype, 472 cmode, offsets, rstr); 473 CeedChk(ierr); 474 return 0; 475 } 476 477 ierr = CeedCalloc(1, rstr); CeedChk(ierr); 478 479 ierr = CeedCalloc(nblk*blksize*elemsize, &blkoffsets); CeedChk(ierr); 480 ierr = CeedPermutePadOffsets(offsets, blkoffsets, nblk, nelem, blksize, 481 elemsize); 482 CeedChk(ierr); 483 484 (*rstr)->ceed = ceed; 485 ceed->refcount++; 486 (*rstr)->refcount = 1; 487 (*rstr)->nelem = nelem; 488 (*rstr)->elemsize = elemsize; 489 (*rstr)->ncomp = ncomp; 490 (*rstr)->compstride = compstride; 491 (*rstr)->lsize = lsize; 492 (*rstr)->nblk = nblk; 493 (*rstr)->blksize = blksize; 494 ierr = ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, 495 (const CeedInt *) blkoffsets, *rstr); CeedChk(ierr); 496 497 if (cmode == CEED_OWN_POINTER) { 498 ierr = CeedFree(&offsets); CeedChk(ierr); 499 } 500 501 return 0; 502 } 503 504 /** 505 @brief Create a blocked strided CeedElemRestriction 506 507 @param ceed A Ceed object where the CeedElemRestriction will be created 508 @param nelem Number of elements described by the restriction 509 @param elemsize Size (number of "nodes") per element 510 @param blksize Number of elements in a block 511 @param ncomp Number of field components per interpolation node 512 (1 for scalar fields) 513 @param lsize The size of the L-vector. This vector may be larger than 514 the elements and fields given by this restriction. 515 @param strides Array for strides between [nodes, components, elements]. 516 Data for node i, component j, element k can be found in 517 the L-vector at index 518 i*strides[0] + j*strides[1] + k*strides[2]. 519 @a CEED_STRIDES_BACKEND may be used with vectors created 520 by a Ceed backend. 521 @param rstr Address of the variable where the newly created 522 CeedElemRestriction will be stored 523 524 @return An error code: 0 - success, otherwise - failure 525 526 @ref User 527 **/ 528 int CeedElemRestrictionCreateBlockedStrided(Ceed ceed, CeedInt nelem, 529 CeedInt elemsize, CeedInt blksize, CeedInt ncomp, CeedInt lsize, 530 const CeedInt strides[3], CeedElemRestriction *rstr) { 531 int ierr; 532 CeedInt nblk = (nelem / blksize) + !!(nelem % blksize); 533 534 if (!ceed->ElemRestrictionCreateBlocked) { 535 Ceed delegate; 536 ierr = CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction"); 537 CeedChk(ierr); 538 539 if (!delegate) 540 // LCOV_EXCL_START 541 return CeedError(ceed, 1, "Backend does not support " 542 "ElemRestrictionCreateBlocked"); 543 // LCOV_EXCL_STOP 544 545 ierr = CeedElemRestrictionCreateBlockedStrided(delegate, nelem, elemsize, 546 blksize, ncomp, lsize, strides, rstr); 547 CeedChk(ierr); 548 return 0; 549 } 550 551 ierr = CeedCalloc(1, rstr); CeedChk(ierr); 552 553 (*rstr)->ceed = ceed; 554 ceed->refcount++; 555 (*rstr)->refcount = 1; 556 (*rstr)->nelem = nelem; 557 (*rstr)->elemsize = elemsize; 558 (*rstr)->ncomp = ncomp; 559 (*rstr)->lsize = lsize; 560 (*rstr)->nblk = nblk; 561 (*rstr)->blksize = blksize; 562 ierr = CeedMalloc(3, &(*rstr)->strides); CeedChk(ierr); 563 for (int i = 0; i<3; i++) 564 (*rstr)->strides[i] = strides[i]; 565 ierr = ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, 566 NULL, *rstr); CeedChk(ierr); 567 568 return 0; 569 } 570 571 /** 572 @brief Create CeedVectors associated with a CeedElemRestriction 573 574 @param rstr CeedElemRestriction 575 @param lvec The address of the L-vector to be created, or NULL 576 @param evec The address of the E-vector to be created, or NULL 577 578 @return An error code: 0 - success, otherwise - failure 579 580 @ref User 581 **/ 582 int CeedElemRestrictionCreateVector(CeedElemRestriction rstr, CeedVector *lvec, 583 CeedVector *evec) { 584 int ierr; 585 CeedInt n, m; 586 m = rstr->lsize; 587 n = rstr->nblk * rstr->blksize * rstr->elemsize * rstr->ncomp; 588 if (lvec) { 589 ierr = CeedVectorCreate(rstr->ceed, m, lvec); CeedChk(ierr); 590 } 591 if (evec) { 592 ierr = CeedVectorCreate(rstr->ceed, n, evec); CeedChk(ierr); 593 } 594 return 0; 595 } 596 597 /** 598 @brief Restrict an L-vector to an E-vector or apply its transpose 599 600 @param rstr CeedElemRestriction 601 @param tmode Apply restriction or transpose 602 @param u Input vector (of size @a lsize when tmode=@ref CEED_NOTRANSPOSE) 603 @param ru Output vector (of shape [@a nelem * @a elemsize] when 604 tmode=@ref CEED_NOTRANSPOSE). Ordering of the e-vector is decided 605 by the backend. 606 @param request Request or @ref CEED_REQUEST_IMMEDIATE 607 608 @return An error code: 0 - success, otherwise - failure 609 610 @ref User 611 **/ 612 int CeedElemRestrictionApply(CeedElemRestriction rstr, CeedTransposeMode tmode, 613 CeedVector u, CeedVector ru, 614 CeedRequest *request) { 615 CeedInt m,n; 616 int ierr; 617 618 if (tmode == CEED_NOTRANSPOSE) { 619 m = rstr->nblk * rstr->blksize * rstr->elemsize * rstr->ncomp; 620 n = rstr->lsize; 621 } else { 622 m = rstr->lsize; 623 n = rstr->nblk * rstr->blksize * rstr->elemsize * rstr->ncomp; 624 } 625 if (n != u->length) 626 // LCOV_EXCL_START 627 return CeedError(rstr->ceed, 2, "Input vector size %d not compatible with " 628 "element restriction (%d, %d)", u->length, m, n); 629 // LCOV_EXCL_STOP 630 if (m != ru->length) 631 // LCOV_EXCL_START 632 return CeedError(rstr->ceed, 2, "Output vector size %d not compatible with " 633 "element restriction (%d, %d)", ru->length, m, n); 634 // LCOV_EXCL_STOP 635 ierr = rstr->Apply(rstr, tmode, u, ru, request); CeedChk(ierr); 636 637 return 0; 638 } 639 640 /** 641 @brief Restrict an L-vector to a block of an E-vector or apply its transpose 642 643 @param rstr CeedElemRestriction 644 @param block Block number to restrict to/from, i.e. block=0 will handle 645 elements [0 : blksize] and block=3 will handle elements 646 [3*blksize : 4*blksize] 647 @param tmode Apply restriction or transpose 648 @param u Input vector (of size @a lsize when tmode=@ref CEED_NOTRANSPOSE) 649 @param ru Output vector (of shape [@a blksize * @a elemsize] when 650 tmode=@ref CEED_NOTRANSPOSE). Ordering of the e-vector is decided 651 by the backend. 652 @param request Request or @ref CEED_REQUEST_IMMEDIATE 653 654 @return An error code: 0 - success, otherwise - failure 655 656 @ref Backend 657 **/ 658 int CeedElemRestrictionApplyBlock(CeedElemRestriction rstr, CeedInt block, 659 CeedTransposeMode tmode, CeedVector u, 660 CeedVector ru, CeedRequest *request) { 661 CeedInt m,n; 662 int ierr; 663 664 if (tmode == CEED_NOTRANSPOSE) { 665 m = rstr->blksize * rstr->elemsize * rstr->ncomp; 666 n = rstr->lsize; 667 } else { 668 m = rstr->lsize; 669 n = rstr->blksize * rstr->elemsize * rstr->ncomp; 670 } 671 if (n != u->length) 672 // LCOV_EXCL_START 673 return CeedError(rstr->ceed, 2, "Input vector size %d not compatible with " 674 "element restriction (%d, %d)", u->length, m, n); 675 // LCOV_EXCL_STOP 676 if (m != ru->length) 677 // LCOV_EXCL_START 678 return CeedError(rstr->ceed, 2, "Output vector size %d not compatible with " 679 "element restriction (%d, %d)", ru->length, m, n); 680 // LCOV_EXCL_STOP 681 if (rstr->blksize*block > rstr->nelem) 682 // LCOV_EXCL_START 683 return CeedError(rstr->ceed, 2, "Cannot retrieve block %d, element %d > " 684 "total elements %d", block, rstr->blksize*block, 685 rstr->nelem); 686 // LCOV_EXCL_STOP 687 ierr = rstr->ApplyBlock(rstr, block, tmode, u, ru, request); 688 CeedChk(ierr); 689 690 return 0; 691 } 692 693 /** 694 @brief Get the L-vector component stride 695 696 @param rstr CeedElemRestriction 697 @param[out] compstride Variable to store component stride 698 699 @return An error code: 0 - success, otherwise - failure 700 701 @ref Backend 702 **/ 703 int CeedElemRestrictionGetCompStride(CeedElemRestriction rstr, 704 CeedInt *compstride) { 705 *compstride = rstr->compstride; 706 return 0; 707 } 708 709 /** 710 @brief Get the total number of elements in the range of a CeedElemRestriction 711 712 @param rstr CeedElemRestriction 713 @param[out] numelem Variable to store number of elements 714 715 @return An error code: 0 - success, otherwise - failure 716 717 @ref Backend 718 **/ 719 int CeedElemRestrictionGetNumElements(CeedElemRestriction rstr, 720 CeedInt *numelem) { 721 *numelem = rstr->nelem; 722 return 0; 723 } 724 725 /** 726 @brief Get the size of elements in the CeedElemRestriction 727 728 @param rstr CeedElemRestriction 729 @param[out] elemsize Variable to store size of elements 730 731 @return An error code: 0 - success, otherwise - failure 732 733 @ref Backend 734 **/ 735 int CeedElemRestrictionGetElementSize(CeedElemRestriction rstr, 736 CeedInt *elemsize) { 737 *elemsize = rstr->elemsize; 738 return 0; 739 } 740 741 /** 742 @brief Get the size of the l-vector for a CeedElemRestriction 743 744 @param rstr CeedElemRestriction 745 @param[out] numnodes Variable to store number of nodes 746 747 @return An error code: 0 - success, otherwise - failure 748 749 @ref Backend 750 **/ 751 int CeedElemRestrictionGetLVectorSize(CeedElemRestriction rstr, 752 CeedInt *lsize) { 753 *lsize = rstr->lsize; 754 return 0; 755 } 756 757 /** 758 @brief Get the number of components in the elements of a 759 CeedElemRestriction 760 761 @param rstr CeedElemRestriction 762 @param[out] numcomp Variable to store number of components 763 764 @return An error code: 0 - success, otherwise - failure 765 766 @ref Backend 767 **/ 768 int CeedElemRestrictionGetNumComponents(CeedElemRestriction rstr, 769 CeedInt *numcomp) { 770 *numcomp = rstr->ncomp; 771 return 0; 772 } 773 774 /** 775 @brief Get the number of blocks in a CeedElemRestriction 776 777 @param rstr CeedElemRestriction 778 @param[out] numblock Variable to store number of blocks 779 780 @return An error code: 0 - success, otherwise - failure 781 782 @ref Backend 783 **/ 784 int CeedElemRestrictionGetNumBlocks(CeedElemRestriction rstr, 785 CeedInt *numblock) { 786 *numblock = rstr->nblk; 787 return 0; 788 } 789 790 /** 791 @brief Get the size of blocks in the CeedElemRestriction 792 793 @param rstr CeedElemRestriction 794 @param[out] blksize Variable to store size of blocks 795 796 @return An error code: 0 - success, otherwise - failure 797 798 @ref Backend 799 **/ 800 int CeedElemRestrictionGetBlockSize(CeedElemRestriction rstr, 801 CeedInt *blksize) { 802 *blksize = rstr->blksize; 803 return 0; 804 } 805 806 /** 807 @brief Get the multiplicity of nodes in a CeedElemRestriction 808 809 @param rstr CeedElemRestriction 810 @param[out] mult Vector to store multiplicity (of size lsize) 811 812 @return An error code: 0 - success, otherwise - failure 813 814 @ref User 815 **/ 816 int CeedElemRestrictionGetMultiplicity(CeedElemRestriction rstr, 817 CeedVector mult) { 818 int ierr; 819 CeedVector evec; 820 821 // Create and set evec 822 ierr = CeedElemRestrictionCreateVector(rstr, NULL, &evec); CeedChk(ierr); 823 ierr = CeedVectorSetValue(evec, 1.0); CeedChk(ierr); 824 ierr = CeedVectorSetValue(mult, 0.0); CeedChk(ierr); 825 826 // Apply to get multiplicity 827 ierr = CeedElemRestrictionApply(rstr, CEED_TRANSPOSE, evec, mult, 828 CEED_REQUEST_IMMEDIATE); CeedChk(ierr); 829 830 // Cleanup 831 ierr = CeedVectorDestroy(&evec); CeedChk(ierr); 832 833 return 0; 834 } 835 836 /** 837 @brief View a CeedElemRestriction 838 839 @param[in] rstr CeedElemRestriction to view 840 @param[in] stream Stream to write; typically stdout/stderr or a file 841 842 @return Error code: 0 - success, otherwise - failure 843 844 @ref User 845 **/ 846 int CeedElemRestrictionView(CeedElemRestriction rstr, FILE *stream) { 847 char stridesstr[500]; 848 if (rstr->strides) 849 sprintf(stridesstr, "[%d, %d, %d]", rstr->strides[0], rstr->strides[1], 850 rstr->strides[2]); 851 else 852 sprintf(stridesstr, "%d", rstr->compstride); 853 854 fprintf(stream, "%sCeedElemRestriction from (%d, %d) to %d elements with %d " 855 "nodes each and %s %s\n", rstr->blksize > 1 ? "Blocked " : "", 856 rstr->lsize, rstr->ncomp, rstr->nelem, rstr->elemsize, 857 rstr->strides ? "strides" : "component stride", stridesstr); 858 return 0; 859 } 860 861 /** 862 @brief Destroy a CeedElemRestriction 863 864 @param rstr CeedElemRestriction to destroy 865 866 @return An error code: 0 - success, otherwise - failure 867 868 @ref User 869 **/ 870 int CeedElemRestrictionDestroy(CeedElemRestriction *rstr) { 871 int ierr; 872 873 if (!*rstr || --(*rstr)->refcount > 0) 874 return 0; 875 if ((*rstr)->numreaders) 876 return CeedError((*rstr)->ceed, 1, "Cannot destroy CeedElemRestriction, " 877 "a process has read access to the offset data"); 878 if ((*rstr)->Destroy) { 879 ierr = (*rstr)->Destroy(*rstr); CeedChk(ierr); 880 } 881 ierr = CeedFree(&(*rstr)->strides); CeedChk(ierr); 882 ierr = CeedDestroy(&(*rstr)->ceed); CeedChk(ierr); 883 ierr = CeedFree(rstr); CeedChk(ierr); 884 return 0; 885 } 886 887 /// @} 888