1 // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors. 2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details. 3 // 4 // SPDX-License-Identifier: BSD-2-Clause 5 // 6 // This file is part of CEED: http://github.com/ceed 7 8 #include <ceed.h> 9 #include <ceed/backend.h> 10 #include <stdbool.h> 11 #include <stdlib.h> 12 #include <string.h> 13 14 #include "ceed-ref.h" 15 16 //------------------------------------------------------------------------------ 17 // Core ElemRestriction Apply Code 18 //------------------------------------------------------------------------------ 19 static inline int CeedElemRestrictionApplyStridedNoTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 20 CeedInt start, CeedInt stop, CeedInt num_elem, CeedInt elem_size, 21 CeedInt v_offset, const CeedScalar *__restrict__ uu, 22 CeedScalar *__restrict__ vv) { 23 // No offsets provided, identity restriction 24 bool has_backend_strides; 25 26 CeedCallBackend(CeedElemRestrictionHasBackendStrides(rstr, &has_backend_strides)); 27 if (has_backend_strides) { 28 // CPU backend strides are {1, elem_size, elem_size*num_comp} 29 // This if branch is left separate to allow better inlining 30 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 31 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 32 CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { 33 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 34 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 35 uu[n + k * elem_size + CeedIntMin(e + j, num_elem - 1) * elem_size * num_comp]; 36 } 37 } 38 } 39 } 40 } else { 41 // User provided strides 42 CeedInt strides[3]; 43 44 CeedCallBackend(CeedElemRestrictionGetStrides(rstr, &strides)); 45 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 46 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 47 CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { 48 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 49 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 50 uu[n * strides[0] + k * strides[1] + CeedIntMin(e + j, num_elem - 1) * strides[2]]; 51 } 52 } 53 } 54 } 55 } 56 return CEED_ERROR_SUCCESS; 57 } 58 59 static inline int CeedElemRestrictionApplyOffsetNoTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 60 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 61 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 62 CeedScalar *__restrict__ vv) { 63 // Default restriction with offsets 64 CeedElemRestriction_Ref *impl; 65 66 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 67 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 68 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 69 CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * block_size; i++) { 70 vv[elem_size * (k * block_size + e * num_comp) + i - v_offset] = uu[impl->offsets[i + e * elem_size] + k * comp_stride]; 71 } 72 } 73 } 74 return CEED_ERROR_SUCCESS; 75 } 76 77 static inline int CeedElemRestrictionApplyOrientedNoTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 78 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 79 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 80 CeedScalar *__restrict__ vv) { 81 // Restriction with orientations 82 CeedElemRestriction_Ref *impl; 83 84 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 85 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 86 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 87 CeedPragmaSIMD for (CeedInt i = 0; i < elem_size * block_size; i++) { 88 vv[elem_size * (k * block_size + e * num_comp) + i - v_offset] = 89 uu[impl->offsets[i + e * elem_size] + k * comp_stride] * (impl->orients[i + e * elem_size] ? -1.0 : 1.0); 90 } 91 } 92 } 93 return CEED_ERROR_SUCCESS; 94 } 95 96 static inline int CeedElemRestrictionApplyCurlOrientedNoTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 97 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 98 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 99 CeedScalar *__restrict__ vv) { 100 // Restriction with tridiagonal transformation 101 CeedElemRestriction_Ref *impl; 102 103 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 104 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 105 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 106 CeedInt n = 0; 107 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 108 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 109 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 110 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size] + 111 uu[impl->offsets[j + (n + 1) * block_size + e * elem_size] + k * comp_stride] * 112 impl->curl_orients[j + (3 * n + 2) * block_size + e * 3 * elem_size]; 113 } 114 for (n = 1; n < elem_size - 1; n++) { 115 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 116 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 117 uu[impl->offsets[j + (n - 1) * block_size + e * elem_size] + k * comp_stride] * 118 impl->curl_orients[j + (3 * n + 0) * block_size + e * 3 * elem_size] + 119 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 120 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size] + 121 uu[impl->offsets[j + (n + 1) * block_size + e * elem_size] + k * comp_stride] * 122 impl->curl_orients[j + (3 * n + 2) * block_size + e * 3 * elem_size]; 123 } 124 } 125 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 126 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 127 uu[impl->offsets[j + (n - 1) * block_size + e * elem_size] + k * comp_stride] * 128 impl->curl_orients[j + (3 * n + 0) * block_size + e * 3 * elem_size] + 129 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 130 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]; 131 } 132 } 133 } 134 return CEED_ERROR_SUCCESS; 135 } 136 137 static inline int CeedElemRestrictionApplyCurlOrientedUnsignedNoTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, 138 const CeedInt block_size, const CeedInt comp_stride, CeedInt start, 139 CeedInt stop, CeedInt num_elem, CeedInt elem_size, 140 CeedInt v_offset, const CeedScalar *__restrict__ uu, 141 CeedScalar *__restrict__ vv) { 142 // Restriction with (unsigned) tridiagonal transformation 143 CeedElemRestriction_Ref *impl; 144 145 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 146 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 147 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 148 CeedInt n = 0; 149 150 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 151 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 152 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 153 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]) + 154 uu[impl->offsets[j + (n + 1) * block_size + e * elem_size] + k * comp_stride] * 155 abs(impl->curl_orients[j + (3 * n + 2) * block_size + e * 3 * elem_size]); 156 } 157 for (n = 1; n < elem_size - 1; n++) { 158 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 159 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 160 uu[impl->offsets[j + (n - 1) * block_size + e * elem_size] + k * comp_stride] * 161 abs(impl->curl_orients[j + (3 * n + 0) * block_size + e * 3 * elem_size]) + 162 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 163 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]) + 164 uu[impl->offsets[j + (n + 1) * block_size + e * elem_size] + k * comp_stride] * 165 abs(impl->curl_orients[j + (3 * n + 2) * block_size + e * 3 * elem_size]); 166 } 167 } 168 CeedPragmaSIMD for (CeedInt j = 0; j < block_size; j++) { 169 vv[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] = 170 uu[impl->offsets[j + (n - 1) * block_size + e * elem_size] + k * comp_stride] * 171 abs(impl->curl_orients[j + (3 * n + 0) * block_size + e * 3 * elem_size]) + 172 uu[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] * 173 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]); 174 } 175 } 176 } 177 return CEED_ERROR_SUCCESS; 178 } 179 180 static inline int CeedElemRestrictionApplyStridedTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 181 CeedInt start, CeedInt stop, CeedInt num_elem, CeedInt elem_size, 182 CeedInt v_offset, const CeedScalar *__restrict__ uu, 183 CeedScalar *__restrict__ vv) { 184 // No offsets provided, identity restriction 185 bool has_backend_strides; 186 187 CeedCallBackend(CeedElemRestrictionHasBackendStrides(rstr, &has_backend_strides)); 188 if (has_backend_strides) { 189 // CPU backend strides are {1, elem_size, elem_size*num_comp} 190 // This if brach is left separate to allow better inlining 191 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 192 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 193 CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { 194 CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(block_size, num_elem - e); j++) { 195 vv[n + k * elem_size + (e + j) * elem_size * num_comp] += uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset]; 196 } 197 } 198 } 199 } 200 } else { 201 // User provided strides 202 CeedInt strides[3]; 203 204 CeedCallBackend(CeedElemRestrictionGetStrides(rstr, &strides)); 205 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 206 CeedPragmaSIMD for (CeedInt k = 0; k < num_comp; k++) { 207 CeedPragmaSIMD for (CeedInt n = 0; n < elem_size; n++) { 208 CeedPragmaSIMD for (CeedInt j = 0; j < CeedIntMin(block_size, num_elem - e); j++) { 209 vv[n * strides[0] + k * strides[1] + (e + j) * strides[2]] += 210 uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset]; 211 } 212 } 213 } 214 } 215 } 216 return CEED_ERROR_SUCCESS; 217 } 218 219 static inline int CeedElemRestrictionApplyOffsetTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 220 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 221 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 222 CeedScalar *__restrict__ vv) { 223 // Default restriction with offsets 224 CeedElemRestriction_Ref *impl; 225 226 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 227 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 228 for (CeedInt k = 0; k < num_comp; k++) { 229 for (CeedInt i = 0; i < elem_size * block_size; i += block_size) { 230 // Iteration bound set to discard padding elements 231 for (CeedInt j = i; j < i + CeedIntMin(block_size, num_elem - e); j++) { 232 CeedScalar uu_val; 233 234 uu_val = uu[elem_size * (k * block_size + e * num_comp) + j - v_offset]; 235 CeedPragmaAtomic vv[impl->offsets[j + e * elem_size] + k * comp_stride] += uu_val; 236 } 237 } 238 } 239 } 240 return CEED_ERROR_SUCCESS; 241 } 242 243 static inline int CeedElemRestrictionApplyOrientedTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 244 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 245 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 246 CeedScalar *__restrict__ vv) { 247 // Restriction with orientations 248 CeedElemRestriction_Ref *impl; 249 250 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 251 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 252 for (CeedInt k = 0; k < num_comp; k++) { 253 for (CeedInt i = 0; i < elem_size * block_size; i += block_size) { 254 // Iteration bound set to discard padding elements 255 for (CeedInt j = i; j < i + CeedIntMin(block_size, num_elem - e); j++) { 256 CeedScalar uu_val; 257 258 uu_val = uu[elem_size * (k * block_size + e * num_comp) + j - v_offset] * (impl->orients[j + e * elem_size] ? -1.0 : 1.0); 259 CeedPragmaAtomic vv[impl->offsets[j + e * elem_size] + k * comp_stride] += uu_val; 260 } 261 } 262 } 263 } 264 return CEED_ERROR_SUCCESS; 265 } 266 267 static inline int CeedElemRestrictionApplyCurlOrientedTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 268 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedInt num_elem, 269 CeedInt elem_size, CeedInt v_offset, const CeedScalar *__restrict__ uu, 270 CeedScalar *__restrict__ vv) { 271 // Restriction with tridiagonal transformation 272 CeedElemRestriction_Ref *impl; 273 274 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 275 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 276 for (CeedInt k = 0; k < num_comp; k++) { 277 // Iteration bound set to discard padding elements 278 const CeedInt block_end = CeedIntMin(block_size, num_elem - e); 279 CeedInt n = 0; 280 281 for (CeedInt j = 0; j < block_end; j++) { 282 CeedScalar uu_val; 283 284 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 285 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size] + 286 uu[e * elem_size * num_comp + (k * elem_size + n + 1) * block_size + j - v_offset] * 287 impl->curl_orients[j + (3 * n + 3) * block_size + e * 3 * elem_size]; 288 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 289 } 290 for (n = 1; n < elem_size - 1; n++) { 291 for (CeedInt j = 0; j < block_end; j++) { 292 CeedScalar uu_val; 293 294 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n - 1) * block_size + j - v_offset] * 295 impl->curl_orients[j + (3 * n - 1) * block_size + e * 3 * elem_size] + 296 uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 297 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size] + 298 uu[e * elem_size * num_comp + (k * elem_size + n + 1) * block_size + j - v_offset] * 299 impl->curl_orients[j + (3 * n + 3) * block_size + e * 3 * elem_size]; 300 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 301 } 302 } 303 for (CeedInt j = 0; j < block_end; j++) { 304 CeedScalar uu_val; 305 306 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n - 1) * block_size + j - v_offset] * 307 impl->curl_orients[j + (3 * n - 1) * block_size + e * 3 * elem_size] + 308 uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 309 impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]; 310 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 311 } 312 } 313 } 314 return CEED_ERROR_SUCCESS; 315 } 316 317 static inline int CeedElemRestrictionApplyCurlOrientedUnsignedTranspose_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, 318 const CeedInt block_size, const CeedInt comp_stride, CeedInt start, 319 CeedInt stop, CeedInt num_elem, CeedInt elem_size, CeedInt v_offset, 320 const CeedScalar *__restrict__ uu, CeedScalar *__restrict__ vv) { 321 // Restriction with (unsigned) tridiagonal transformation 322 CeedElemRestriction_Ref *impl; 323 324 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 325 for (CeedInt e = start * block_size; e < stop * block_size; e += block_size) { 326 for (CeedInt k = 0; k < num_comp; k++) { 327 // Iteration bound set to discard padding elements 328 const CeedInt block_end = CeedIntMin(block_size, num_elem - e); 329 CeedInt n = 0; 330 331 for (CeedInt j = 0; j < block_end; j++) { 332 CeedScalar uu_val; 333 334 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 335 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]) + 336 uu[e * elem_size * num_comp + (k * elem_size + n + 1) * block_size + j - v_offset] * 337 abs(impl->curl_orients[j + (3 * n + 3) * block_size + e * 3 * elem_size]); 338 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 339 } 340 for (n = 1; n < elem_size - 1; n++) { 341 for (CeedInt j = 0; j < block_end; j++) { 342 CeedScalar uu_val; 343 344 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n - 1) * block_size + j - v_offset] * 345 abs(impl->curl_orients[j + (3 * n - 1) * block_size + e * 3 * elem_size]) + 346 uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 347 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]) + 348 uu[e * elem_size * num_comp + (k * elem_size + n + 1) * block_size + j - v_offset] * 349 abs(impl->curl_orients[j + (3 * n + 3) * block_size + e * 3 * elem_size]); 350 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 351 } 352 } 353 for (CeedInt j = 0; j < block_end; j++) { 354 CeedScalar uu_val; 355 356 uu_val = uu[e * elem_size * num_comp + (k * elem_size + n - 1) * block_size + j - v_offset] * 357 abs(impl->curl_orients[j + (3 * n - 1) * block_size + e * 3 * elem_size]) + 358 uu[e * elem_size * num_comp + (k * elem_size + n) * block_size + j - v_offset] * 359 abs(impl->curl_orients[j + (3 * n + 1) * block_size + e * 3 * elem_size]); 360 CeedPragmaAtomic vv[impl->offsets[j + n * block_size + e * elem_size] + k * comp_stride] += uu_val; 361 } 362 } 363 } 364 return CEED_ERROR_SUCCESS; 365 } 366 367 static inline int CeedElemRestrictionApplyAtPointsInElement_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, CeedInt start, CeedInt stop, 368 CeedTransposeMode t_mode, const CeedScalar *__restrict__ uu, 369 CeedScalar *__restrict__ vv) { 370 CeedInt num_points, l_vec_offset, e_vec_offset = 0; 371 CeedElemRestriction_Ref *impl; 372 373 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 374 for (CeedInt e = start; e < stop; e++) { 375 l_vec_offset = impl->offsets[e]; 376 CeedCallBackend(CeedElemRestrictionGetNumPointsInElement(rstr, e, &num_points)); 377 if (t_mode == CEED_NOTRANSPOSE) { 378 for (CeedInt i = 0; i < num_points; i++) { 379 for (CeedInt j = 0; j < num_comp; j++) vv[j * num_points + i + e_vec_offset] = uu[impl->offsets[i + l_vec_offset] * num_comp + j]; 380 } 381 } else { 382 for (CeedInt i = 0; i < num_points; i++) { 383 for (CeedInt j = 0; j < num_comp; j++) vv[impl->offsets[i + l_vec_offset] * num_comp + j] = uu[j * num_points + i + e_vec_offset]; 384 } 385 } 386 e_vec_offset += num_points * num_comp; 387 } 388 return CEED_ERROR_SUCCESS; 389 } 390 391 static inline int CeedElemRestrictionApply_Ref_Core(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, 392 const CeedInt comp_stride, CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, 393 bool use_orients, CeedVector u, CeedVector v, CeedRequest *request) { 394 CeedInt num_elem, elem_size, v_offset; 395 CeedRestrictionType rstr_type; 396 const CeedScalar *uu; 397 CeedScalar *vv; 398 399 CeedCallBackend(CeedElemRestrictionGetNumElements(rstr, &num_elem)); 400 CeedCallBackend(CeedElemRestrictionGetElementSize(rstr, &elem_size)); 401 v_offset = start * block_size * elem_size * num_comp; 402 CeedCallBackend(CeedElemRestrictionGetType(rstr, &rstr_type)); 403 CeedCallBackend(CeedVectorGetArrayRead(u, CEED_MEM_HOST, &uu)); 404 405 if (t_mode == CEED_TRANSPOSE) { 406 // Sum into for transpose mode, E-vector to L-vector 407 CeedCallBackend(CeedVectorGetArray(v, CEED_MEM_HOST, &vv)); 408 } else { 409 // Overwrite for notranspose mode, L-vector to E-vector 410 CeedCallBackend(CeedVectorGetArrayWrite(v, CEED_MEM_HOST, &vv)); 411 } 412 if (t_mode == CEED_TRANSPOSE) { 413 // Restriction from E-vector to L-vector 414 // Performing v += r^T * u 415 // uu has shape [elem_size, num_comp, num_elem], row-major 416 // vv has shape [nnodes, num_comp] 417 // Sum into for transpose mode 418 switch (rstr_type) { 419 case CEED_RESTRICTION_STRIDED: 420 CeedCallBackend( 421 CeedElemRestrictionApplyStridedTranspose_Ref_Core(rstr, num_comp, block_size, start, stop, num_elem, elem_size, v_offset, uu, vv)); 422 break; 423 case CEED_RESTRICTION_STANDARD: 424 CeedCallBackend(CeedElemRestrictionApplyOffsetTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, elem_size, 425 v_offset, uu, vv)); 426 break; 427 case CEED_RESTRICTION_ORIENTED: 428 if (use_signs) { 429 CeedCallBackend(CeedElemRestrictionApplyOrientedTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 430 elem_size, v_offset, uu, vv)); 431 } else { 432 CeedCallBackend(CeedElemRestrictionApplyOffsetTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, elem_size, 433 v_offset, uu, vv)); 434 } 435 break; 436 case CEED_RESTRICTION_CURL_ORIENTED: 437 if (use_signs && use_orients) { 438 CeedCallBackend(CeedElemRestrictionApplyCurlOrientedTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 439 elem_size, v_offset, uu, vv)); 440 } else if (use_orients) { 441 CeedCallBackend(CeedElemRestrictionApplyCurlOrientedUnsignedTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, 442 num_elem, elem_size, v_offset, uu, vv)); 443 } else { 444 CeedCallBackend(CeedElemRestrictionApplyOffsetTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, elem_size, 445 v_offset, uu, vv)); 446 } 447 break; 448 case CEED_RESTRICTION_POINTS: 449 CeedCallBackend(CeedElemRestrictionApplyAtPointsInElement_Ref_Core(rstr, num_comp, start, stop, t_mode, uu, vv)); 450 break; 451 } 452 } else { 453 // Restriction from L-vector to E-vector 454 // Perform: v = r * u 455 // vv has shape [elem_size, num_comp, num_elem], row-major 456 // uu has shape [nnodes, num_comp] 457 // Overwrite for notranspose mode 458 switch (rstr_type) { 459 case CEED_RESTRICTION_STRIDED: 460 CeedCallBackend( 461 CeedElemRestrictionApplyStridedNoTranspose_Ref_Core(rstr, num_comp, block_size, start, stop, num_elem, elem_size, v_offset, uu, vv)); 462 break; 463 case CEED_RESTRICTION_STANDARD: 464 CeedCallBackend(CeedElemRestrictionApplyOffsetNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, elem_size, 465 v_offset, uu, vv)); 466 break; 467 case CEED_RESTRICTION_ORIENTED: 468 if (use_signs) { 469 CeedCallBackend(CeedElemRestrictionApplyOrientedNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 470 elem_size, v_offset, uu, vv)); 471 } else { 472 CeedCallBackend(CeedElemRestrictionApplyOffsetNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 473 elem_size, v_offset, uu, vv)); 474 } 475 break; 476 case CEED_RESTRICTION_CURL_ORIENTED: 477 if (use_signs && use_orients) { 478 CeedCallBackend(CeedElemRestrictionApplyCurlOrientedNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 479 elem_size, v_offset, uu, vv)); 480 } else if (use_orients) { 481 CeedCallBackend(CeedElemRestrictionApplyCurlOrientedUnsignedNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, 482 num_elem, elem_size, v_offset, uu, vv)); 483 } else { 484 CeedCallBackend(CeedElemRestrictionApplyOffsetNoTranspose_Ref_Core(rstr, num_comp, block_size, comp_stride, start, stop, num_elem, 485 elem_size, v_offset, uu, vv)); 486 } 487 break; 488 case CEED_RESTRICTION_POINTS: 489 CeedCallBackend(CeedElemRestrictionApplyAtPointsInElement_Ref_Core(rstr, num_comp, start, stop, t_mode, uu, vv)); 490 break; 491 } 492 } 493 CeedCallBackend(CeedVectorRestoreArrayRead(u, &uu)); 494 CeedCallBackend(CeedVectorRestoreArray(v, &vv)); 495 if (request != CEED_REQUEST_IMMEDIATE && request != CEED_REQUEST_ORDERED) *request = NULL; 496 return CEED_ERROR_SUCCESS; 497 } 498 499 //------------------------------------------------------------------------------ 500 // ElemRestriction Apply - Common Sizes 501 //------------------------------------------------------------------------------ 502 static int CeedElemRestrictionApply_Ref_110(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 503 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 504 CeedVector v, CeedRequest *request) { 505 return CeedElemRestrictionApply_Ref_Core(rstr, 1, 1, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 506 } 507 508 static int CeedElemRestrictionApply_Ref_111(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 509 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 510 CeedVector v, CeedRequest *request) { 511 return CeedElemRestrictionApply_Ref_Core(rstr, 1, 1, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 512 } 513 514 static int CeedElemRestrictionApply_Ref_180(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 515 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 516 CeedVector v, CeedRequest *request) { 517 return CeedElemRestrictionApply_Ref_Core(rstr, 1, 8, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 518 } 519 520 static int CeedElemRestrictionApply_Ref_181(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 521 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 522 CeedVector v, CeedRequest *request) { 523 return CeedElemRestrictionApply_Ref_Core(rstr, 1, 8, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 524 } 525 526 static int CeedElemRestrictionApply_Ref_310(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 527 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 528 CeedVector v, CeedRequest *request) { 529 return CeedElemRestrictionApply_Ref_Core(rstr, 3, 1, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 530 } 531 532 static int CeedElemRestrictionApply_Ref_311(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 533 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 534 CeedVector v, CeedRequest *request) { 535 return CeedElemRestrictionApply_Ref_Core(rstr, 3, 1, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 536 } 537 538 static int CeedElemRestrictionApply_Ref_380(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 539 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 540 CeedVector v, CeedRequest *request) { 541 return CeedElemRestrictionApply_Ref_Core(rstr, 3, 8, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 542 } 543 544 static int CeedElemRestrictionApply_Ref_381(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 545 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 546 CeedVector v, CeedRequest *request) { 547 return CeedElemRestrictionApply_Ref_Core(rstr, 3, 8, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 548 } 549 550 // LCOV_EXCL_START 551 static int CeedElemRestrictionApply_Ref_510(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 552 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 553 CeedVector v, CeedRequest *request) { 554 return CeedElemRestrictionApply_Ref_Core(rstr, 5, 1, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 555 } 556 // LCOV_EXCL_STOP 557 558 static int CeedElemRestrictionApply_Ref_511(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 559 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 560 CeedVector v, CeedRequest *request) { 561 return CeedElemRestrictionApply_Ref_Core(rstr, 5, 1, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 562 } 563 564 // LCOV_EXCL_START 565 static int CeedElemRestrictionApply_Ref_580(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 566 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 567 CeedVector v, CeedRequest *request) { 568 return CeedElemRestrictionApply_Ref_Core(rstr, 5, 8, comp_stride, start, stop, t_mode, use_signs, use_orients, u, v, request); 569 } 570 // LCOV_EXCL_STOP 571 572 static int CeedElemRestrictionApply_Ref_581(CeedElemRestriction rstr, const CeedInt num_comp, const CeedInt block_size, const CeedInt comp_stride, 573 CeedInt start, CeedInt stop, CeedTransposeMode t_mode, bool use_signs, bool use_orients, CeedVector u, 574 CeedVector v, CeedRequest *request) { 575 return CeedElemRestrictionApply_Ref_Core(rstr, 5, 8, 1, start, stop, t_mode, use_signs, use_orients, u, v, request); 576 } 577 578 //------------------------------------------------------------------------------ 579 // ElemRestriction Apply 580 //------------------------------------------------------------------------------ 581 static int CeedElemRestrictionApply_Ref(CeedElemRestriction rstr, CeedTransposeMode t_mode, CeedVector u, CeedVector v, CeedRequest *request) { 582 CeedInt num_block, block_size, num_comp, comp_stride; 583 CeedElemRestriction_Ref *impl; 584 585 CeedCallBackend(CeedElemRestrictionGetNumBlocks(rstr, &num_block)); 586 CeedCallBackend(CeedElemRestrictionGetBlockSize(rstr, &block_size)); 587 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 588 CeedCallBackend(CeedElemRestrictionGetCompStride(rstr, &comp_stride)); 589 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 590 CeedCallBackend(impl->Apply(rstr, num_comp, block_size, comp_stride, 0, num_block, t_mode, true, true, u, v, request)); 591 return CEED_ERROR_SUCCESS; 592 } 593 594 //------------------------------------------------------------------------------ 595 // ElemRestriction Apply Unsigned 596 //------------------------------------------------------------------------------ 597 static int CeedElemRestrictionApplyUnsigned_Ref(CeedElemRestriction rstr, CeedTransposeMode t_mode, CeedVector u, CeedVector v, 598 CeedRequest *request) { 599 CeedInt num_block, block_size, num_comp, comp_stride; 600 CeedElemRestriction_Ref *impl; 601 602 CeedCallBackend(CeedElemRestrictionGetNumBlocks(rstr, &num_block)); 603 CeedCallBackend(CeedElemRestrictionGetBlockSize(rstr, &block_size)); 604 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 605 CeedCallBackend(CeedElemRestrictionGetCompStride(rstr, &comp_stride)); 606 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 607 CeedCallBackend(impl->Apply(rstr, num_comp, block_size, comp_stride, 0, num_block, t_mode, false, true, u, v, request)); 608 return CEED_ERROR_SUCCESS; 609 } 610 611 //------------------------------------------------------------------------------ 612 // ElemRestriction Apply Unoriented 613 //------------------------------------------------------------------------------ 614 static int CeedElemRestrictionApplyUnoriented_Ref(CeedElemRestriction rstr, CeedTransposeMode t_mode, CeedVector u, CeedVector v, 615 CeedRequest *request) { 616 CeedInt num_block, block_size, num_comp, comp_stride; 617 CeedElemRestriction_Ref *impl; 618 619 CeedCallBackend(CeedElemRestrictionGetNumBlocks(rstr, &num_block)); 620 CeedCallBackend(CeedElemRestrictionGetBlockSize(rstr, &block_size)); 621 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 622 CeedCallBackend(CeedElemRestrictionGetCompStride(rstr, &comp_stride)); 623 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 624 CeedCallBackend(impl->Apply(rstr, num_comp, block_size, comp_stride, 0, num_block, t_mode, false, false, u, v, request)); 625 return CEED_ERROR_SUCCESS; 626 } 627 628 //------------------------------------------------------------------------------ 629 // ElemRestriction Apply Points 630 //------------------------------------------------------------------------------ 631 static int CeedElemRestrictionApplyAtPointsInElement_Ref(CeedElemRestriction rstr, CeedInt elem, CeedTransposeMode t_mode, CeedVector u, CeedVector v, 632 CeedRequest *request) { 633 CeedInt num_comp; 634 CeedElemRestriction_Ref *impl; 635 636 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 637 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 638 return impl->Apply(rstr, num_comp, 0, 1, elem, elem + 1, t_mode, false, false, u, v, request); 639 } 640 641 //------------------------------------------------------------------------------ 642 // ElemRestriction Apply Block 643 //------------------------------------------------------------------------------ 644 static int CeedElemRestrictionApplyBlock_Ref(CeedElemRestriction rstr, CeedInt block, CeedTransposeMode t_mode, CeedVector u, CeedVector v, 645 CeedRequest *request) { 646 CeedInt block_size, num_comp, comp_stride; 647 CeedElemRestriction_Ref *impl; 648 649 CeedCallBackend(CeedElemRestrictionGetBlockSize(rstr, &block_size)); 650 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 651 CeedCallBackend(CeedElemRestrictionGetCompStride(rstr, &comp_stride)); 652 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 653 CeedCallBackend(impl->Apply(rstr, num_comp, block_size, comp_stride, block, block + 1, t_mode, true, true, u, v, request)); 654 return CEED_ERROR_SUCCESS; 655 } 656 657 //------------------------------------------------------------------------------ 658 // ElemRestriction Get Offsets 659 //------------------------------------------------------------------------------ 660 static int CeedElemRestrictionGetOffsets_Ref(CeedElemRestriction rstr, CeedMemType mem_type, const CeedInt **offsets) { 661 Ceed ceed; 662 CeedElemRestriction_Ref *impl; 663 664 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 665 CeedCallBackend(CeedElemRestrictionGetCeed(rstr, &ceed)); 666 667 CeedCheck(mem_type == CEED_MEM_HOST, ceed, CEED_ERROR_BACKEND, "Can only provide to HOST memory"); 668 669 *offsets = impl->offsets; 670 return CEED_ERROR_SUCCESS; 671 } 672 673 //------------------------------------------------------------------------------ 674 // ElemRestriction Get Orientations 675 //------------------------------------------------------------------------------ 676 static int CeedElemRestrictionGetOrientations_Ref(CeedElemRestriction rstr, CeedMemType mem_type, const bool **orients) { 677 Ceed ceed; 678 CeedElemRestriction_Ref *impl; 679 680 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 681 CeedCallBackend(CeedElemRestrictionGetCeed(rstr, &ceed)); 682 683 CeedCheck(mem_type == CEED_MEM_HOST, ceed, CEED_ERROR_BACKEND, "Can only provide to HOST memory"); 684 685 *orients = impl->orients; 686 return CEED_ERROR_SUCCESS; 687 } 688 689 //------------------------------------------------------------------------------ 690 // ElemRestriction Get Curl-Conforming Orientations 691 //------------------------------------------------------------------------------ 692 static int CeedElemRestrictionGetCurlOrientations_Ref(CeedElemRestriction rstr, CeedMemType mem_type, const CeedInt8 **curl_orients) { 693 Ceed ceed; 694 CeedElemRestriction_Ref *impl; 695 696 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 697 CeedCallBackend(CeedElemRestrictionGetCeed(rstr, &ceed)); 698 699 CeedCheck(mem_type == CEED_MEM_HOST, ceed, CEED_ERROR_BACKEND, "Can only provide to HOST memory"); 700 701 *curl_orients = impl->curl_orients; 702 return CEED_ERROR_SUCCESS; 703 } 704 705 //------------------------------------------------------------------------------ 706 // ElemRestriction Destroy 707 //------------------------------------------------------------------------------ 708 static int CeedElemRestrictionDestroy_Ref(CeedElemRestriction rstr) { 709 CeedElemRestriction_Ref *impl; 710 711 CeedCallBackend(CeedElemRestrictionGetData(rstr, &impl)); 712 CeedCallBackend(CeedFree(&impl->offsets_allocated)); 713 CeedCallBackend(CeedFree(&impl->orients_allocated)); 714 CeedCallBackend(CeedFree(&impl->curl_orients_allocated)); 715 CeedCallBackend(CeedFree(&impl)); 716 return CEED_ERROR_SUCCESS; 717 } 718 719 //------------------------------------------------------------------------------ 720 // ElemRestriction Create 721 //------------------------------------------------------------------------------ 722 int CeedElemRestrictionCreate_Ref(CeedMemType mem_type, CeedCopyMode copy_mode, const CeedInt *offsets, const bool *orients, 723 const CeedInt8 *curl_orients, CeedElemRestriction rstr) { 724 Ceed ceed; 725 CeedInt num_elem, elem_size, num_block, block_size, num_comp, comp_stride, num_points = 0, num_offsets; 726 CeedRestrictionType rstr_type; 727 CeedElemRestriction_Ref *impl; 728 729 CeedCallBackend(CeedElemRestrictionGetCeed(rstr, &ceed)); 730 CeedCallBackend(CeedElemRestrictionGetNumElements(rstr, &num_elem)); 731 CeedCallBackend(CeedElemRestrictionGetElementSize(rstr, &elem_size)); 732 CeedCallBackend(CeedElemRestrictionGetNumBlocks(rstr, &num_block)); 733 CeedCallBackend(CeedElemRestrictionGetBlockSize(rstr, &block_size)); 734 CeedCallBackend(CeedElemRestrictionGetNumComponents(rstr, &num_comp)); 735 CeedCallBackend(CeedElemRestrictionGetCompStride(rstr, &comp_stride)); 736 CeedInt layout[3] = {1, elem_size, elem_size * num_comp}; 737 738 CeedCheck(mem_type == CEED_MEM_HOST, ceed, CEED_ERROR_BACKEND, "Only MemType = HOST supported"); 739 740 CeedCallBackend(CeedCalloc(1, &impl)); 741 CeedCallBackend(CeedElemRestrictionSetData(rstr, impl)); 742 CeedCallBackend(CeedElemRestrictionSetELayout(rstr, layout)); 743 744 // Offsets data 745 CeedCallBackend(CeedElemRestrictionGetType(rstr, &rstr_type)); 746 if (rstr_type != CEED_RESTRICTION_STRIDED) { 747 const char *resource; 748 749 // Check indices for ref or memcheck backends 750 { 751 Ceed current = ceed, parent = NULL; 752 753 CeedCallBackend(CeedGetParent(current, &parent)); 754 while (current != parent) { 755 current = parent; 756 CeedCallBackend(CeedGetParent(current, &parent)); 757 } 758 CeedCallBackend(CeedGetResource(parent, &resource)); 759 } 760 if (!strcmp(resource, "/cpu/self/ref/serial") || !strcmp(resource, "/cpu/self/ref/blocked") || !strcmp(resource, "/cpu/self/memcheck/serial") || 761 !strcmp(resource, "/cpu/self/memcheck/blocked")) { 762 CeedSize l_size; 763 764 CeedCallBackend(CeedElemRestrictionGetLVectorSize(rstr, &l_size)); 765 for (CeedInt i = 0; i < num_elem * elem_size; i++) { 766 CeedCheck(offsets[i] >= 0 && offsets[i] + (num_comp - 1) * comp_stride < l_size, ceed, CEED_ERROR_BACKEND, 767 "Restriction offset %" CeedInt_FMT " (%" CeedInt_FMT ") out of range [0, %" CeedInt_FMT "]", i, offsets[i], l_size); 768 } 769 } 770 771 // Copy data 772 if (rstr_type == CEED_RESTRICTION_POINTS) CeedCallBackend(CeedElemRestrictionGetNumPoints(rstr, &num_points)); 773 num_offsets = rstr_type == CEED_RESTRICTION_POINTS ? (num_elem + 1 + num_points) : (num_elem * elem_size); 774 switch (copy_mode) { 775 case CEED_COPY_VALUES: 776 CeedCallBackend(CeedMalloc(num_offsets, &impl->offsets_allocated)); 777 memcpy(impl->offsets_allocated, offsets, num_offsets * sizeof(offsets[0])); 778 impl->offsets = impl->offsets_allocated; 779 break; 780 case CEED_OWN_POINTER: 781 impl->offsets_allocated = (CeedInt *)offsets; 782 impl->offsets = impl->offsets_allocated; 783 break; 784 case CEED_USE_POINTER: 785 impl->offsets = offsets; 786 } 787 788 // Orientation data 789 if (rstr_type == CEED_RESTRICTION_ORIENTED) { 790 CeedCheck(orients != NULL, ceed, CEED_ERROR_BACKEND, "No orients array provided for oriented restriction"); 791 switch (copy_mode) { 792 case CEED_COPY_VALUES: 793 CeedCallBackend(CeedMalloc(num_offsets, &impl->orients_allocated)); 794 memcpy(impl->orients_allocated, orients, num_offsets * sizeof(orients[0])); 795 impl->orients = impl->orients_allocated; 796 break; 797 case CEED_OWN_POINTER: 798 impl->orients_allocated = (bool *)orients; 799 impl->orients = impl->orients_allocated; 800 break; 801 case CEED_USE_POINTER: 802 impl->orients = orients; 803 } 804 } else if (rstr_type == CEED_RESTRICTION_CURL_ORIENTED) { 805 CeedCheck(curl_orients != NULL, ceed, CEED_ERROR_BACKEND, "No curl_orients array provided for oriented restriction"); 806 switch (copy_mode) { 807 case CEED_COPY_VALUES: 808 CeedCallBackend(CeedMalloc(3 * num_offsets, &impl->curl_orients_allocated)); 809 memcpy(impl->curl_orients_allocated, curl_orients, 3 * num_offsets * sizeof(curl_orients[0])); 810 impl->curl_orients = impl->curl_orients_allocated; 811 break; 812 case CEED_OWN_POINTER: 813 impl->curl_orients_allocated = (CeedInt8 *)curl_orients; 814 impl->curl_orients = impl->curl_orients_allocated; 815 break; 816 case CEED_USE_POINTER: 817 impl->curl_orients = curl_orients; 818 } 819 } 820 } 821 822 // Set apply function based upon num_comp, block_size, and comp_stride 823 CeedInt index = -1; 824 825 if (block_size < 10) index = 100 * num_comp + 10 * block_size + (comp_stride == 1); 826 switch (index) { 827 case 110: 828 impl->Apply = CeedElemRestrictionApply_Ref_110; 829 break; 830 case 111: 831 impl->Apply = CeedElemRestrictionApply_Ref_111; 832 break; 833 case 180: 834 impl->Apply = CeedElemRestrictionApply_Ref_180; 835 break; 836 case 181: 837 impl->Apply = CeedElemRestrictionApply_Ref_181; 838 break; 839 case 310: 840 impl->Apply = CeedElemRestrictionApply_Ref_310; 841 break; 842 case 311: 843 impl->Apply = CeedElemRestrictionApply_Ref_311; 844 break; 845 case 380: 846 impl->Apply = CeedElemRestrictionApply_Ref_380; 847 break; 848 case 381: 849 impl->Apply = CeedElemRestrictionApply_Ref_381; 850 break; 851 // LCOV_EXCL_START 852 case 510: 853 impl->Apply = CeedElemRestrictionApply_Ref_510; 854 break; 855 // LCOV_EXCL_STOP 856 case 511: 857 impl->Apply = CeedElemRestrictionApply_Ref_511; 858 break; 859 // LCOV_EXCL_START 860 case 580: 861 impl->Apply = CeedElemRestrictionApply_Ref_580; 862 break; 863 // LCOV_EXCL_STOP 864 case 581: 865 impl->Apply = CeedElemRestrictionApply_Ref_581; 866 break; 867 default: 868 impl->Apply = CeedElemRestrictionApply_Ref_Core; 869 break; 870 } 871 872 // Register backend functions 873 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "Apply", CeedElemRestrictionApply_Ref)); 874 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "ApplyUnsigned", CeedElemRestrictionApplyUnsigned_Ref)); 875 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "ApplyUnoriented", CeedElemRestrictionApplyUnoriented_Ref)); 876 if (rstr_type == CEED_RESTRICTION_POINTS) { 877 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "ApplyAtPointsInElement", CeedElemRestrictionApplyAtPointsInElement_Ref)); 878 } 879 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "ApplyBlock", CeedElemRestrictionApplyBlock_Ref)); 880 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "GetOffsets", CeedElemRestrictionGetOffsets_Ref)); 881 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "GetOrientations", CeedElemRestrictionGetOrientations_Ref)); 882 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "GetCurlOrientations", CeedElemRestrictionGetCurlOrientations_Ref)); 883 CeedCallBackend(CeedSetBackendFunction(ceed, "ElemRestriction", rstr, "Destroy", CeedElemRestrictionDestroy_Ref)); 884 return CEED_ERROR_SUCCESS; 885 } 886 887 //------------------------------------------------------------------------------ 888