From 0534fe31af6cf9aa1f88f25323138d938f5803e4 Mon Sep 17 00:00:00 2001 From: Sebastian Grimberg Date: Mon, 24 Apr 2023 16:19:44 -0700 Subject: [PATCH] WIP: enum CeedRestrictionType for CeedElemRestriction type --- backends/blocked/ceed-blocked-operator.c | 49 ++-- backends/cuda-ref/ceed-cuda-restriction.c | 6 +- backends/hip-ref/ceed-hip-ref-restriction.c | 6 +- backends/magma/ceed-magma-restriction.c | 6 +- backends/opt/ceed-opt-operator.c | 50 ++-- backends/ref/ceed-ref-restriction.c | 300 ++++++++++---------- backends/ref/ceed-ref.c | 4 - backends/ref/ceed-ref.h | 7 +- include/ceed-impl.h | 10 +- include/ceed/backend.h | 14 + interface/ceed-elemrestriction.c | 64 +++-- interface/ceed.c | 4 - 12 files changed, 277 insertions(+), 243 deletions(-) diff --git a/backends/blocked/ceed-blocked-operator.c b/backends/blocked/ceed-blocked-operator.c index cc2f2081dc..ec909a59a2 100644 --- a/backends/blocked/ceed-blocked-operator.c +++ b/backends/blocked/ceed-blocked-operator.c @@ -49,37 +49,40 @@ static int CeedOperatorSetupFields_Blocked(CeedQFunction qf, CeedOperator op, bo CeedCallBackend(CeedElemRestrictionGetElementSize(r, &elem_size)); CeedCallBackend(CeedElemRestrictionGetLVectorSize(r, &l_size)); CeedCallBackend(CeedElemRestrictionGetNumComponents(r, &num_comp)); - - bool strided; - CeedCallBackend(CeedElemRestrictionIsStrided(r, &strided)); - if (strided) { - CeedInt strides[3]; - CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); - CeedCallBackend( - CeedElemRestrictionCreateBlockedStrided(ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, &blk_restr[i + start_e])); - } else { - const CeedInt *offsets = NULL; - const bool *orients = NULL; - const CeedInt *curl_orients = NULL; - CeedCallBackend(CeedElemRestrictionGetOffsets(r, CEED_MEM_HOST, &offsets)); - CeedCallBackend(CeedElemRestrictionGetOrientations(r, CEED_MEM_HOST, &orients)); - CeedCallBackend(CeedElemRestrictionGetCurlOrientations(r, CEED_MEM_HOST, &curl_orients)); - CeedCallBackend(CeedElemRestrictionGetCompStride(r, &comp_stride)); - if (!orients && !curl_orients) { + CeedCallBackend(CeedElemRestrictionGetCompStride(r, &comp_stride)); + + CeedRestrictionType rstr_type; + CeedCallBackend(CeedElemRestrictionGetType(r, &rstr_type)); + const CeedInt *offsets = NULL; + CeedCallBackend(CeedElemRestrictionGetOffsets(r, CEED_MEM_HOST, &offsets)); + switch (rstr_type) { + case CEED_RESTRICTION_DEFAULT: { CeedCallBackend(CeedElemRestrictionCreateBlocked(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, &blk_restr[i + start_e])); - } else if (!curl_orients) { + } break; + case CEED_RESTRICTION_ORIENTED: { + const bool *orients = NULL; + CeedCallBackend(CeedElemRestrictionGetOrientations(r, CEED_MEM_HOST, &orients)); CeedCallBackend(CeedElemRestrictionCreateBlockedOriented(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, orients, &blk_restr[i + start_e])); - } else { + CeedCallBackend(CeedElemRestrictionRestoreOrientations(r, &orients)); + } break; + case CEED_RESTRICTION_CURL_ORIENTED: { + const CeedInt *curl_orients = NULL; + CeedCallBackend(CeedElemRestrictionGetCurlOrientations(r, CEED_MEM_HOST, &curl_orients)); CeedCallBackend(CeedElemRestrictionCreateBlockedCurlOriented(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, curl_orients, &blk_restr[i + start_e])); - } - CeedCallBackend(CeedElemRestrictionRestoreOffsets(r, &offsets)); - CeedCallBackend(CeedElemRestrictionRestoreOrientations(r, &orients)); - CeedCallBackend(CeedElemRestrictionRestoreCurlOrientations(r, &curl_orients)); + CeedCallBackend(CeedElemRestrictionRestoreCurlOrientations(r, &curl_orients)); + } break; + case CEED_RESTRICTION_STRIDED: { + CeedInt strides[3]; + CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + CeedCallBackend( + CeedElemRestrictionCreateBlockedStrided(ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, &blk_restr[i + start_e])); + } break; } + CeedCallBackend(CeedElemRestrictionRestoreOffsets(r, &offsets)); CeedCallBackend(CeedElemRestrictionCreateVector(blk_restr[i + start_e], NULL, &e_vecs_full[i + start_e])); } diff --git a/backends/cuda-ref/ceed-cuda-restriction.c b/backends/cuda-ref/ceed-cuda-restriction.c index 5ce2f92498..eec735881c 100644 --- a/backends/cuda-ref/ceed-cuda-restriction.c +++ b/backends/cuda-ref/ceed-cuda-restriction.c @@ -221,7 +221,8 @@ static int CeedElemRestrictionOffset_Cuda(const CeedElemRestriction r, const Cee //------------------------------------------------------------------------------ // Create restriction //------------------------------------------------------------------------------ -int CeedElemRestrictionCreate_Cuda(CeedMemType m_type, CeedCopyMode copy_mode, const CeedInt *indices, CeedElemRestriction r) { +int CeedElemRestrictionCreate_Cuda(CeedMemType m_type, CeedCopyMode copy_mode, const CeedInt *indices, const bool *orients, + const CeedInt *curl_orients, CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); CeedElemRestriction_Cuda *impl; @@ -343,7 +344,8 @@ int CeedElemRestrictionCreate_Cuda(CeedMemType m_type, CeedCopyMode copy_mode, c //------------------------------------------------------------------------------ // Blocked not supported //------------------------------------------------------------------------------ -int CeedElemRestrictionCreateBlocked_Cuda(const CeedMemType m_type, const CeedCopyMode copy_mode, const CeedInt *indices, CeedElemRestriction r) { +int CeedElemRestrictionCreateBlocked_Cuda(const CeedMemType m_type, const CeedCopyMode copy_mode, const CeedInt *indices, const bool *orients, + const CeedInt *curl_orients, CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); return CeedError(ceed, CEED_ERROR_BACKEND, "Backend does not implement blocked restrictions"); diff --git a/backends/hip-ref/ceed-hip-ref-restriction.c b/backends/hip-ref/ceed-hip-ref-restriction.c index d72975031d..a6243aa36a 100644 --- a/backends/hip-ref/ceed-hip-ref-restriction.c +++ b/backends/hip-ref/ceed-hip-ref-restriction.c @@ -219,7 +219,8 @@ static int CeedElemRestrictionOffset_Hip(const CeedElemRestriction r, const Ceed //------------------------------------------------------------------------------ // Create restriction //------------------------------------------------------------------------------ -int CeedElemRestrictionCreate_Hip(CeedMemType mtype, CeedCopyMode cmode, const CeedInt *indices, CeedElemRestriction r) { +int CeedElemRestrictionCreate_Hip(CeedMemType mtype, CeedCopyMode cmode, const CeedInt *indices, const bool *orients, const CeedInt *curl_orients, + CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); CeedElemRestriction_Hip *impl; @@ -341,7 +342,8 @@ int CeedElemRestrictionCreate_Hip(CeedMemType mtype, CeedCopyMode cmode, const C //------------------------------------------------------------------------------ // Blocked not supported //------------------------------------------------------------------------------ -int CeedElemRestrictionCreateBlocked_Hip(const CeedMemType mtype, const CeedCopyMode cmode, const CeedInt *indices, CeedElemRestriction r) { +int CeedElemRestrictionCreateBlocked_Hip(const CeedMemType mtype, const CeedCopyMode cmode, const CeedInt *indices, const bool *orients, + const CeedInt *curl_orients, CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); return CeedError(ceed, CEED_ERROR_BACKEND, "Backend does not implement blocked restrictions"); diff --git a/backends/magma/ceed-magma-restriction.c b/backends/magma/ceed-magma-restriction.c index 32038e0988..d89c79396f 100644 --- a/backends/magma/ceed-magma-restriction.c +++ b/backends/magma/ceed-magma-restriction.c @@ -160,7 +160,8 @@ static int CeedElemRestrictionDestroy_Magma(CeedElemRestriction r) { return CEED_ERROR_SUCCESS; } -int CeedElemRestrictionCreate_Magma(CeedMemType mtype, CeedCopyMode cmode, const CeedInt *offsets, CeedElemRestriction r) { +int CeedElemRestrictionCreate_Magma(CeedMemType mtype, CeedCopyMode cmode, const CeedInt *offsets, const bool *orients, const CeedInt *curl_orients, + CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); @@ -278,7 +279,8 @@ int CeedElemRestrictionCreate_Magma(CeedMemType mtype, CeedCopyMode cmode, const return CEED_ERROR_SUCCESS; } -int CeedElemRestrictionCreateBlocked_Magma(const CeedMemType mtype, const CeedCopyMode cmode, const CeedInt *offsets, const CeedElemRestriction r) { +int CeedElemRestrictionCreateBlocked_Magma(const CeedMemType mtype, const CeedCopyMode cmode, const CeedInt *offsets, const bool *orients, + const CeedInt *curl_orients, const CeedElemRestriction r) { Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); return CeedError(ceed, CEED_ERROR_BACKEND, "Backend does not implement blocked restrictions"); diff --git a/backends/opt/ceed-opt-operator.c b/backends/opt/ceed-opt-operator.c index 901101cd11..c2504803bd 100644 --- a/backends/opt/ceed-opt-operator.c +++ b/backends/opt/ceed-opt-operator.c @@ -42,7 +42,6 @@ static int CeedOperatorSetupFields_Opt(CeedQFunction qf, CeedOperator op, bool i if (eval_mode != CEED_EVAL_WEIGHT) { CeedCallBackend(CeedOperatorFieldGetElemRestriction(op_fields[i], &r)); - Ceed ceed; CeedCallBackend(CeedElemRestrictionGetCeed(r, &ceed)); CeedSize l_size; CeedInt num_elem, elem_size, comp_stride; @@ -50,37 +49,40 @@ static int CeedOperatorSetupFields_Opt(CeedQFunction qf, CeedOperator op, bool i CeedCallBackend(CeedElemRestrictionGetElementSize(r, &elem_size)); CeedCallBackend(CeedElemRestrictionGetLVectorSize(r, &l_size)); CeedCallBackend(CeedElemRestrictionGetNumComponents(r, &num_comp)); - - bool strided; - CeedCallBackend(CeedElemRestrictionIsStrided(r, &strided)); - if (strided) { - CeedInt strides[3]; - CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); - CeedCallBackend( - CeedElemRestrictionCreateBlockedStrided(ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, &blk_restr[i + start_e])); - } else { - const CeedInt *offsets = NULL; - const bool *orients = NULL; - const CeedInt *curl_orients = NULL; - CeedCallBackend(CeedElemRestrictionGetOffsets(r, CEED_MEM_HOST, &offsets)); - CeedCallBackend(CeedElemRestrictionGetOrientations(r, CEED_MEM_HOST, &orients)); - CeedCallBackend(CeedElemRestrictionGetCurlOrientations(r, CEED_MEM_HOST, &curl_orients)); - CeedCallBackend(CeedElemRestrictionGetCompStride(r, &comp_stride)); - if (!orients && !curl_orients) { + CeedCallBackend(CeedElemRestrictionGetCompStride(r, &comp_stride)); + + CeedRestrictionType rstr_type; + CeedCallBackend(CeedElemRestrictionGetType(r, &rstr_type)); + const CeedInt *offsets = NULL; + CeedCallBackend(CeedElemRestrictionGetOffsets(r, CEED_MEM_HOST, &offsets)); + switch (rstr_type) { + case CEED_RESTRICTION_DEFAULT: { CeedCallBackend(CeedElemRestrictionCreateBlocked(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, &blk_restr[i + start_e])); - } else if (!curl_orients) { + } break; + case CEED_RESTRICTION_ORIENTED: { + const bool *orients = NULL; + CeedCallBackend(CeedElemRestrictionGetOrientations(r, CEED_MEM_HOST, &orients)); CeedCallBackend(CeedElemRestrictionCreateBlockedOriented(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, orients, &blk_restr[i + start_e])); - } else { + CeedCallBackend(CeedElemRestrictionRestoreOrientations(r, &orients)); + } break; + case CEED_RESTRICTION_CURL_ORIENTED: { + const CeedInt *curl_orients = NULL; + CeedCallBackend(CeedElemRestrictionGetCurlOrientations(r, CEED_MEM_HOST, &curl_orients)); CeedCallBackend(CeedElemRestrictionCreateBlockedCurlOriented(ceed, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, CEED_MEM_HOST, CEED_COPY_VALUES, offsets, curl_orients, &blk_restr[i + start_e])); - } - CeedCallBackend(CeedElemRestrictionRestoreOffsets(r, &offsets)); - CeedCallBackend(CeedElemRestrictionRestoreOrientations(r, &orients)); - CeedCallBackend(CeedElemRestrictionRestoreCurlOrientations(r, &curl_orients)); + CeedCallBackend(CeedElemRestrictionRestoreCurlOrientations(r, &curl_orients)); + } break; + case CEED_RESTRICTION_STRIDED: { + CeedInt strides[3]; + CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + CeedCallBackend( + CeedElemRestrictionCreateBlockedStrided(ceed, num_elem, elem_size, blk_size, num_comp, l_size, strides, &blk_restr[i + start_e])); + } break; } + CeedCallBackend(CeedElemRestrictionRestoreOffsets(r, &offsets)); CeedCallBackend(CeedElemRestrictionCreateVector(blk_restr[i + start_e], NULL, &e_vecs_full[i + start_e])); } diff --git a/backends/ref/ceed-ref-restriction.c b/backends/ref/ceed-ref-restriction.c index ccf3d6bc71..56a745c098 100644 --- a/backends/ref/ceed-ref-restriction.c +++ b/backends/ref/ceed-ref-restriction.c @@ -26,9 +26,8 @@ static inline int CeedElemRestrictionApply_Ref_Core(CeedElemRestriction r, const CeedCallBackend(CeedElemRestrictionGetNumElements(r, &num_elem)); CeedCallBackend(CeedElemRestrictionGetElementSize(r, &elem_size)); v_offset = start * blk_size * elem_size * num_comp; - bool is_oriented, is_curl_oriented; - is_oriented = (impl->orients != NULL); - is_curl_oriented = (impl->curl_orients != NULL); + CeedRestrictionType rstr_type; + CeedCallBackend(CeedElemRestrictionGetType(r, &rstr_type)); CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_HOST, &uu)); if (t_mode == CEED_TRANSPOSE) { @@ -40,55 +39,69 @@ static inline int CeedElemRestrictionApply_Ref_Core(CeedElemRestriction r, const } // Restriction from L-vector to E-vector // Perform: v = r * u + // vv has shape [elem_size, num_comp, num_elem], row-major + // uu has shape [nnodes, num_comp] if (t_mode == CEED_NOTRANSPOSE) { - // No offsets provided, Identity Restriction - if (!impl->offsets) { - bool has_backend_strides; - CeedCallBackend(CeedElemRestrictionHasBackendStrides(r, &has_backend_strides)); - if (has_backend_strides) { - // CPU backend strides are {1, elem_size, elem_size*num_comp} - // This if branch is left separate to allow better inlining - for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { - CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { - CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { - CeedPragmaSIMD for (CeedInt j = 0; j < blk_size; j++) { - vv[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset] = - uu[n + k * elem_size + CeedIntMin(e + j, num_elem - 1) * elem_size * num_comp]; + switch (rstr_type) { + case CEED_RESTRICTION_STRIDED: { + // No offsets provided, Identity Restriction + bool has_backend_strides; + CeedCallBackend(CeedElemRestrictionHasBackendStrides(r, &has_backend_strides)); + if (has_backend_strides) { + // CPU backend strides are {1, elem_size, elem_size*num_comp} + // This if branch is left separate to allow better inlining + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { + CeedPragmaSIMD for (CeedInt j = 0; j < blk_size; j++) { + vv[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset] = + uu[n + k * elem_size + CeedIntMin(e + j, num_elem - 1) * elem_size * num_comp]; + } + } + } + } + } else { + // User provided strides + CeedInt strides[3]; + CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { + CeedPragmaSIMD for (CeedInt j = 0; j < blk_size; j++) { + vv[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset] = + uu[n * strides[0] + k * strides[1] + CeedIntMin(e + j, num_elem - 1) * strides[2]]; + } } } } } - } else { - // User provided strides - CeedInt strides[3]; - CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + } break; + case CEED_RESTRICTION_DEFAULT: + // Default restriction with offsets for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { - CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { - CeedPragmaSIMD for (CeedInt j = 0; j < blk_size; j++) { - vv[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset] = - uu[n * strides[0] + k * strides[1] + CeedIntMin(e + j, num_elem - 1) * strides[2]]; - } + CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * blk_size; i++) { + vv[elem_size * (k * blk_size + num_comp * e) + i - v_offset] = uu[impl->offsets[i + elem_size * e] + k * comp_stride]; } } } - } - } else { - // Offsets provided, standard or blocked restriction - // vv has shape [elem_size, num_comp, num_elem], row-major - // uu has shape [nnodes, num_comp] - for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { - CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { - CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * blk_size; i++) { - if (!is_oriented && !is_curl_oriented) { - // Unsigned restriction - vv[elem_size * (k * blk_size + num_comp * e) + i - v_offset] = uu[impl->offsets[i + elem_size * e] + k * comp_stride]; - } else if (!is_curl_oriented) { - // Signed restriction + break; + case CEED_RESTRICTION_ORIENTED: + // Restriction with orientations + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * blk_size; i++) { vv[elem_size * (k * blk_size + num_comp * e) + i - v_offset] = uu[impl->offsets[i + elem_size * e] + k * comp_stride] * (impl->orients[i + elem_size * e] ? -1.0 : 1.0); - } else { - // Restriction with tridiagonal transformation + } + } + } + break; + case CEED_RESTRICTION_CURL_ORIENTED: + // Restriction with tridiagonal transformation + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * blk_size; i++) { CeedInt ii = i % elem_size; if (ii == 0) { vv[elem_size * (k * blk_size + num_comp * e) + i - v_offset] = @@ -107,61 +120,81 @@ static inline int CeedElemRestrictionApply_Ref_Core(CeedElemRestriction r, const } } } - } + break; } } else { // Restriction from E-vector to L-vector // Performing v += r^T * u - // No offsets provided, Identity Restriction - if (!impl->offsets) { - bool has_backend_strides; - CeedCallBackend(CeedElemRestrictionHasBackendStrides(r, &has_backend_strides)); - if (has_backend_strides) { - // CPU backend strides are {1, elem_size, elem_size*num_comp} - // This if brach is left separate to allow better inlining - for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { - CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { - CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { - CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(blk_size, num_elem - e); j++) { - vv[n + k * elem_size + (e + j) * elem_size * num_comp] += - uu[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset]; + // uu has shape [elem_size, num_comp, num_elem] + // vv has shape [nnodes, num_comp] + switch (rstr_type) { + case CEED_RESTRICTION_STRIDED: { + // No offsets provided, Identity Restriction + bool has_backend_strides; + CeedCallBackend(CeedElemRestrictionHasBackendStrides(r, &has_backend_strides)); + if (has_backend_strides) { + // CPU backend strides are {1, elem_size, elem_size*num_comp} + // This if brach is left separate to allow better inlining + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { + CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(blk_size, num_elem - e); j++) { + vv[n + k * elem_size + (e + j) * elem_size * num_comp] += + uu[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset]; + } + } + } + } + } else { + // User provided strides + CeedInt strides[3]; + CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { + CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { + CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(blk_size, num_elem - e); j++) { + vv[n * strides[0] + k * strides[1] + (e + j) * strides[2]] += + uu[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset]; + } } } } } - } else { - // User provided strides - CeedInt strides[3]; - CeedCallBackend(CeedElemRestrictionGetStrides(r, &strides)); + } break; + case CEED_RESTRICTION_DEFAULT: + // Default restriction with offsets for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { - CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { - CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { - CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(blk_size, num_elem - e); j++) { - vv[n * strides[0] + k * strides[1] + (e + j) * strides[2]] += - uu[e * elem_size * num_comp + (k * elem_size + n) * blk_size + j - v_offset]; + for (CeedInt k = 0; k < num_comp; k++) { + for (CeedInt i = 0; i < elem_size * blk_size; i += blk_size) { + // Iteration bound set to discard padding elements + for (CeedInt j = i; j < i + CeedIntMin(blk_size, num_elem - e); j++) { + vv[impl->offsets[j + e * elem_size] + k * comp_stride] += uu[elem_size * (k * blk_size + num_comp * e) + j - v_offset]; } } } } - } - } else { - // Offsets provided, standard or blocked restriction - // uu has shape [elem_size, num_comp, num_elem] - // vv has shape [nnodes, num_comp] - for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { - for (CeedInt k = 0; k < num_comp; k++) { - for (CeedInt i = 0; i < elem_size * blk_size; i += blk_size) { - // Iteration bound set to discard padding elements - for (CeedInt j = i; j < i + CeedIntMin(blk_size, num_elem - e); j++) { - if (!is_oriented && !is_curl_oriented) { - // Unsigned restriction - vv[impl->offsets[j + e * elem_size] + k * comp_stride] += uu[elem_size * (k * blk_size + num_comp * e) + j - v_offset]; - } else if (!is_curl_oriented) { - // Signed restriction + break; + case CEED_RESTRICTION_ORIENTED: + // Restriction with orientations + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + for (CeedInt k = 0; k < num_comp; k++) { + for (CeedInt i = 0; i < elem_size * blk_size; i += blk_size) { + // Iteration bound set to discard padding elements + for (CeedInt j = i; j < i + CeedIntMin(blk_size, num_elem - e); j++) { vv[impl->offsets[j + e * elem_size] + k * comp_stride] += uu[elem_size * (k * blk_size + num_comp * e) + j - v_offset] * (impl->orients[j + e * elem_size] ? -1.0 : 1.0); - } else { - // Restriction with tridiagonal transformation + } + } + } + } + break; + case CEED_RESTRICTION_CURL_ORIENTED: + // Restriction with tridiagonal transformation + for (CeedInt e = start * blk_size; e < stop * blk_size; e += blk_size) { + for (CeedInt k = 0; k < num_comp; k++) { + for (CeedInt i = 0; i < elem_size * blk_size; i += blk_size) { + // Iteration bound set to discard padding elements + for (CeedInt j = i; j < i + CeedIntMin(blk_size, num_elem - e); j++) { CeedInt jj = j % elem_size; if (jj == 0) { vv[impl->offsets[j + e * elem_size] + k * comp_stride] += @@ -181,7 +214,7 @@ static inline int CeedElemRestrictionApply_Ref_Core(CeedElemRestriction r, const } } } - } + break; } } CeedCallBackend(CeedVectorRestoreArrayRead(u, &uu)); @@ -357,7 +390,8 @@ static int CeedElemRestrictionDestroy_Ref(CeedElemRestriction r) { //------------------------------------------------------------------------------ // ElemRestriction Create //------------------------------------------------------------------------------ -int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, CeedElemRestriction r) { +int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, + const CeedInt *curl_orients, CeedElemRestriction r) { CeedElemRestriction_Ref *impl; CeedInt num_elem, elem_size, num_blk, blk_size, num_comp, comp_stride; CeedCallBackend(CeedElemRestrictionGetNumElements(r, &num_elem)); @@ -373,9 +407,9 @@ int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, CeedCallBackend(CeedCalloc(1, &impl)); // Offsets data - bool is_strided; - CeedCallBackend(CeedElemRestrictionIsStrided(r, &is_strided)); - if (!is_strided) { + CeedRestrictionType rstr_type; + CeedCallBackend(CeedElemRestrictionGetType(r, &rstr_type)); + if (rstr_type != CEED_RESTRICTION_STRIDED) { // Check indices for ref or memcheck backends Ceed parent_ceed = ceed, curr_ceed = NULL; while (parent_ceed != curr_ceed) { @@ -409,6 +443,37 @@ int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, case CEED_USE_POINTER: impl->offsets = offsets; } + + // Orientation data + if (rstr_type == CEED_RESTRICTION_ORIENTED) { + switch (copy_mode) { + case CEED_COPY_VALUES: + CeedCallBackend(CeedMalloc(num_elem * elem_size, &impl->orients_allocated)); + memcpy(impl->orients_allocated, orients, num_elem * elem_size * sizeof(orients[0])); + impl->orients = impl->orients_allocated; + break; + case CEED_OWN_POINTER: + impl->orients_allocated = (bool *)orients; + impl->orients = impl->orients_allocated; + break; + case CEED_USE_POINTER: + impl->orients = orients; + } + } else if (rstr_type == CEED_RESTRICTION_CURL_ORIENTED) { + switch (copy_mode) { + case CEED_COPY_VALUES: + CeedCallBackend(CeedMalloc(num_elem * 3 * elem_size, &impl->curl_orients_allocated)); + memcpy(impl->curl_orients_allocated, curl_orients, num_elem * 3 * elem_size * sizeof(curl_orients[0])); + impl->curl_orients = impl->curl_orients_allocated; + break; + case CEED_OWN_POINTER: + impl->curl_orients_allocated = (CeedInt *)curl_orients; + impl->curl_orients = impl->curl_orients_allocated; + break; + case CEED_USE_POINTER: + impl->curl_orients = curl_orients; + } + } } CeedCallBackend(CeedElemRestrictionSetData(r, impl)); @@ -474,68 +539,3 @@ int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, } //------------------------------------------------------------------------------ -// ElemRestriction Create Oriented -//------------------------------------------------------------------------------ -int CeedElemRestrictionCreateOriented_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, - CeedElemRestriction r) { - // Set up for normal restriction with explicit offsets. This sets up dispatch to - // CeedElemRestrictionApply_Ref_* and manages the impl->offsets array copy/allocation. - CeedCallBackend(CeedElemRestrictionCreate_Ref(mem_type, copy_mode, offsets, r)); - - CeedElemRestriction_Ref *impl; - CeedInt num_elem, elem_size; - CeedCallBackend(CeedElemRestrictionGetData(r, &impl)); - CeedCallBackend(CeedElemRestrictionGetNumElements(r, &num_elem)); - CeedCallBackend(CeedElemRestrictionGetElementSize(r, &elem_size)); - - // Orientation data - switch (copy_mode) { - case CEED_COPY_VALUES: - CeedCallBackend(CeedMalloc(num_elem * elem_size, &impl->orients_allocated)); - memcpy(impl->orients_allocated, orients, num_elem * elem_size * sizeof(orients[0])); - impl->orients = impl->orients_allocated; - break; - case CEED_OWN_POINTER: - impl->orients_allocated = (bool *)orients; - impl->orients = impl->orients_allocated; - break; - case CEED_USE_POINTER: - impl->orients = orients; - } - - return CEED_ERROR_SUCCESS; -} - -//------------------------------------------------------------------------------ -// ElemRestriction Create Curl-Conforming Oriented -//------------------------------------------------------------------------------ -int CeedElemRestrictionCreateCurlOriented_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const CeedInt *curl_orients, - CeedElemRestriction r) { - // Set up for normal restriction with explicit offsets. This sets up dispatch to - // CeedElemRestrictionApply_Ref_* and manages the impl->offsets array copy/allocation. - CeedCallBackend(CeedElemRestrictionCreate_Ref(mem_type, copy_mode, offsets, r)); - - CeedElemRestriction_Ref *impl; - CeedInt num_elem, elem_size; - CeedCallBackend(CeedElemRestrictionGetData(r, &impl)); - CeedCallBackend(CeedElemRestrictionGetNumElements(r, &num_elem)); - CeedCallBackend(CeedElemRestrictionGetElementSize(r, &elem_size)); - - // Orientation data - switch (copy_mode) { - case CEED_COPY_VALUES: - CeedCallBackend(CeedMalloc(num_elem * 3 * elem_size, &impl->curl_orients_allocated)); - memcpy(impl->curl_orients_allocated, curl_orients, num_elem * 3 * elem_size * sizeof(curl_orients[0])); - impl->curl_orients = impl->curl_orients_allocated; - break; - case CEED_OWN_POINTER: - impl->curl_orients_allocated = (CeedInt *)curl_orients; - impl->curl_orients = impl->curl_orients_allocated; - break; - case CEED_USE_POINTER: - impl->curl_orients = curl_orients; - } - - return CEED_ERROR_SUCCESS; -} -//------------------------------------------------------------------------------ diff --git a/backends/ref/ceed-ref.c b/backends/ref/ceed-ref.c index 1c21825968..4e5e4f540e 100644 --- a/backends/ref/ceed-ref.c +++ b/backends/ref/ceed-ref.c @@ -26,11 +26,7 @@ static int CeedInit_Ref(const char *resource, Ceed ceed) { CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "BasisCreateHcurl", CeedBasisCreateHcurl_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "TensorContractCreate", CeedTensorContractCreate_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreate", CeedElemRestrictionCreate_Ref)); - CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreateOriented", CeedElemRestrictionCreateOriented_Ref)); - CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreateCurlOriented", CeedElemRestrictionCreateCurlOriented_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreateBlocked", CeedElemRestrictionCreate_Ref)); - CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreateBlockedOriented", CeedElemRestrictionCreateOriented_Ref)); - CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "ElemRestrictionCreateBlockedCurlOriented", CeedElemRestrictionCreateCurlOriented_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "QFunctionCreate", CeedQFunctionCreate_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "QFunctionContextCreate", CeedQFunctionContextCreate_Ref)); CeedCallBackend(CeedSetBackendFunction(ceed, "Ceed", ceed, "OperatorCreate", CeedOperatorCreate_Ref)); diff --git a/backends/ref/ceed-ref.h b/backends/ref/ceed-ref.h index 5b88d1eab8..869ec20731 100644 --- a/backends/ref/ceed-ref.h +++ b/backends/ref/ceed-ref.h @@ -61,11 +61,8 @@ typedef struct { CEED_INTERN int CeedVectorCreate_Ref(CeedSize n, CeedVector vec); -CEED_INTERN int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, CeedElemRestriction r); -CEED_INTERN int CeedElemRestrictionCreateOriented_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, - CeedElemRestriction r); -CEED_INTERN int CeedElemRestrictionCreateCurlOriented_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, - const CeedInt *curl_orients, CeedElemRestriction r); +CEED_INTERN int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, + const CeedInt *curl_orients, CeedElemRestriction r); CEED_INTERN int CeedBasisCreateTensorH1_Ref(CeedInt dim, CeedInt P_1d, CeedInt Q_1d, const CeedScalar *interp_1d, const CeedScalar *grad_1d, const CeedScalar *q_ref_1d, const CeedScalar *q_weight_1d, CeedBasis basis); diff --git a/include/ceed-impl.h b/include/ceed-impl.h index dfb3d83742..75af1592d7 100644 --- a/include/ceed-impl.h +++ b/include/ceed-impl.h @@ -97,12 +97,8 @@ struct Ceed_private { int (*GetPreferredMemType)(CeedMemType *); int (*Destroy)(Ceed); int (*VectorCreate)(CeedSize, CeedVector); - int (*ElemRestrictionCreate)(CeedMemType, CeedCopyMode, const CeedInt *, CeedElemRestriction); - int (*ElemRestrictionCreateOriented)(CeedMemType, CeedCopyMode, const CeedInt *, const bool *, CeedElemRestriction); - int (*ElemRestrictionCreateCurlOriented)(CeedMemType, CeedCopyMode, const CeedInt *, const CeedInt *, CeedElemRestriction); - int (*ElemRestrictionCreateBlocked)(CeedMemType, CeedCopyMode, const CeedInt *, CeedElemRestriction); - int (*ElemRestrictionCreateBlockedOriented)(CeedMemType, CeedCopyMode, const CeedInt *, const bool *, CeedElemRestriction); - int (*ElemRestrictionCreateBlockedCurlOriented)(CeedMemType, CeedCopyMode, const CeedInt *, const CeedInt *, CeedElemRestriction); + int (*ElemRestrictionCreate)(CeedMemType, CeedCopyMode, const CeedInt *, const bool *, const CeedInt *, CeedElemRestriction); + int (*ElemRestrictionCreateBlocked)(CeedMemType, CeedCopyMode, const CeedInt *, const bool *, const CeedInt *, CeedElemRestriction); int (*BasisCreateTensorH1)(CeedInt, CeedInt, CeedInt, const CeedScalar *, const CeedScalar *, const CeedScalar *, const CeedScalar *, CeedBasis); int (*BasisCreateH1)(CeedElemTopology, CeedInt, CeedInt, CeedInt, const CeedScalar *, const CeedScalar *, const CeedScalar *, const CeedScalar *, CeedBasis); @@ -169,6 +165,8 @@ struct CeedElemRestriction_private { CeedInt num_blk; /* number of blocks of elements */ CeedInt *strides; /* strides between [nodes, components, elements] */ CeedInt layout[3]; /* E-vector layout [nodes, components, elements] */ + CeedRestrictionType + rstr_type; /* initialized in element restriction constructor for default, oriented, curl-oriented, or strided element restriction */ uint64_t num_readers; /* number of instances of offset read only access */ void *data; /* place for the backend to store any data */ }; diff --git a/include/ceed/backend.h b/include/ceed/backend.h index 04976b3dd1..7c3c154bd6 100644 --- a/include/ceed/backend.h +++ b/include/ceed/backend.h @@ -164,6 +164,20 @@ CEED_EXTERN int CeedVectorGetData(CeedVector vec, void *data); CEED_EXTERN int CeedVectorSetData(CeedVector vec, void *data); CEED_EXTERN int CeedVectorReference(CeedVector vec); +/// Type of element restriction; +/// @ingroup CeedElemRestriction +typedef enum { + /// Default element restriction with offsets + CEED_RESTRICTION_DEFAULT = 1, + /// Oriented element restriction + CEED_RESTRICTION_ORIENTED = 2, + /// Curl-oriented element restriction + CEED_RESTRICTION_CURL_ORIENTED = 3, + /// Strided element restriction + CEED_RESTRICTION_STRIDED = 4, +} CeedRestrictionType; + +CEED_EXTERN int CeedElemRestrictionGetType(CeedElemRestriction rstr, CeedRestrictionType *rstr_type); CEED_EXTERN int CeedElemRestrictionIsStrided(CeedElemRestriction rstr, bool *is_strided); CEED_EXTERN int CeedElemRestrictionGetStrides(CeedElemRestriction rstr, CeedInt (*strides)[3]); CEED_EXTERN int CeedElemRestrictionHasBackendStrides(CeedElemRestriction rstr, bool *has_backend_strides); diff --git a/interface/ceed-elemrestriction.c b/interface/ceed-elemrestriction.c index dae29651aa..bcb65bba4b 100644 --- a/interface/ceed-elemrestriction.c +++ b/interface/ceed-elemrestriction.c @@ -80,18 +80,33 @@ int CeedPermutePadOrientations(const bool *orients, bool *blk_orients, CeedInt n /// @addtogroup CeedElemRestrictionBackend /// @{ +/** + @brief Get the type of a CeedElemRestriction + + @param[in] rstr CeedElemRestriction + @param[out] rstr_type Variable to store restriction type + + @return An error code: 0 - success, otherwise - failure + + @ref Backend +**/ +int CeedElemRestrictionGetType(CeedElemRestriction rstr, CeedRestrictionType *rstr_type) { + *rstr_type = rstr->rstr_type; + return CEED_ERROR_SUCCESS; +} + /** @brief Get the strided status of a CeedElemRestriction - @param[in] rstr CeedElemRestriction - @param[out] is_strided Variable to store strided status, 1 if strided else 0 + @param[in] rstr CeedElemRestriction + @param[out] is_strided Variable to store strided status, 1 if strided else 0 @return An error code: 0 - success, otherwise - failure @ref Backend **/ int CeedElemRestrictionIsStrided(CeedElemRestriction rstr, bool *is_strided) { - *is_strided = rstr->strides ? true : false; + *is_strided = (rstr->rstr_type == CEED_RESTRICTION_STRIDED); return CEED_ERROR_SUCCESS; } @@ -417,7 +432,8 @@ int CeedElemRestrictionCreate(Ceed ceed, CeedInt num_elem, CeedInt elem_size, Ce (*rstr)->l_size = l_size; (*rstr)->num_blk = num_elem; (*rstr)->blk_size = 1; - CeedCall(ceed->ElemRestrictionCreate(mem_type, copy_mode, offsets, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_DEFAULT; + CeedCall(ceed->ElemRestrictionCreate(mem_type, copy_mode, offsets, NULL, NULL, *rstr)); return CEED_ERROR_SUCCESS; } @@ -445,11 +461,11 @@ unknowns corresponding to element i, where 0 <= i < @a num_elem. All offsets mus int CeedElemRestrictionCreateOriented(Ceed ceed, CeedInt num_elem, CeedInt elem_size, CeedInt num_comp, CeedInt comp_stride, CeedSize l_size, CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, CeedElemRestriction *rstr) { - if (!ceed->ElemRestrictionCreateOriented) { + if (!ceed->ElemRestrictionCreate) { Ceed delegate; CeedCall(CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction")); - CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not implement ElemRestrictionCreateOriented"); + CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not implement ElemRestrictionCreate"); CeedCall( CeedElemRestrictionCreateOriented(delegate, num_elem, elem_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, orients, rstr)); return CEED_ERROR_SUCCESS; @@ -470,7 +486,8 @@ int CeedElemRestrictionCreateOriented(Ceed ceed, CeedInt num_elem, CeedInt elem_ (*rstr)->l_size = l_size; (*rstr)->num_blk = num_elem; (*rstr)->blk_size = 1; - CeedCall(ceed->ElemRestrictionCreateOriented(mem_type, copy_mode, offsets, orients, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_ORIENTED; + CeedCall(ceed->ElemRestrictionCreate(mem_type, copy_mode, offsets, orients, NULL, *rstr)); return CEED_ERROR_SUCCESS; } @@ -501,11 +518,11 @@ to resolve face orientation issues for 3D meshes (https://dl.acm.org/doi/pdf/10. int CeedElemRestrictionCreateCurlOriented(Ceed ceed, CeedInt num_elem, CeedInt elem_size, CeedInt num_comp, CeedInt comp_stride, CeedSize l_size, CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const CeedInt *curl_orients, CeedElemRestriction *rstr) { - if (!ceed->ElemRestrictionCreateCurlOriented) { + if (!ceed->ElemRestrictionCreate) { Ceed delegate; CeedCall(CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction")); - CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not implement ElemRestrictionCreateCurlOriented"); + CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not implement ElemRestrictionCreate"); CeedCall(CeedElemRestrictionCreateCurlOriented(delegate, num_elem, elem_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, curl_orients, rstr)); return CEED_ERROR_SUCCESS; @@ -526,7 +543,8 @@ int CeedElemRestrictionCreateCurlOriented(Ceed ceed, CeedInt num_elem, CeedInt e (*rstr)->l_size = l_size; (*rstr)->num_blk = num_elem; (*rstr)->blk_size = 1; - CeedCall(ceed->ElemRestrictionCreateCurlOriented(mem_type, copy_mode, offsets, curl_orients, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_CURL_ORIENTED; + CeedCall(ceed->ElemRestrictionCreate(mem_type, copy_mode, offsets, NULL, curl_orients, *rstr)); return CEED_ERROR_SUCCESS; } @@ -569,9 +587,10 @@ int CeedElemRestrictionCreateStrided(Ceed ceed, CeedInt num_elem, CeedInt elem_s (*rstr)->l_size = l_size; (*rstr)->num_blk = num_elem; (*rstr)->blk_size = 1; + (*rstr)->rstr_type = CEED_RESTRICTION_STRIDED; CeedCall(CeedMalloc(3, &(*rstr)->strides)); for (CeedInt i = 0; i < 3; i++) (*rstr)->strides[i] = strides[i]; - CeedCall(ceed->ElemRestrictionCreate(CEED_MEM_HOST, CEED_OWN_POINTER, NULL, *rstr)); + CeedCall(ceed->ElemRestrictionCreate(CEED_MEM_HOST, CEED_OWN_POINTER, NULL, NULL, NULL, *rstr)); return CEED_ERROR_SUCCESS; } @@ -633,7 +652,8 @@ int CeedElemRestrictionCreateBlocked(Ceed ceed, CeedInt num_elem, CeedInt elem_s (*rstr)->l_size = l_size; (*rstr)->num_blk = num_blk; (*rstr)->blk_size = blk_size; - CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_DEFAULT; + CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, NULL, NULL, *rstr)); if (copy_mode == CEED_OWN_POINTER) { CeedCall(CeedFree(&offsets)); } @@ -671,11 +691,11 @@ int CeedElemRestrictionCreateBlockedOriented(Ceed ceed, CeedInt num_elem, CeedIn bool *blk_orients; CeedInt num_blk = (num_elem / blk_size) + !!(num_elem % blk_size); - if (!ceed->ElemRestrictionCreateBlockedOriented) { + if (!ceed->ElemRestrictionCreateBlocked) { Ceed delegate; CeedCall(CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction")); - CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not support ElemRestrictionCreateBlockedOriented"); + CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not support ElemRestrictionCreateBlocked"); CeedCall(CeedElemRestrictionCreateBlockedOriented(delegate, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, orients, rstr)); return CEED_ERROR_SUCCESS; @@ -703,8 +723,8 @@ int CeedElemRestrictionCreateBlockedOriented(Ceed ceed, CeedInt num_elem, CeedIn (*rstr)->l_size = l_size; (*rstr)->num_blk = num_blk; (*rstr)->blk_size = blk_size; - CeedCall( - ceed->ElemRestrictionCreateBlockedOriented(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, (const bool *)blk_orients, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_ORIENTED; + CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, (const bool *)blk_orients, NULL, *rstr)); if (copy_mode == CEED_OWN_POINTER) { CeedCall(CeedFree(&offsets)); } @@ -745,11 +765,11 @@ int CeedElemRestrictionCreateBlockedCurlOriented(Ceed ceed, CeedInt num_elem, Ce CeedInt *blk_curl_orients; CeedInt num_blk = (num_elem / blk_size) + !!(num_elem % blk_size); - if (!ceed->ElemRestrictionCreateBlockedCurlOriented) { + if (!ceed->ElemRestrictionCreateBlocked) { Ceed delegate; CeedCall(CeedGetObjectDelegate(ceed, &delegate, "ElemRestriction")); - CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not support ElemRestrictionCreateBlockedCurlOriented"); + CeedCheck(delegate, ceed, CEED_ERROR_UNSUPPORTED, "Backend does not support ElemRestrictionCreateBlocked"); CeedCall(CeedElemRestrictionCreateBlockedCurlOriented(delegate, num_elem, elem_size, blk_size, num_comp, comp_stride, l_size, mem_type, copy_mode, offsets, curl_orients, rstr)); return CEED_ERROR_SUCCESS; @@ -777,8 +797,9 @@ int CeedElemRestrictionCreateBlockedCurlOriented(Ceed ceed, CeedInt num_elem, Ce (*rstr)->l_size = l_size; (*rstr)->num_blk = num_blk; (*rstr)->blk_size = blk_size; - CeedCall(ceed->ElemRestrictionCreateBlockedCurlOriented(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, - (const CeedInt *)blk_curl_orients, *rstr)); + (*rstr)->rstr_type = CEED_RESTRICTION_CURL_ORIENTED; + CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, (const CeedInt *)blk_offsets, NULL, (const CeedInt *)blk_curl_orients, + *rstr)); if (copy_mode == CEED_OWN_POINTER) { CeedCall(CeedFree(&offsets)); } @@ -830,9 +851,10 @@ int CeedElemRestrictionCreateBlockedStrided(Ceed ceed, CeedInt num_elem, CeedInt (*rstr)->l_size = l_size; (*rstr)->num_blk = num_blk; (*rstr)->blk_size = blk_size; + (*rstr)->rstr_type = CEED_RESTRICTION_STRIDED; CeedCall(CeedMalloc(3, &(*rstr)->strides)); for (CeedInt i = 0; i < 3; i++) (*rstr)->strides[i] = strides[i]; - CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, NULL, *rstr)); + CeedCall(ceed->ElemRestrictionCreateBlocked(CEED_MEM_HOST, CEED_OWN_POINTER, NULL, NULL, NULL, *rstr)); return CEED_ERROR_SUCCESS; } diff --git a/interface/ceed.c b/interface/ceed.c index eebc46925e..f34ab10d3f 100644 --- a/interface/ceed.c +++ b/interface/ceed.c @@ -785,11 +785,7 @@ int CeedInit(const char *resource, Ceed *ceed) { CEED_FTABLE_ENTRY(Ceed, Destroy), CEED_FTABLE_ENTRY(Ceed, VectorCreate), CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreate), - CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreateOriented), - CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreateCurlOriented), CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreateBlocked), - CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreateBlockedOriented), - CEED_FTABLE_ENTRY(Ceed, ElemRestrictionCreateBlockedCurlOriented), CEED_FTABLE_ENTRY(Ceed, BasisCreateTensorH1), CEED_FTABLE_ENTRY(Ceed, BasisCreateH1), CEED_FTABLE_ENTRY(Ceed, BasisCreateHdiv),