10036de2cSjeremylt /// @file 20036de2cSjeremylt /// Test creation, use, and destruction of a blocked element restriction with multiple components in the lvector, and libCEED owns ind pointer 30036de2cSjeremylt /// \test Test creation, use, and destruction of a blocked element restriction with multiple components in the lvector, and libCEED owns ind pointer 40036de2cSjeremylt #include <ceed.h> 50acb07cdSJeremy L Thompson #include <ceed/backend.h> 60acb07cdSJeremy L Thompson #include <stdlib.h> 70acb07cdSJeremy L Thompson #include <string.h> 80036de2cSjeremylt 90036de2cSjeremylt int main(int argc, char **argv) { 100036de2cSjeremylt Ceed ceed; 110036de2cSjeremylt CeedVector x, y; 12d1d35e2fSjeremylt CeedInt num_elem = 8; 130acb07cdSJeremy L Thompson CeedInt elem_size = 2; 140acb07cdSJeremy L Thompson CeedInt num_blk = 2; 15d1d35e2fSjeremylt CeedInt blk_size = 5; 16d1d35e2fSjeremylt CeedInt num_comp = 3; 170acb07cdSJeremy L Thompson CeedInt ind[elem_size * num_elem]; 180acb07cdSJeremy L Thompson CeedInt *ceed_ind = malloc(sizeof(CeedInt) * elem_size * num_elem); 19d1d35e2fSjeremylt CeedScalar a[num_comp * (num_elem + 1)]; 200acb07cdSJeremy L Thompson const CeedScalar *xx, *yy; 210acb07cdSJeremy L Thompson CeedInt layout[3]; 220036de2cSjeremylt CeedElemRestriction r; 230036de2cSjeremylt 240036de2cSjeremylt CeedInit(argv[1], &ceed); 250036de2cSjeremylt 260acb07cdSJeremy L Thompson CeedVectorCreate(ceed, num_comp * (num_elem + 1), &x); 270acb07cdSJeremy L Thompson for (CeedInt i = 0; i < num_elem + 1; i++) { 28d1d35e2fSjeremylt a[i + 0 * (num_elem + 1)] = 10 + i; 29d1d35e2fSjeremylt a[i + 1 * (num_elem + 1)] = 20 + i; 30d1d35e2fSjeremylt a[i + 2 * (num_elem + 1)] = 30 + i; 310036de2cSjeremylt } 320036de2cSjeremylt CeedVectorSetArray(x, CEED_MEM_HOST, CEED_USE_POINTER, a); 330acb07cdSJeremy L Thompson 34d1d35e2fSjeremylt for (CeedInt i = 0; i < num_elem; i++) { 350036de2cSjeremylt ind[2 * i + 0] = i; 360036de2cSjeremylt ind[2 * i + 1] = i + 1; 370036de2cSjeremylt } 380acb07cdSJeremy L Thompson memcpy(ceed_ind, ind, sizeof(CeedInt) * elem_size * num_elem); 39*2b730f8bSJeremy L Thompson CeedElemRestrictionCreateBlocked(ceed, num_elem, elem_size, blk_size, num_comp, num_elem + 1, num_comp * (num_elem + 1), CEED_MEM_HOST, 400acb07cdSJeremy L Thompson CEED_OWN_POINTER, ceed_ind, &r); 410acb07cdSJeremy L Thompson CeedVectorCreate(ceed, num_comp * num_blk * blk_size * elem_size, &y); 420036de2cSjeremylt 430036de2cSjeremylt // NoTranspose 440036de2cSjeremylt CeedElemRestrictionApply(r, CEED_NOTRANSPOSE, x, y, CEED_REQUEST_IMMEDIATE); 450acb07cdSJeremy L Thompson CeedVectorGetArrayRead(y, CEED_MEM_HOST, &yy); 460acb07cdSJeremy L Thompson CeedElemRestrictionGetELayout(r, &layout); 47*2b730f8bSJeremy L Thompson for (CeedInt i = 0; i < elem_size; i++) { // Node 48*2b730f8bSJeremy L Thompson for (CeedInt j = 0; j < num_comp; j++) { // Component 490acb07cdSJeremy L Thompson for (CeedInt k = 0; k < num_elem; k++) { // Element 500acb07cdSJeremy L Thompson CeedInt block = k / blk_size; 510acb07cdSJeremy L Thompson CeedInt elem = k % blk_size; 52*2b730f8bSJeremy L Thompson CeedInt index = (i * blk_size + elem) * layout[0] + j * layout[1] * blk_size + block * layout[2] * blk_size; 530acb07cdSJeremy L Thompson if (yy[index] != a[ind[k * elem_size + i] + j * (num_elem + 1)]) 540acb07cdSJeremy L Thompson // LCOV_EXCL_START 55*2b730f8bSJeremy L Thompson printf("Error in restricted array y[%" CeedInt_FMT "][%" CeedInt_FMT "][%" CeedInt_FMT "] = %f\n", i, j, k, (double)yy[index]); 560acb07cdSJeremy L Thompson // LCOV_EXCL_STOP 570acb07cdSJeremy L Thompson } 58*2b730f8bSJeremy L Thompson } 59*2b730f8bSJeremy L Thompson } 600acb07cdSJeremy L Thompson CeedVectorRestoreArrayRead(y, &yy); 610036de2cSjeremylt 620036de2cSjeremylt // Transpose 630acb07cdSJeremy L Thompson CeedVectorSetValue(x, 0); 640036de2cSjeremylt CeedElemRestrictionApply(r, CEED_TRANSPOSE, y, x, CEED_REQUEST_IMMEDIATE); 650acb07cdSJeremy L Thompson CeedVectorGetArrayRead(x, CEED_MEM_HOST, &xx); 660acb07cdSJeremy L Thompson for (CeedInt i = 0; i < num_elem + 1; i++) { 670acb07cdSJeremy L Thompson for (CeedInt j = 0; j < num_comp; j++) { 68*2b730f8bSJeremy L Thompson if (xx[i + j * (num_elem + 1)] != ((j + 1) * 10 + i) * (i > 0 && i < num_elem ? 2.0 : 1.0)) { 690acb07cdSJeremy L Thompson // LCOV_EXCL_START 70*2b730f8bSJeremy L Thompson printf("Error in restricted array x[%" CeedInt_FMT "][%" CeedInt_FMT "] = %f\n", j, i, (double)xx[i + j * (num_elem + 1)]); 710acb07cdSJeremy L Thompson // LCOV_EXCL_STOP 720acb07cdSJeremy L Thompson } 730acb07cdSJeremy L Thompson } 74*2b730f8bSJeremy L Thompson } 750acb07cdSJeremy L Thompson CeedVectorRestoreArrayRead(x, &xx); 760036de2cSjeremylt 770036de2cSjeremylt CeedVectorDestroy(&x); 780036de2cSjeremylt CeedVectorDestroy(&y); 790036de2cSjeremylt CeedElemRestrictionDestroy(&r); 800036de2cSjeremylt CeedDestroy(&ceed); 810036de2cSjeremylt return 0; 820036de2cSjeremylt } 83