1d7b241e6Sjeremylt // Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2d7b241e6Sjeremylt // the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3d7b241e6Sjeremylt // reserved. See files LICENSE and NOTICE for details. 4d7b241e6Sjeremylt // 5d7b241e6Sjeremylt // This file is part of CEED, a collection of benchmarks, miniapps, software 6d7b241e6Sjeremylt // libraries and APIs for efficient high-order finite element and spectral 7d7b241e6Sjeremylt // element discretizations for exascale applications. For more information and 8d7b241e6Sjeremylt // source code availability see http://github.com/ceed. 9d7b241e6Sjeremylt // 10d7b241e6Sjeremylt // The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11d7b241e6Sjeremylt // a collaborative effort of two U.S. Department of Energy organizations (Office 12d7b241e6Sjeremylt // of Science and the National Nuclear Security Administration) responsible for 13d7b241e6Sjeremylt // the planning and preparation of a capable exascale ecosystem, including 14d7b241e6Sjeremylt // software, applications, hardware, advanced system engineering and early 15d7b241e6Sjeremylt // testbed platforms, in support of the nation's exascale computing imperative. 16d7b241e6Sjeremylt 17d7b241e6Sjeremylt #include <ceed-impl.h> 18d863ab9bSjeremylt #include <ceed-backend.h> 19d7b241e6Sjeremylt #include <string.h> 203bd813ffSjeremylt #include <math.h> 21d7b241e6Sjeremylt 22dfdf5a53Sjeremylt /// @file 237a982d89SJeremy L. Thompson /// Implementation of CeedOperator interfaces 247a982d89SJeremy L. Thompson 257a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 267a982d89SJeremy L. Thompson /// CeedOperator Library Internal Functions 277a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 287a982d89SJeremy L. Thompson /// @addtogroup CeedOperatorDeveloper 297a982d89SJeremy L. Thompson /// @{ 307a982d89SJeremy L. Thompson 317a982d89SJeremy L. Thompson /** 327a982d89SJeremy L. Thompson @brief Duplicate a CeedOperator with a reference Ceed to fallback for advanced 337a982d89SJeremy L. Thompson CeedOperator functionality 347a982d89SJeremy L. Thompson 357a982d89SJeremy L. Thompson @param op CeedOperator to create fallback for 367a982d89SJeremy L. Thompson 377a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 387a982d89SJeremy L. Thompson 397a982d89SJeremy L. Thompson @ref Developer 407a982d89SJeremy L. Thompson **/ 417a982d89SJeremy L. Thompson int CeedOperatorCreateFallback(CeedOperator op) { 427a982d89SJeremy L. Thompson int ierr; 437a982d89SJeremy L. Thompson 447a982d89SJeremy L. Thompson // Fallback Ceed 457a982d89SJeremy L. Thompson const char *resource, *fallbackresource; 467a982d89SJeremy L. Thompson ierr = CeedGetResource(op->ceed, &resource); CeedChk(ierr); 477a982d89SJeremy L. Thompson ierr = CeedGetOperatorFallbackResource(op->ceed, &fallbackresource); 487a982d89SJeremy L. Thompson CeedChk(ierr); 497a982d89SJeremy L. Thompson if (!strcmp(resource, fallbackresource)) 507a982d89SJeremy L. Thompson // LCOV_EXCL_START 517a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Backend %s cannot create an operator" 527a982d89SJeremy L. Thompson "fallback to resource %s", resource, fallbackresource); 537a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 547a982d89SJeremy L. Thompson 557a982d89SJeremy L. Thompson // Fallback Ceed 567a982d89SJeremy L. Thompson Ceed ceedref; 577a982d89SJeremy L. Thompson if (!op->ceed->opfallbackceed) { 587a982d89SJeremy L. Thompson ierr = CeedInit(fallbackresource, &ceedref); CeedChk(ierr); 597a982d89SJeremy L. Thompson ceedref->opfallbackparent = op->ceed; 607a982d89SJeremy L. Thompson op->ceed->opfallbackceed = ceedref; 617a982d89SJeremy L. Thompson } 627a982d89SJeremy L. Thompson ceedref = op->ceed->opfallbackceed; 637a982d89SJeremy L. Thompson 647a982d89SJeremy L. Thompson // Clone Op 657a982d89SJeremy L. Thompson CeedOperator opref; 667a982d89SJeremy L. Thompson ierr = CeedCalloc(1, &opref); CeedChk(ierr); 677a982d89SJeremy L. Thompson memcpy(opref, op, sizeof(*opref)); CeedChk(ierr); 687a982d89SJeremy L. Thompson opref->data = NULL; 697a982d89SJeremy L. Thompson opref->setupdone = 0; 707a982d89SJeremy L. Thompson opref->ceed = ceedref; 717a982d89SJeremy L. Thompson ierr = ceedref->OperatorCreate(opref); CeedChk(ierr); 727a982d89SJeremy L. Thompson op->opfallback = opref; 737a982d89SJeremy L. Thompson 747a982d89SJeremy L. Thompson // Clone QF 757a982d89SJeremy L. Thompson CeedQFunction qfref; 767a982d89SJeremy L. Thompson ierr = CeedCalloc(1, &qfref); CeedChk(ierr); 777a982d89SJeremy L. Thompson memcpy(qfref, (op->qf), sizeof(*qfref)); CeedChk(ierr); 787a982d89SJeremy L. Thompson qfref->data = NULL; 797a982d89SJeremy L. Thompson qfref->ceed = ceedref; 807a982d89SJeremy L. Thompson ierr = ceedref->QFunctionCreate(qfref); CeedChk(ierr); 817a982d89SJeremy L. Thompson opref->qf = qfref; 827a982d89SJeremy L. Thompson op->qffallback = qfref; 837a982d89SJeremy L. Thompson 847a982d89SJeremy L. Thompson return 0; 857a982d89SJeremy L. Thompson } 867a982d89SJeremy L. Thompson 877a982d89SJeremy L. Thompson /** 887a982d89SJeremy L. Thompson @brief Check if a CeedOperator is ready to be used. 897a982d89SJeremy L. Thompson 907a982d89SJeremy L. Thompson @param[in] ceed Ceed object for error handling 917a982d89SJeremy L. Thompson @param[in] op CeedOperator to check 927a982d89SJeremy L. Thompson 937a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 947a982d89SJeremy L. Thompson 957a982d89SJeremy L. Thompson @ref Developer 967a982d89SJeremy L. Thompson **/ 977a982d89SJeremy L. Thompson static int CeedOperatorCheckReady(Ceed ceed, CeedOperator op) { 987a982d89SJeremy L. Thompson CeedQFunction qf = op->qf; 997a982d89SJeremy L. Thompson 1007a982d89SJeremy L. Thompson if (op->composite) { 1017a982d89SJeremy L. Thompson if (!op->numsub) 1027a982d89SJeremy L. Thompson // LCOV_EXCL_START 1037a982d89SJeremy L. Thompson return CeedError(ceed, 1, "No suboperators set"); 1047a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 1057a982d89SJeremy L. Thompson } else { 1067a982d89SJeremy L. Thompson if (op->nfields == 0) 1077a982d89SJeremy L. Thompson // LCOV_EXCL_START 1087a982d89SJeremy L. Thompson return CeedError(ceed, 1, "No operator fields set"); 1097a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 1107a982d89SJeremy L. Thompson if (op->nfields < qf->numinputfields + qf->numoutputfields) 1117a982d89SJeremy L. Thompson // LCOV_EXCL_START 1127a982d89SJeremy L. Thompson return CeedError(ceed, 1, "Not all operator fields set"); 1137a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 1147a982d89SJeremy L. Thompson if (!op->hasrestriction) 1157a982d89SJeremy L. Thompson // LCOV_EXCL_START 1167a982d89SJeremy L. Thompson return CeedError(ceed, 1,"At least one restriction required"); 1177a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 1187a982d89SJeremy L. Thompson if (op->numqpoints == 0) 1197a982d89SJeremy L. Thompson // LCOV_EXCL_START 1207a982d89SJeremy L. Thompson return CeedError(ceed, 1,"At least one non-collocated basis required"); 1217a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 1227a982d89SJeremy L. Thompson } 1237a982d89SJeremy L. Thompson 1247a982d89SJeremy L. Thompson return 0; 1257a982d89SJeremy L. Thompson } 1267a982d89SJeremy L. Thompson 1277a982d89SJeremy L. Thompson /** 1287a982d89SJeremy L. Thompson @brief View a field of a CeedOperator 1297a982d89SJeremy L. Thompson 1307a982d89SJeremy L. Thompson @param[in] field Operator field to view 1314c4400c7SValeria Barra @param[in] qffield QFunction field (carries field name) 1327a982d89SJeremy L. Thompson @param[in] fieldnumber Number of field being viewed 1334c4400c7SValeria Barra @param[in] sub true indicates sub-operator, which increases indentation; false for top-level operator 1344c4400c7SValeria Barra @param[in] in true for an input field; false for output field 1357a982d89SJeremy L. Thompson @param[in] stream Stream to view to, e.g., stdout 1367a982d89SJeremy L. Thompson 1377a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 1387a982d89SJeremy L. Thompson 1397a982d89SJeremy L. Thompson @ref Utility 1407a982d89SJeremy L. Thompson **/ 1417a982d89SJeremy L. Thompson static int CeedOperatorFieldView(CeedOperatorField field, 1427a982d89SJeremy L. Thompson CeedQFunctionField qffield, 1437a982d89SJeremy L. Thompson CeedInt fieldnumber, bool sub, bool in, 1447a982d89SJeremy L. Thompson FILE *stream) { 1457a982d89SJeremy L. Thompson const char *pre = sub ? " " : ""; 1467a982d89SJeremy L. Thompson const char *inout = in ? "Input" : "Output"; 1477a982d89SJeremy L. Thompson 1487a982d89SJeremy L. Thompson fprintf(stream, "%s %s Field [%d]:\n" 1497a982d89SJeremy L. Thompson "%s Name: \"%s\"\n", 1507a982d89SJeremy L. Thompson pre, inout, fieldnumber, pre, qffield->fieldname); 1517a982d89SJeremy L. Thompson 1527a982d89SJeremy L. Thompson if (field->basis == CEED_BASIS_COLLOCATED) 1537a982d89SJeremy L. Thompson fprintf(stream, "%s Collocated basis\n", pre); 1547a982d89SJeremy L. Thompson 1557a982d89SJeremy L. Thompson if (field->vec == CEED_VECTOR_ACTIVE) 1567a982d89SJeremy L. Thompson fprintf(stream, "%s Active vector\n", pre); 1577a982d89SJeremy L. Thompson else if (field->vec == CEED_VECTOR_NONE) 1587a982d89SJeremy L. Thompson fprintf(stream, "%s No vector\n", pre); 1597a982d89SJeremy L. Thompson 1607a982d89SJeremy L. Thompson return 0; 1617a982d89SJeremy L. Thompson } 1627a982d89SJeremy L. Thompson 1637a982d89SJeremy L. Thompson /** 1647a982d89SJeremy L. Thompson @brief View a single CeedOperator 1657a982d89SJeremy L. Thompson 1667a982d89SJeremy L. Thompson @param[in] op CeedOperator to view 1677a982d89SJeremy L. Thompson @param[in] sub Boolean flag for sub-operator 1687a982d89SJeremy L. Thompson @param[in] stream Stream to write; typically stdout/stderr or a file 1697a982d89SJeremy L. Thompson 1707a982d89SJeremy L. Thompson @return Error code: 0 - success, otherwise - failure 1717a982d89SJeremy L. Thompson 1727a982d89SJeremy L. Thompson @ref Utility 1737a982d89SJeremy L. Thompson **/ 1747a982d89SJeremy L. Thompson int CeedOperatorSingleView(CeedOperator op, bool sub, FILE *stream) { 1757a982d89SJeremy L. Thompson int ierr; 1767a982d89SJeremy L. Thompson const char *pre = sub ? " " : ""; 1777a982d89SJeremy L. Thompson 1787a982d89SJeremy L. Thompson CeedInt totalfields; 1797a982d89SJeremy L. Thompson ierr = CeedOperatorGetNumArgs(op, &totalfields); CeedChk(ierr); 1807a982d89SJeremy L. Thompson 1817a982d89SJeremy L. Thompson fprintf(stream, "%s %d Field%s\n", pre, totalfields, totalfields>1 ? "s" : ""); 1827a982d89SJeremy L. Thompson 1837a982d89SJeremy L. Thompson fprintf(stream, "%s %d Input Field%s:\n", pre, op->qf->numinputfields, 1847a982d89SJeremy L. Thompson op->qf->numinputfields>1 ? "s" : ""); 1857a982d89SJeremy L. Thompson for (CeedInt i=0; i<op->qf->numinputfields; i++) { 1867a982d89SJeremy L. Thompson ierr = CeedOperatorFieldView(op->inputfields[i], op->qf->inputfields[i], 1877a982d89SJeremy L. Thompson i, sub, 1, stream); CeedChk(ierr); 1887a982d89SJeremy L. Thompson } 1897a982d89SJeremy L. Thompson 1907a982d89SJeremy L. Thompson fprintf(stream, "%s %d Output Field%s:\n", pre, op->qf->numoutputfields, 1917a982d89SJeremy L. Thompson op->qf->numoutputfields>1 ? "s" : ""); 1927a982d89SJeremy L. Thompson for (CeedInt i=0; i<op->qf->numoutputfields; i++) { 1937a982d89SJeremy L. Thompson ierr = CeedOperatorFieldView(op->outputfields[i], op->qf->inputfields[i], 1947a982d89SJeremy L. Thompson i, sub, 0, stream); CeedChk(ierr); 1957a982d89SJeremy L. Thompson } 1967a982d89SJeremy L. Thompson 1977a982d89SJeremy L. Thompson return 0; 1987a982d89SJeremy L. Thompson } 1997a982d89SJeremy L. Thompson 2007a982d89SJeremy L. Thompson /// @} 2017a982d89SJeremy L. Thompson 2027a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 2037a982d89SJeremy L. Thompson /// CeedOperator Backend API 2047a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 2057a982d89SJeremy L. Thompson /// @addtogroup CeedOperatorBackend 2067a982d89SJeremy L. Thompson /// @{ 2077a982d89SJeremy L. Thompson 2087a982d89SJeremy L. Thompson /** 2097a982d89SJeremy L. Thompson @brief Get the Ceed associated with a CeedOperator 2107a982d89SJeremy L. Thompson 2117a982d89SJeremy L. Thompson @param op CeedOperator 2127a982d89SJeremy L. Thompson @param[out] ceed Variable to store Ceed 2137a982d89SJeremy L. Thompson 2147a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 2157a982d89SJeremy L. Thompson 2167a982d89SJeremy L. Thompson @ref Backend 2177a982d89SJeremy L. Thompson **/ 2187a982d89SJeremy L. Thompson 2197a982d89SJeremy L. Thompson int CeedOperatorGetCeed(CeedOperator op, Ceed *ceed) { 2207a982d89SJeremy L. Thompson *ceed = op->ceed; 2217a982d89SJeremy L. Thompson return 0; 2227a982d89SJeremy L. Thompson } 2237a982d89SJeremy L. Thompson 2247a982d89SJeremy L. Thompson /** 2257a982d89SJeremy L. Thompson @brief Get the number of elements associated with a CeedOperator 2267a982d89SJeremy L. Thompson 2277a982d89SJeremy L. Thompson @param op CeedOperator 2287a982d89SJeremy L. Thompson @param[out] numelem Variable to store number of elements 2297a982d89SJeremy L. Thompson 2307a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 2317a982d89SJeremy L. Thompson 2327a982d89SJeremy L. Thompson @ref Backend 2337a982d89SJeremy L. Thompson **/ 2347a982d89SJeremy L. Thompson 2357a982d89SJeremy L. Thompson int CeedOperatorGetNumElements(CeedOperator op, CeedInt *numelem) { 2367a982d89SJeremy L. Thompson if (op->composite) 2377a982d89SJeremy L. Thompson // LCOV_EXCL_START 2387a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not defined for composite operator"); 2397a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 2407a982d89SJeremy L. Thompson 2417a982d89SJeremy L. Thompson *numelem = op->numelements; 2427a982d89SJeremy L. Thompson return 0; 2437a982d89SJeremy L. Thompson } 2447a982d89SJeremy L. Thompson 2457a982d89SJeremy L. Thompson /** 2467a982d89SJeremy L. Thompson @brief Get the number of quadrature points associated with a CeedOperator 2477a982d89SJeremy L. Thompson 2487a982d89SJeremy L. Thompson @param op CeedOperator 2497a982d89SJeremy L. Thompson @param[out] numqpts Variable to store vector number of quadrature points 2507a982d89SJeremy L. Thompson 2517a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 2527a982d89SJeremy L. Thompson 2537a982d89SJeremy L. Thompson @ref Backend 2547a982d89SJeremy L. Thompson **/ 2557a982d89SJeremy L. Thompson 2567a982d89SJeremy L. Thompson int CeedOperatorGetNumQuadraturePoints(CeedOperator op, CeedInt *numqpts) { 2577a982d89SJeremy L. Thompson if (op->composite) 2587a982d89SJeremy L. Thompson // LCOV_EXCL_START 2597a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not defined for composite operator"); 2607a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 2617a982d89SJeremy L. Thompson 2627a982d89SJeremy L. Thompson *numqpts = op->numqpoints; 2637a982d89SJeremy L. Thompson return 0; 2647a982d89SJeremy L. Thompson } 2657a982d89SJeremy L. Thompson 2667a982d89SJeremy L. Thompson /** 2677a982d89SJeremy L. Thompson @brief Get the number of arguments associated with a CeedOperator 2687a982d89SJeremy L. Thompson 2697a982d89SJeremy L. Thompson @param op CeedOperator 2707a982d89SJeremy L. Thompson @param[out] numargs Variable to store vector number of arguments 2717a982d89SJeremy L. Thompson 2727a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 2737a982d89SJeremy L. Thompson 2747a982d89SJeremy L. Thompson @ref Backend 2757a982d89SJeremy L. Thompson **/ 2767a982d89SJeremy L. Thompson 2777a982d89SJeremy L. Thompson int CeedOperatorGetNumArgs(CeedOperator op, CeedInt *numargs) { 2787a982d89SJeremy L. Thompson if (op->composite) 2797a982d89SJeremy L. Thompson // LCOV_EXCL_START 2807a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not defined for composite operators"); 2817a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 2827a982d89SJeremy L. Thompson 2837a982d89SJeremy L. Thompson *numargs = op->nfields; 2847a982d89SJeremy L. Thompson return 0; 2857a982d89SJeremy L. Thompson } 2867a982d89SJeremy L. Thompson 2877a982d89SJeremy L. Thompson /** 2887a982d89SJeremy L. Thompson @brief Get the setup status of a CeedOperator 2897a982d89SJeremy L. Thompson 2907a982d89SJeremy L. Thompson @param op CeedOperator 2917a982d89SJeremy L. Thompson @param[out] setupdone Variable to store setup status 2927a982d89SJeremy L. Thompson 2937a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 2947a982d89SJeremy L. Thompson 2957a982d89SJeremy L. Thompson @ref Backend 2967a982d89SJeremy L. Thompson **/ 2977a982d89SJeremy L. Thompson 2987a982d89SJeremy L. Thompson int CeedOperatorGetSetupStatus(CeedOperator op, bool *setupdone) { 2997a982d89SJeremy L. Thompson *setupdone = op->setupdone; 3007a982d89SJeremy L. Thompson return 0; 3017a982d89SJeremy L. Thompson } 3027a982d89SJeremy L. Thompson 3037a982d89SJeremy L. Thompson /** 3047a982d89SJeremy L. Thompson @brief Get the QFunction associated with a CeedOperator 3057a982d89SJeremy L. Thompson 3067a982d89SJeremy L. Thompson @param op CeedOperator 3077a982d89SJeremy L. Thompson @param[out] qf Variable to store QFunction 3087a982d89SJeremy L. Thompson 3097a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 3107a982d89SJeremy L. Thompson 3117a982d89SJeremy L. Thompson @ref Backend 3127a982d89SJeremy L. Thompson **/ 3137a982d89SJeremy L. Thompson 3147a982d89SJeremy L. Thompson int CeedOperatorGetQFunction(CeedOperator op, CeedQFunction *qf) { 3157a982d89SJeremy L. Thompson if (op->composite) 3167a982d89SJeremy L. Thompson // LCOV_EXCL_START 3177a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not defined for composite operator"); 3187a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 3197a982d89SJeremy L. Thompson 3207a982d89SJeremy L. Thompson *qf = op->qf; 3217a982d89SJeremy L. Thompson return 0; 3227a982d89SJeremy L. Thompson } 3237a982d89SJeremy L. Thompson 3247a982d89SJeremy L. Thompson /** 325*c04a41a7SJeremy L Thompson @brief Get a boolean value indicating if the CeedOperator is composite 326*c04a41a7SJeremy L Thompson 327*c04a41a7SJeremy L Thompson @param op CeedOperator 328*c04a41a7SJeremy L Thompson @param[out] isComposite Variable to store composite status 329*c04a41a7SJeremy L Thompson 330*c04a41a7SJeremy L Thompson @return An error code: 0 - success, otherwise - failure 331*c04a41a7SJeremy L Thompson 332*c04a41a7SJeremy L Thompson @ref Backend 333*c04a41a7SJeremy L Thompson **/ 334*c04a41a7SJeremy L Thompson 335*c04a41a7SJeremy L Thompson int CeedOperatorGetCompositeStatus(CeedOperator op, bool *isComposite) { 336*c04a41a7SJeremy L Thompson *isComposite = op->composite; 337*c04a41a7SJeremy L Thompson return 0; 338*c04a41a7SJeremy L Thompson } 339*c04a41a7SJeremy L Thompson 340*c04a41a7SJeremy L Thompson /** 3417a982d89SJeremy L. Thompson @brief Get the number of suboperators associated with a CeedOperator 3427a982d89SJeremy L. Thompson 3437a982d89SJeremy L. Thompson @param op CeedOperator 3447a982d89SJeremy L. Thompson @param[out] numsub Variable to store number of suboperators 3457a982d89SJeremy L. Thompson 3467a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 3477a982d89SJeremy L. Thompson 3487a982d89SJeremy L. Thompson @ref Backend 3497a982d89SJeremy L. Thompson **/ 3507a982d89SJeremy L. Thompson 3517a982d89SJeremy L. Thompson int CeedOperatorGetNumSub(CeedOperator op, CeedInt *numsub) { 3527a982d89SJeremy L. Thompson if (!op->composite) 3537a982d89SJeremy L. Thompson // LCOV_EXCL_START 3547a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not a composite operator"); 3557a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 3567a982d89SJeremy L. Thompson 3577a982d89SJeremy L. Thompson *numsub = op->numsub; 3587a982d89SJeremy L. Thompson return 0; 3597a982d89SJeremy L. Thompson } 3607a982d89SJeremy L. Thompson 3617a982d89SJeremy L. Thompson /** 3627a982d89SJeremy L. Thompson @brief Get the list of suboperators associated with a CeedOperator 3637a982d89SJeremy L. Thompson 3647a982d89SJeremy L. Thompson @param op CeedOperator 3657a982d89SJeremy L. Thompson @param[out] suboperators Variable to store list of suboperators 3667a982d89SJeremy L. Thompson 3677a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 3687a982d89SJeremy L. Thompson 3697a982d89SJeremy L. Thompson @ref Backend 3707a982d89SJeremy L. Thompson **/ 3717a982d89SJeremy L. Thompson 3727a982d89SJeremy L. Thompson int CeedOperatorGetSubList(CeedOperator op, CeedOperator **suboperators) { 3737a982d89SJeremy L. Thompson if (!op->composite) 3747a982d89SJeremy L. Thompson // LCOV_EXCL_START 3757a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not a composite operator"); 3767a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 3777a982d89SJeremy L. Thompson 3787a982d89SJeremy L. Thompson *suboperators = op->suboperators; 3797a982d89SJeremy L. Thompson return 0; 3807a982d89SJeremy L. Thompson } 3817a982d89SJeremy L. Thompson 3827a982d89SJeremy L. Thompson /** 3837a982d89SJeremy L. Thompson @brief Get the backend data of a CeedOperator 3847a982d89SJeremy L. Thompson 3857a982d89SJeremy L. Thompson @param op CeedOperator 3867a982d89SJeremy L. Thompson @param[out] data Variable to store data 3877a982d89SJeremy L. Thompson 3887a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 3897a982d89SJeremy L. Thompson 3907a982d89SJeremy L. Thompson @ref Backend 3917a982d89SJeremy L. Thompson **/ 3927a982d89SJeremy L. Thompson 3937a982d89SJeremy L. Thompson int CeedOperatorGetData(CeedOperator op, void **data) { 3947a982d89SJeremy L. Thompson *data = op->data; 3957a982d89SJeremy L. Thompson return 0; 3967a982d89SJeremy L. Thompson } 3977a982d89SJeremy L. Thompson 3987a982d89SJeremy L. Thompson /** 3997a982d89SJeremy L. Thompson @brief Set the backend data of a CeedOperator 4007a982d89SJeremy L. Thompson 4017a982d89SJeremy L. Thompson @param[out] op CeedOperator 4027a982d89SJeremy L. Thompson @param data Data to set 4037a982d89SJeremy L. Thompson 4047a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4057a982d89SJeremy L. Thompson 4067a982d89SJeremy L. Thompson @ref Backend 4077a982d89SJeremy L. Thompson **/ 4087a982d89SJeremy L. Thompson 4097a982d89SJeremy L. Thompson int CeedOperatorSetData(CeedOperator op, void **data) { 4107a982d89SJeremy L. Thompson op->data = *data; 4117a982d89SJeremy L. Thompson return 0; 4127a982d89SJeremy L. Thompson } 4137a982d89SJeremy L. Thompson 4147a982d89SJeremy L. Thompson /** 4157a982d89SJeremy L. Thompson @brief Set the setup flag of a CeedOperator to True 4167a982d89SJeremy L. Thompson 4177a982d89SJeremy L. Thompson @param op CeedOperator 4187a982d89SJeremy L. Thompson 4197a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4207a982d89SJeremy L. Thompson 4217a982d89SJeremy L. Thompson @ref Backend 4227a982d89SJeremy L. Thompson **/ 4237a982d89SJeremy L. Thompson 4247a982d89SJeremy L. Thompson int CeedOperatorSetSetupDone(CeedOperator op) { 4257a982d89SJeremy L. Thompson op->setupdone = 1; 4267a982d89SJeremy L. Thompson return 0; 4277a982d89SJeremy L. Thompson } 4287a982d89SJeremy L. Thompson 4297a982d89SJeremy L. Thompson /** 4307a982d89SJeremy L. Thompson @brief Get the CeedOperatorFields of a CeedOperator 4317a982d89SJeremy L. Thompson 4327a982d89SJeremy L. Thompson @param op CeedOperator 4337a982d89SJeremy L. Thompson @param[out] inputfields Variable to store inputfields 4347a982d89SJeremy L. Thompson @param[out] outputfields Variable to store outputfields 4357a982d89SJeremy L. Thompson 4367a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4377a982d89SJeremy L. Thompson 4387a982d89SJeremy L. Thompson @ref Backend 4397a982d89SJeremy L. Thompson **/ 4407a982d89SJeremy L. Thompson 4417a982d89SJeremy L. Thompson int CeedOperatorGetFields(CeedOperator op, CeedOperatorField **inputfields, 4427a982d89SJeremy L. Thompson CeedOperatorField **outputfields) { 4437a982d89SJeremy L. Thompson if (op->composite) 4447a982d89SJeremy L. Thompson // LCOV_EXCL_START 4457a982d89SJeremy L. Thompson return CeedError(op->ceed, 1, "Not defined for composite operator"); 4467a982d89SJeremy L. Thompson // LCOV_EXCL_STOP 4477a982d89SJeremy L. Thompson 4487a982d89SJeremy L. Thompson if (inputfields) *inputfields = op->inputfields; 4497a982d89SJeremy L. Thompson if (outputfields) *outputfields = op->outputfields; 4507a982d89SJeremy L. Thompson return 0; 4517a982d89SJeremy L. Thompson } 4527a982d89SJeremy L. Thompson 4537a982d89SJeremy L. Thompson /** 4547a982d89SJeremy L. Thompson @brief Get the CeedElemRestriction of a CeedOperatorField 4557a982d89SJeremy L. Thompson 4567a982d89SJeremy L. Thompson @param opfield CeedOperatorField 4577a982d89SJeremy L. Thompson @param[out] rstr Variable to store CeedElemRestriction 4587a982d89SJeremy L. Thompson 4597a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4607a982d89SJeremy L. Thompson 4617a982d89SJeremy L. Thompson @ref Backend 4627a982d89SJeremy L. Thompson **/ 4637a982d89SJeremy L. Thompson 4647a982d89SJeremy L. Thompson int CeedOperatorFieldGetElemRestriction(CeedOperatorField opfield, 4657a982d89SJeremy L. Thompson CeedElemRestriction *rstr) { 4667a982d89SJeremy L. Thompson *rstr = opfield->Erestrict; 4677a982d89SJeremy L. Thompson return 0; 4687a982d89SJeremy L. Thompson } 4697a982d89SJeremy L. Thompson 4707a982d89SJeremy L. Thompson /** 4717a982d89SJeremy L. Thompson @brief Get the CeedBasis of a CeedOperatorField 4727a982d89SJeremy L. Thompson 4737a982d89SJeremy L. Thompson @param opfield CeedOperatorField 4747a982d89SJeremy L. Thompson @param[out] basis Variable to store CeedBasis 4757a982d89SJeremy L. Thompson 4767a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4777a982d89SJeremy L. Thompson 4787a982d89SJeremy L. Thompson @ref Backend 4797a982d89SJeremy L. Thompson **/ 4807a982d89SJeremy L. Thompson 4817a982d89SJeremy L. Thompson int CeedOperatorFieldGetBasis(CeedOperatorField opfield, CeedBasis *basis) { 4827a982d89SJeremy L. Thompson *basis = opfield->basis; 4837a982d89SJeremy L. Thompson return 0; 4847a982d89SJeremy L. Thompson } 4857a982d89SJeremy L. Thompson 4867a982d89SJeremy L. Thompson /** 4877a982d89SJeremy L. Thompson @brief Get the CeedVector of a CeedOperatorField 4887a982d89SJeremy L. Thompson 4897a982d89SJeremy L. Thompson @param opfield CeedOperatorField 4907a982d89SJeremy L. Thompson @param[out] vec Variable to store CeedVector 4917a982d89SJeremy L. Thompson 4927a982d89SJeremy L. Thompson @return An error code: 0 - success, otherwise - failure 4937a982d89SJeremy L. Thompson 4947a982d89SJeremy L. Thompson @ref Backend 4957a982d89SJeremy L. Thompson **/ 4967a982d89SJeremy L. Thompson 4977a982d89SJeremy L. Thompson int CeedOperatorFieldGetVector(CeedOperatorField opfield, CeedVector *vec) { 4987a982d89SJeremy L. Thompson *vec = opfield->vec; 4997a982d89SJeremy L. Thompson return 0; 5007a982d89SJeremy L. Thompson } 5017a982d89SJeremy L. Thompson 5027a982d89SJeremy L. Thompson /// @} 5037a982d89SJeremy L. Thompson 5047a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 5057a982d89SJeremy L. Thompson /// CeedOperator Public API 5067a982d89SJeremy L. Thompson /// ---------------------------------------------------------------------------- 5077a982d89SJeremy L. Thompson /// @addtogroup CeedOperatorUser 508dfdf5a53Sjeremylt /// @{ 509d7b241e6Sjeremylt 510d7b241e6Sjeremylt /** 5110219ea01SJeremy L Thompson @brief Create a CeedOperator and associate a CeedQFunction. A CeedBasis and 5120219ea01SJeremy L Thompson CeedElemRestriction can be associated with CeedQFunction fields with 5130219ea01SJeremy L Thompson \ref CeedOperatorSetField. 514d7b241e6Sjeremylt 515b11c1e72Sjeremylt @param ceed A Ceed object where the CeedOperator will be created 516d7b241e6Sjeremylt @param qf QFunction defining the action of the operator at quadrature points 51734138859Sjeremylt @param dqf QFunction defining the action of the Jacobian of @a qf (or 5184cc79fe7SJed Brown @ref CEED_QFUNCTION_NONE) 519d7b241e6Sjeremylt @param dqfT QFunction defining the action of the transpose of the Jacobian 5204cc79fe7SJed Brown of @a qf (or @ref CEED_QFUNCTION_NONE) 521b11c1e72Sjeremylt @param[out] op Address of the variable where the newly created 522b11c1e72Sjeremylt CeedOperator will be stored 523b11c1e72Sjeremylt 524b11c1e72Sjeremylt @return An error code: 0 - success, otherwise - failure 525dfdf5a53Sjeremylt 5267a982d89SJeremy L. Thompson @ref User 527d7b241e6Sjeremylt */ 528d7b241e6Sjeremylt int CeedOperatorCreate(Ceed ceed, CeedQFunction qf, CeedQFunction dqf, 529d7b241e6Sjeremylt CeedQFunction dqfT, CeedOperator *op) { 530d7b241e6Sjeremylt int ierr; 531d7b241e6Sjeremylt 5325fe0d4faSjeremylt if (!ceed->OperatorCreate) { 5335fe0d4faSjeremylt Ceed delegate; 534aefd8378Sjeremylt ierr = CeedGetObjectDelegate(ceed, &delegate, "Operator"); CeedChk(ierr); 5355fe0d4faSjeremylt 5365fe0d4faSjeremylt if (!delegate) 537c042f62fSJeremy L Thompson // LCOV_EXCL_START 5385fe0d4faSjeremylt return CeedError(ceed, 1, "Backend does not support OperatorCreate"); 539c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 5405fe0d4faSjeremylt 5415fe0d4faSjeremylt ierr = CeedOperatorCreate(delegate, qf, dqf, dqfT, op); CeedChk(ierr); 5425fe0d4faSjeremylt return 0; 5435fe0d4faSjeremylt } 5445fe0d4faSjeremylt 545b3b7035fSJeremy L Thompson if (!qf || qf == CEED_QFUNCTION_NONE) 546b3b7035fSJeremy L Thompson // LCOV_EXCL_START 547b3b7035fSJeremy L Thompson return CeedError(ceed, 1, "Operator must have a valid QFunction."); 548b3b7035fSJeremy L Thompson // LCOV_EXCL_STOP 549d7b241e6Sjeremylt ierr = CeedCalloc(1, op); CeedChk(ierr); 550d7b241e6Sjeremylt (*op)->ceed = ceed; 551d7b241e6Sjeremylt ceed->refcount++; 552d7b241e6Sjeremylt (*op)->refcount = 1; 553d7b241e6Sjeremylt (*op)->qf = qf; 554d7b241e6Sjeremylt qf->refcount++; 555442e7f0bSjeremylt if (dqf && dqf != CEED_QFUNCTION_NONE) { 556d7b241e6Sjeremylt (*op)->dqf = dqf; 557442e7f0bSjeremylt dqf->refcount++; 558442e7f0bSjeremylt } 559442e7f0bSjeremylt if (dqfT && dqfT != CEED_QFUNCTION_NONE) { 560d7b241e6Sjeremylt (*op)->dqfT = dqfT; 561442e7f0bSjeremylt dqfT->refcount++; 562442e7f0bSjeremylt } 563fe2413ffSjeremylt ierr = CeedCalloc(16, &(*op)->inputfields); CeedChk(ierr); 564fe2413ffSjeremylt ierr = CeedCalloc(16, &(*op)->outputfields); CeedChk(ierr); 565d7b241e6Sjeremylt ierr = ceed->OperatorCreate(*op); CeedChk(ierr); 566d7b241e6Sjeremylt return 0; 567d7b241e6Sjeremylt } 568d7b241e6Sjeremylt 569d7b241e6Sjeremylt /** 57052d6035fSJeremy L Thompson @brief Create an operator that composes the action of several operators 57152d6035fSJeremy L Thompson 57252d6035fSJeremy L Thompson @param ceed A Ceed object where the CeedOperator will be created 57352d6035fSJeremy L Thompson @param[out] op Address of the variable where the newly created 57452d6035fSJeremy L Thompson Composite CeedOperator will be stored 57552d6035fSJeremy L Thompson 57652d6035fSJeremy L Thompson @return An error code: 0 - success, otherwise - failure 57752d6035fSJeremy L Thompson 5787a982d89SJeremy L. Thompson @ref User 57952d6035fSJeremy L Thompson */ 58052d6035fSJeremy L Thompson int CeedCompositeOperatorCreate(Ceed ceed, CeedOperator *op) { 58152d6035fSJeremy L Thompson int ierr; 58252d6035fSJeremy L Thompson 58352d6035fSJeremy L Thompson if (!ceed->CompositeOperatorCreate) { 58452d6035fSJeremy L Thompson Ceed delegate; 585aefd8378Sjeremylt ierr = CeedGetObjectDelegate(ceed, &delegate, "Operator"); CeedChk(ierr); 58652d6035fSJeremy L Thompson 587250756a7Sjeremylt if (delegate) { 58852d6035fSJeremy L Thompson ierr = CeedCompositeOperatorCreate(delegate, op); CeedChk(ierr); 58952d6035fSJeremy L Thompson return 0; 59052d6035fSJeremy L Thompson } 591250756a7Sjeremylt } 59252d6035fSJeremy L Thompson 59352d6035fSJeremy L Thompson ierr = CeedCalloc(1, op); CeedChk(ierr); 59452d6035fSJeremy L Thompson (*op)->ceed = ceed; 59552d6035fSJeremy L Thompson ceed->refcount++; 59652d6035fSJeremy L Thompson (*op)->composite = true; 59752d6035fSJeremy L Thompson ierr = CeedCalloc(16, &(*op)->suboperators); CeedChk(ierr); 598250756a7Sjeremylt 599250756a7Sjeremylt if (ceed->CompositeOperatorCreate) { 60052d6035fSJeremy L Thompson ierr = ceed->CompositeOperatorCreate(*op); CeedChk(ierr); 601250756a7Sjeremylt } 60252d6035fSJeremy L Thompson return 0; 60352d6035fSJeremy L Thompson } 60452d6035fSJeremy L Thompson 60552d6035fSJeremy L Thompson /** 606b11c1e72Sjeremylt @brief Provide a field to a CeedOperator for use by its CeedQFunction 607d7b241e6Sjeremylt 608d7b241e6Sjeremylt This function is used to specify both active and passive fields to a 609d7b241e6Sjeremylt CeedOperator. For passive fields, a vector @arg v must be provided. Passive 610d7b241e6Sjeremylt fields can inputs or outputs (updated in-place when operator is applied). 611d7b241e6Sjeremylt 612d7b241e6Sjeremylt Active fields must be specified using this function, but their data (in a 613d7b241e6Sjeremylt CeedVector) is passed in CeedOperatorApply(). There can be at most one active 614d7b241e6Sjeremylt input and at most one active output. 615d7b241e6Sjeremylt 6168c91a0c9SJeremy L Thompson @param op CeedOperator on which to provide the field 6178795c945Sjeremylt @param fieldname Name of the field (to be matched with the name used by 6188795c945Sjeremylt CeedQFunction) 619b11c1e72Sjeremylt @param r CeedElemRestriction 6204cc79fe7SJed Brown @param b CeedBasis in which the field resides or @ref CEED_BASIS_COLLOCATED 621b11c1e72Sjeremylt if collocated with quadrature points 6224cc79fe7SJed Brown @param v CeedVector to be used by CeedOperator or @ref CEED_VECTOR_ACTIVE 6234cc79fe7SJed Brown if field is active or @ref CEED_VECTOR_NONE if using 6244cc79fe7SJed Brown @ref CEED_EVAL_WEIGHT in the QFunction 625b11c1e72Sjeremylt 626b11c1e72Sjeremylt @return An error code: 0 - success, otherwise - failure 627dfdf5a53Sjeremylt 6287a982d89SJeremy L. Thompson @ref User 629b11c1e72Sjeremylt **/ 630d7b241e6Sjeremylt int CeedOperatorSetField(CeedOperator op, const char *fieldname, 631a8d32208Sjeremylt CeedElemRestriction r, CeedBasis b, CeedVector v) { 632d7b241e6Sjeremylt int ierr; 63352d6035fSJeremy L Thompson if (op->composite) 634c042f62fSJeremy L Thompson // LCOV_EXCL_START 63552d6035fSJeremy L Thompson return CeedError(op->ceed, 1, "Cannot add field to composite operator."); 636c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 6378b067b84SJed Brown if (!r) 638c042f62fSJeremy L Thompson // LCOV_EXCL_START 639c042f62fSJeremy L Thompson return CeedError(op->ceed, 1, 640c042f62fSJeremy L Thompson "ElemRestriction r for field \"%s\" must be non-NULL.", 641c042f62fSJeremy L Thompson fieldname); 642c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 6438b067b84SJed Brown if (!b) 644c042f62fSJeremy L Thompson // LCOV_EXCL_START 645c042f62fSJeremy L Thompson return CeedError(op->ceed, 1, "Basis b for field \"%s\" must be non-NULL.", 646c042f62fSJeremy L Thompson fieldname); 647c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 6488b067b84SJed Brown if (!v) 649c042f62fSJeremy L Thompson // LCOV_EXCL_START 650c042f62fSJeremy L Thompson return CeedError(op->ceed, 1, "Vector v for field \"%s\" must be non-NULL.", 651c042f62fSJeremy L Thompson fieldname); 652c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 65352d6035fSJeremy L Thompson 654d7b241e6Sjeremylt CeedInt numelements; 655d7b241e6Sjeremylt ierr = CeedElemRestrictionGetNumElements(r, &numelements); CeedChk(ierr); 65615910d16Sjeremylt if (r != CEED_ELEMRESTRICTION_NONE && op->hasrestriction && 65715910d16Sjeremylt op->numelements != numelements) 658c042f62fSJeremy L Thompson // LCOV_EXCL_START 659d7b241e6Sjeremylt return CeedError(op->ceed, 1, 6601d102b48SJeremy L Thompson "ElemRestriction with %d elements incompatible with prior " 6611d102b48SJeremy L Thompson "%d elements", numelements, op->numelements); 662c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 66315910d16Sjeremylt if (r != CEED_ELEMRESTRICTION_NONE) { 664d7b241e6Sjeremylt op->numelements = numelements; 6652cb0afc5Sjeremylt op->hasrestriction = true; // Restriction set, but numelements may be 0 66615910d16Sjeremylt } 667d7b241e6Sjeremylt 668783c99b3SValeria Barra if (b != CEED_BASIS_COLLOCATED) { 669d7b241e6Sjeremylt CeedInt numqpoints; 670d7b241e6Sjeremylt ierr = CeedBasisGetNumQuadraturePoints(b, &numqpoints); CeedChk(ierr); 671d7b241e6Sjeremylt if (op->numqpoints && op->numqpoints != numqpoints) 672c042f62fSJeremy L Thompson // LCOV_EXCL_START 6731d102b48SJeremy L Thompson return CeedError(op->ceed, 1, "Basis with %d quadrature points " 6741d102b48SJeremy L Thompson "incompatible with prior %d points", numqpoints, 6751d102b48SJeremy L Thompson op->numqpoints); 676c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 677d7b241e6Sjeremylt op->numqpoints = numqpoints; 678d7b241e6Sjeremylt } 67915910d16Sjeremylt CeedQFunctionField qfield; 680d1bcdac9Sjeremylt CeedOperatorField *ofield; 681d7b241e6Sjeremylt for (CeedInt i=0; i<op->qf->numinputfields; i++) { 682fe2413ffSjeremylt if (!strcmp(fieldname, (*op->qf->inputfields[i]).fieldname)) { 68315910d16Sjeremylt qfield = op->qf->inputfields[i]; 684d7b241e6Sjeremylt ofield = &op->inputfields[i]; 685d7b241e6Sjeremylt goto found; 686d7b241e6Sjeremylt } 687d7b241e6Sjeremylt } 688d7b241e6Sjeremylt for (CeedInt i=0; i<op->qf->numoutputfields; i++) { 689fe2413ffSjeremylt if (!strcmp(fieldname, (*op->qf->outputfields[i]).fieldname)) { 69015910d16Sjeremylt qfield = op->qf->inputfields[i]; 691d7b241e6Sjeremylt ofield = &op->outputfields[i]; 692d7b241e6Sjeremylt goto found; 693d7b241e6Sjeremylt } 694d7b241e6Sjeremylt } 695c042f62fSJeremy L Thompson // LCOV_EXCL_START 696d7b241e6Sjeremylt return CeedError(op->ceed, 1, "QFunction has no knowledge of field '%s'", 697d7b241e6Sjeremylt fieldname); 698c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 699d7b241e6Sjeremylt found: 70015910d16Sjeremylt if (r == CEED_ELEMRESTRICTION_NONE && qfield->emode != CEED_EVAL_WEIGHT) 701b0d62198Sjeremylt // LCOV_EXCL_START 70215910d16Sjeremylt return CeedError(op->ceed, 1, "CEED_ELEMRESTRICTION_NONE can only be used " 70315910d16Sjeremylt "for a field with eval mode CEED_EVAL_WEIGHT"); 704b0d62198Sjeremylt // LCOV_EXCL_STOP 705fe2413ffSjeremylt ierr = CeedCalloc(1, ofield); CeedChk(ierr); 706fe2413ffSjeremylt (*ofield)->Erestrict = r; 70771352a93Sjeremylt r->refcount += 1; 708fe2413ffSjeremylt (*ofield)->basis = b; 70971352a93Sjeremylt if (b != CEED_BASIS_COLLOCATED) 71071352a93Sjeremylt b->refcount += 1; 711fe2413ffSjeremylt (*ofield)->vec = v; 71271352a93Sjeremylt if (v != CEED_VECTOR_ACTIVE && v != CEED_VECTOR_NONE) 71371352a93Sjeremylt v->refcount += 1; 714d7b241e6Sjeremylt op->nfields += 1; 715d7b241e6Sjeremylt return 0; 716d7b241e6Sjeremylt } 717d7b241e6Sjeremylt 718d7b241e6Sjeremylt /** 71952d6035fSJeremy L Thompson @brief Add a sub-operator to a composite CeedOperator 720288c0443SJeremy L Thompson 72134138859Sjeremylt @param[out] compositeop Composite CeedOperator 72234138859Sjeremylt @param subop Sub-operator CeedOperator 72352d6035fSJeremy L Thompson 72452d6035fSJeremy L Thompson @return An error code: 0 - success, otherwise - failure 72552d6035fSJeremy L Thompson 7267a982d89SJeremy L. Thompson @ref User 72752d6035fSJeremy L Thompson */ 728692c2638Sjeremylt int CeedCompositeOperatorAddSub(CeedOperator compositeop, CeedOperator subop) { 72952d6035fSJeremy L Thompson if (!compositeop->composite) 730c042f62fSJeremy L Thompson // LCOV_EXCL_START 7311d102b48SJeremy L Thompson return CeedError(compositeop->ceed, 1, "CeedOperator is not a composite " 7321d102b48SJeremy L Thompson "operator"); 733c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 73452d6035fSJeremy L Thompson 73552d6035fSJeremy L Thompson if (compositeop->numsub == CEED_COMPOSITE_MAX) 736c042f62fSJeremy L Thompson // LCOV_EXCL_START 7371d102b48SJeremy L Thompson return CeedError(compositeop->ceed, 1, "Cannot add additional suboperators"); 738c042f62fSJeremy L Thompson // LCOV_EXCL_STOP 73952d6035fSJeremy L Thompson 74052d6035fSJeremy L Thompson compositeop->suboperators[compositeop->numsub] = subop; 74152d6035fSJeremy L Thompson subop->refcount++; 74252d6035fSJeremy L Thompson compositeop->numsub++; 74352d6035fSJeremy L Thompson return 0; 74452d6035fSJeremy L Thompson } 74552d6035fSJeremy L Thompson 74652d6035fSJeremy L Thompson /** 7471d102b48SJeremy L Thompson @brief Assemble a linear CeedQFunction associated with a CeedOperator 7481d102b48SJeremy L Thompson 7491d102b48SJeremy L Thompson This returns a CeedVector containing a matrix at each quadrature point 7501d102b48SJeremy L Thompson providing the action of the CeedQFunction associated with the CeedOperator. 7511d102b48SJeremy L Thompson The vector 'assembled' is of shape 7521d102b48SJeremy L Thompson [num_elements, num_input_fields, num_output_fields, num_quad_points] 7531d102b48SJeremy L Thompson and contains column-major matrices representing the action of the 7541d102b48SJeremy L Thompson CeedQFunction for a corresponding quadrature point on an element. Inputs and 7551d102b48SJeremy L Thompson outputs are in the order provided by the user when adding CeedOperator fields. 7564cc79fe7SJed Brown For example, a CeedQFunction with inputs 'u' and 'gradu' and outputs 'gradv' and 7571d102b48SJeremy L Thompson 'v', provided in that order, would result in an assembled QFunction that 7581d102b48SJeremy L Thompson consists of (1 + dim) x (dim + 1) matrices at each quadrature point acting 7591d102b48SJeremy L Thompson on the input [u, du_0, du_1] and producing the output [dv_0, dv_1, v]. 7601d102b48SJeremy L Thompson 7611d102b48SJeremy L Thompson @param op CeedOperator to assemble CeedQFunction 7621d102b48SJeremy L Thompson @param[out] assembled CeedVector to store assembled CeedQFunction at 7631d102b48SJeremy L Thompson quadrature points 7641d102b48SJeremy L Thompson @param[out] rstr CeedElemRestriction for CeedVector containing assembled 7651d102b48SJeremy L Thompson CeedQFunction 7661d102b48SJeremy L Thompson @param request Address of CeedRequest for non-blocking completion, else 7674cc79fe7SJed Brown @ref CEED_REQUEST_IMMEDIATE 7681d102b48SJeremy L Thompson 7691d102b48SJeremy L Thompson @return An error code: 0 - success, otherwise - failure 7701d102b48SJeremy L Thompson 7717a982d89SJeremy L. Thompson @ref User 7721d102b48SJeremy L Thompson **/ 7731d102b48SJeremy L Thompson int CeedOperatorAssembleLinearQFunction(CeedOperator op, CeedVector *assembled, 7747f823360Sjeremylt CeedElemRestriction *rstr, 7757f823360Sjeremylt CeedRequest *request) { 7761d102b48SJeremy L Thompson int ierr; 7771d102b48SJeremy L Thompson Ceed ceed = op->ceed; 778250756a7Sjeremylt ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 7791d102b48SJeremy L Thompson 7807172caadSjeremylt // Backend version 7815107b09fSJeremy L Thompson if (op->AssembleLinearQFunction) { 7821d102b48SJeremy L Thompson ierr = op->AssembleLinearQFunction(op, assembled, rstr, request); 7831d102b48SJeremy L Thompson CeedChk(ierr); 7845107b09fSJeremy L Thompson } else { 7855107b09fSJeremy L Thompson // Fallback to reference Ceed 7865107b09fSJeremy L Thompson if (!op->opfallback) { 7875107b09fSJeremy L Thompson ierr = CeedOperatorCreateFallback(op); CeedChk(ierr); 7885107b09fSJeremy L Thompson } 7895107b09fSJeremy L Thompson // Assemble 7905107b09fSJeremy L Thompson ierr = op->opfallback->AssembleLinearQFunction(op->opfallback, assembled, 7915107b09fSJeremy L Thompson rstr, request); CeedChk(ierr); 7925107b09fSJeremy L Thompson } 793250756a7Sjeremylt 7941d102b48SJeremy L Thompson return 0; 7951d102b48SJeremy L Thompson } 7961d102b48SJeremy L Thompson 7971d102b48SJeremy L Thompson /** 798d965c7a7SJeremy L Thompson @brief Assemble the diagonal of a square linear CeedOperator 799b7ec98d8SJeremy L Thompson 800b7ec98d8SJeremy L Thompson This returns a CeedVector containing the diagonal of a linear CeedOperator. 801b7ec98d8SJeremy L Thompson 802*c04a41a7SJeremy L Thompson Note: Currently only non-composite CeedOperators with a single field and 803*c04a41a7SJeremy L Thompson composite CeedOperators with single field sub-operators are supported. 804d965c7a7SJeremy L Thompson 805b7ec98d8SJeremy L Thompson @param op CeedOperator to assemble CeedQFunction 806b7ec98d8SJeremy L Thompson @param[out] assembled CeedVector to store assembled CeedOperator diagonal 807b7ec98d8SJeremy L Thompson @param request Address of CeedRequest for non-blocking completion, else 8084cc79fe7SJed Brown @ref CEED_REQUEST_IMMEDIATE 809b7ec98d8SJeremy L Thompson 810b7ec98d8SJeremy L Thompson @return An error code: 0 - success, otherwise - failure 811b7ec98d8SJeremy L Thompson 8127a982d89SJeremy L. Thompson @ref User 813b7ec98d8SJeremy L Thompson **/ 8147f823360Sjeremylt int CeedOperatorAssembleLinearDiagonal(CeedOperator op, CeedVector *assembled, 8157f823360Sjeremylt CeedRequest *request) { 816b7ec98d8SJeremy L Thompson int ierr; 817b7ec98d8SJeremy L Thompson Ceed ceed = op->ceed; 818250756a7Sjeremylt ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 819b7ec98d8SJeremy L Thompson 820b7ec98d8SJeremy L Thompson // Use backend version, if available 821b7ec98d8SJeremy L Thompson if (op->AssembleLinearDiagonal) { 822b7ec98d8SJeremy L Thompson ierr = op->AssembleLinearDiagonal(op, assembled, request); CeedChk(ierr); 8237172caadSjeremylt } else { 8247172caadSjeremylt // Fallback to reference Ceed 8257172caadSjeremylt if (!op->opfallback) { 8267172caadSjeremylt ierr = CeedOperatorCreateFallback(op); CeedChk(ierr); 827b7ec98d8SJeremy L Thompson } 8287172caadSjeremylt // Assemble 8297172caadSjeremylt ierr = op->opfallback->AssembleLinearDiagonal(op->opfallback, assembled, 8307172caadSjeremylt request); CeedChk(ierr); 831b7ec98d8SJeremy L Thompson } 832b7ec98d8SJeremy L Thompson 833b7ec98d8SJeremy L Thompson return 0; 834b7ec98d8SJeremy L Thompson } 835b7ec98d8SJeremy L Thompson 836b7ec98d8SJeremy L Thompson /** 837d965c7a7SJeremy L Thompson @brief Assemble the point block diagonal of a square linear CeedOperator 838d965c7a7SJeremy L Thompson 839d965c7a7SJeremy L Thompson This returns a CeedVector containing the point block diagonal of a linear 840d965c7a7SJeremy L Thompson CeedOperator. 841d965c7a7SJeremy L Thompson 842*c04a41a7SJeremy L Thompson Note: Currently only non-composite CeedOperators with a single field and 843*c04a41a7SJeremy L Thompson composite CeedOperators with single field sub-operators are supported. 844d965c7a7SJeremy L Thompson 845d965c7a7SJeremy L Thompson @param op CeedOperator to assemble CeedQFunction 846d965c7a7SJeremy L Thompson @param[out] assembled CeedVector to store assembled CeedOperator point block 847d965c7a7SJeremy L Thompson diagonal, provided in row-major form with an 848d965c7a7SJeremy L Thompson @a ncomp * @a ncomp block at each node. The dimensions 849d965c7a7SJeremy L Thompson of this vector are derived from the active vector 850d965c7a7SJeremy L Thompson for the CeedOperator. The array has shape 851d965c7a7SJeremy L Thompson [nodes, component out, component in]. 852d965c7a7SJeremy L Thompson @param request Address of CeedRequest for non-blocking completion, else 853d965c7a7SJeremy L Thompson CEED_REQUEST_IMMEDIATE 854d965c7a7SJeremy L Thompson 855d965c7a7SJeremy L Thompson @return An error code: 0 - success, otherwise - failure 856d965c7a7SJeremy L Thompson 857d965c7a7SJeremy L Thompson @ref User 858d965c7a7SJeremy L Thompson **/ 859d965c7a7SJeremy L Thompson int CeedOperatorAssembleLinearPointBlockDiagonal(CeedOperator op, 860d965c7a7SJeremy L Thompson CeedVector *assembled, 861d965c7a7SJeremy L Thompson CeedRequest *request) { 862d965c7a7SJeremy L Thompson int ierr; 863d965c7a7SJeremy L Thompson Ceed ceed = op->ceed; 864d965c7a7SJeremy L Thompson ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 865d965c7a7SJeremy L Thompson 866d965c7a7SJeremy L Thompson // Use backend version, if available 867d965c7a7SJeremy L Thompson if (op->AssembleLinearPointBlockDiagonal) { 868d965c7a7SJeremy L Thompson ierr = op->AssembleLinearPointBlockDiagonal(op, assembled, request); 869d965c7a7SJeremy L Thompson CeedChk(ierr); 870d965c7a7SJeremy L Thompson } else { 871d965c7a7SJeremy L Thompson // Fallback to reference Ceed 872d965c7a7SJeremy L Thompson if (!op->opfallback) { 873d965c7a7SJeremy L Thompson ierr = CeedOperatorCreateFallback(op); CeedChk(ierr); 874d965c7a7SJeremy L Thompson } 875d965c7a7SJeremy L Thompson // Assemble 876d965c7a7SJeremy L Thompson ierr = op->opfallback->AssembleLinearPointBlockDiagonal(op->opfallback, 877d965c7a7SJeremy L Thompson assembled, request); CeedChk(ierr); 878d965c7a7SJeremy L Thompson } 879d965c7a7SJeremy L Thompson 880d965c7a7SJeremy L Thompson return 0; 881d965c7a7SJeremy L Thompson } 882d965c7a7SJeremy L Thompson 883d965c7a7SJeremy L Thompson /** 8843bd813ffSjeremylt @brief Build a FDM based approximate inverse for each element for a 8853bd813ffSjeremylt CeedOperator 8863bd813ffSjeremylt 8873bd813ffSjeremylt This returns a CeedOperator and CeedVector to apply a Fast Diagonalization 8883bd813ffSjeremylt Method based approximate inverse. This function obtains the simultaneous 8893bd813ffSjeremylt diagonalization for the 1D mass and Laplacian operators, 8903bd813ffSjeremylt M = V^T V, K = V^T S V. 8913bd813ffSjeremylt The assembled QFunction is used to modify the eigenvalues from simultaneous 8923bd813ffSjeremylt diagonalization and obtain an approximate inverse of the form 8933bd813ffSjeremylt V^T S^hat V. The CeedOperator must be linear and non-composite. The 8943bd813ffSjeremylt associated CeedQFunction must therefore also be linear. 8953bd813ffSjeremylt 8963bd813ffSjeremylt @param op CeedOperator to create element inverses 8973bd813ffSjeremylt @param[out] fdminv CeedOperator to apply the action of a FDM based inverse 8983bd813ffSjeremylt for each element 8993bd813ffSjeremylt @param request Address of CeedRequest for non-blocking completion, else 9004cc79fe7SJed Brown @ref CEED_REQUEST_IMMEDIATE 9013bd813ffSjeremylt 9023bd813ffSjeremylt @return An error code: 0 - success, otherwise - failure 9033bd813ffSjeremylt 9047a982d89SJeremy L. Thompson @ref Backend 9053bd813ffSjeremylt **/ 9063bd813ffSjeremylt int CeedOperatorCreateFDMElementInverse(CeedOperator op, CeedOperator *fdminv, 9073bd813ffSjeremylt CeedRequest *request) { 9083bd813ffSjeremylt int ierr; 9093bd813ffSjeremylt Ceed ceed = op->ceed; 910713f43c3Sjeremylt ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 9113bd813ffSjeremylt 912713f43c3Sjeremylt // Use backend version, if available 913713f43c3Sjeremylt if (op->CreateFDMElementInverse) { 914713f43c3Sjeremylt ierr = op->CreateFDMElementInverse(op, fdminv, request); CeedChk(ierr); 915713f43c3Sjeremylt } else { 916713f43c3Sjeremylt // Fallback to reference Ceed 917713f43c3Sjeremylt if (!op->opfallback) { 918713f43c3Sjeremylt ierr = CeedOperatorCreateFallback(op); CeedChk(ierr); 9193bd813ffSjeremylt } 920713f43c3Sjeremylt // Assemble 921713f43c3Sjeremylt ierr = op->opfallback->CreateFDMElementInverse(op->opfallback, fdminv, 9223bd813ffSjeremylt request); CeedChk(ierr); 9233bd813ffSjeremylt } 9243bd813ffSjeremylt 9253bd813ffSjeremylt return 0; 9263bd813ffSjeremylt } 9273bd813ffSjeremylt 9287a982d89SJeremy L. Thompson /** 9297a982d89SJeremy L. Thompson @brief View a CeedOperator 9307a982d89SJeremy L. Thompson 9317a982d89SJeremy L. Thompson @param[in] op CeedOperator to view 9327a982d89SJeremy L. Thompson @param[in] stream Stream to write; typically stdout/stderr or a file 9337a982d89SJeremy L. Thompson 9347a982d89SJeremy L. Thompson @return Error code: 0 - success, otherwise - failure 9357a982d89SJeremy L. Thompson 9367a982d89SJeremy L. Thompson @ref User 9377a982d89SJeremy L. Thompson **/ 9387a982d89SJeremy L. Thompson int CeedOperatorView(CeedOperator op, FILE *stream) { 9397a982d89SJeremy L. Thompson int ierr; 9407a982d89SJeremy L. Thompson 9417a982d89SJeremy L. Thompson if (op->composite) { 9427a982d89SJeremy L. Thompson fprintf(stream, "Composite CeedOperator\n"); 9437a982d89SJeremy L. Thompson 9447a982d89SJeremy L. Thompson for (CeedInt i=0; i<op->numsub; i++) { 9457a982d89SJeremy L. Thompson fprintf(stream, " SubOperator [%d]:\n", i); 9467a982d89SJeremy L. Thompson ierr = CeedOperatorSingleView(op->suboperators[i], 1, stream); 9477a982d89SJeremy L. Thompson CeedChk(ierr); 9487a982d89SJeremy L. Thompson } 9497a982d89SJeremy L. Thompson } else { 9507a982d89SJeremy L. Thompson fprintf(stream, "CeedOperator\n"); 9517a982d89SJeremy L. Thompson ierr = CeedOperatorSingleView(op, 0, stream); CeedChk(ierr); 9527a982d89SJeremy L. Thompson } 9537a982d89SJeremy L. Thompson 9547a982d89SJeremy L. Thompson return 0; 9557a982d89SJeremy L. Thompson } 9563bd813ffSjeremylt 9573bd813ffSjeremylt /** 9583bd813ffSjeremylt @brief Apply CeedOperator to a vector 959d7b241e6Sjeremylt 960d7b241e6Sjeremylt This computes the action of the operator on the specified (active) input, 961d7b241e6Sjeremylt yielding its (active) output. All inputs and outputs must be specified using 962d7b241e6Sjeremylt CeedOperatorSetField(). 963d7b241e6Sjeremylt 964d7b241e6Sjeremylt @param op CeedOperator to apply 9654cc79fe7SJed Brown @param[in] in CeedVector containing input state or @ref CEED_VECTOR_NONE if 96634138859Sjeremylt there are no active inputs 967b11c1e72Sjeremylt @param[out] out CeedVector to store result of applying operator (must be 9684cc79fe7SJed Brown distinct from @a in) or @ref CEED_VECTOR_NONE if there are no 96934138859Sjeremylt active outputs 970d7b241e6Sjeremylt @param request Address of CeedRequest for non-blocking completion, else 9714cc79fe7SJed Brown @ref CEED_REQUEST_IMMEDIATE 972b11c1e72Sjeremylt 973b11c1e72Sjeremylt @return An error code: 0 - success, otherwise - failure 974dfdf5a53Sjeremylt 9757a982d89SJeremy L. Thompson @ref User 976b11c1e72Sjeremylt **/ 977692c2638Sjeremylt int CeedOperatorApply(CeedOperator op, CeedVector in, CeedVector out, 978692c2638Sjeremylt CeedRequest *request) { 979d7b241e6Sjeremylt int ierr; 980d7b241e6Sjeremylt Ceed ceed = op->ceed; 981250756a7Sjeremylt ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 982d7b241e6Sjeremylt 983250756a7Sjeremylt if (op->numelements) { 984250756a7Sjeremylt // Standard Operator 985cae8b89aSjeremylt if (op->Apply) { 986250756a7Sjeremylt ierr = op->Apply(op, in, out, request); CeedChk(ierr); 987cae8b89aSjeremylt } else { 988cae8b89aSjeremylt // Zero all output vectors 989250756a7Sjeremylt CeedQFunction qf = op->qf; 990cae8b89aSjeremylt for (CeedInt i=0; i<qf->numoutputfields; i++) { 991cae8b89aSjeremylt CeedVector vec = op->outputfields[i]->vec; 992cae8b89aSjeremylt if (vec == CEED_VECTOR_ACTIVE) 993cae8b89aSjeremylt vec = out; 994cae8b89aSjeremylt if (vec != CEED_VECTOR_NONE) { 995cae8b89aSjeremylt ierr = CeedVectorSetValue(vec, 0.0); CeedChk(ierr); 996cae8b89aSjeremylt } 997cae8b89aSjeremylt } 998250756a7Sjeremylt // Apply 999250756a7Sjeremylt ierr = op->ApplyAdd(op, in, out, request); CeedChk(ierr); 1000250756a7Sjeremylt } 1001250756a7Sjeremylt } else if (op->composite) { 1002250756a7Sjeremylt // Composite Operator 1003250756a7Sjeremylt if (op->ApplyComposite) { 1004250756a7Sjeremylt ierr = op->ApplyComposite(op, in, out, request); CeedChk(ierr); 1005250756a7Sjeremylt } else { 1006250756a7Sjeremylt CeedInt numsub; 1007250756a7Sjeremylt ierr = CeedOperatorGetNumSub(op, &numsub); CeedChk(ierr); 1008250756a7Sjeremylt CeedOperator *suboperators; 1009250756a7Sjeremylt ierr = CeedOperatorGetSubList(op, &suboperators); CeedChk(ierr); 1010250756a7Sjeremylt 1011250756a7Sjeremylt // Zero all output vectors 1012250756a7Sjeremylt if (out != CEED_VECTOR_NONE) { 1013cae8b89aSjeremylt ierr = CeedVectorSetValue(out, 0.0); CeedChk(ierr); 1014cae8b89aSjeremylt } 1015250756a7Sjeremylt for (CeedInt i=0; i<numsub; i++) { 1016250756a7Sjeremylt for (CeedInt j=0; j<suboperators[i]->qf->numoutputfields; j++) { 1017250756a7Sjeremylt CeedVector vec = suboperators[i]->outputfields[j]->vec; 1018250756a7Sjeremylt if (vec != CEED_VECTOR_ACTIVE && vec != CEED_VECTOR_NONE) { 1019250756a7Sjeremylt ierr = CeedVectorSetValue(vec, 0.0); CeedChk(ierr); 1020250756a7Sjeremylt } 1021250756a7Sjeremylt } 1022250756a7Sjeremylt } 1023250756a7Sjeremylt // Apply 1024250756a7Sjeremylt for (CeedInt i=0; i<op->numsub; i++) { 1025250756a7Sjeremylt ierr = CeedOperatorApplyAdd(op->suboperators[i], in, out, request); 1026cae8b89aSjeremylt CeedChk(ierr); 1027cae8b89aSjeremylt } 1028cae8b89aSjeremylt } 1029250756a7Sjeremylt } 1030250756a7Sjeremylt 1031cae8b89aSjeremylt return 0; 1032cae8b89aSjeremylt } 1033cae8b89aSjeremylt 1034cae8b89aSjeremylt /** 1035cae8b89aSjeremylt @brief Apply CeedOperator to a vector and add result to output vector 1036cae8b89aSjeremylt 1037cae8b89aSjeremylt This computes the action of the operator on the specified (active) input, 1038cae8b89aSjeremylt yielding its (active) output. All inputs and outputs must be specified using 1039cae8b89aSjeremylt CeedOperatorSetField(). 1040cae8b89aSjeremylt 1041cae8b89aSjeremylt @param op CeedOperator to apply 1042cae8b89aSjeremylt @param[in] in CeedVector containing input state or NULL if there are no 1043cae8b89aSjeremylt active inputs 1044cae8b89aSjeremylt @param[out] out CeedVector to sum in result of applying operator (must be 1045cae8b89aSjeremylt distinct from @a in) or NULL if there are no active outputs 1046cae8b89aSjeremylt @param request Address of CeedRequest for non-blocking completion, else 10474cc79fe7SJed Brown @ref CEED_REQUEST_IMMEDIATE 1048cae8b89aSjeremylt 1049cae8b89aSjeremylt @return An error code: 0 - success, otherwise - failure 1050cae8b89aSjeremylt 10517a982d89SJeremy L. Thompson @ref User 1052cae8b89aSjeremylt **/ 1053cae8b89aSjeremylt int CeedOperatorApplyAdd(CeedOperator op, CeedVector in, CeedVector out, 1054cae8b89aSjeremylt CeedRequest *request) { 1055cae8b89aSjeremylt int ierr; 1056cae8b89aSjeremylt Ceed ceed = op->ceed; 1057250756a7Sjeremylt ierr = CeedOperatorCheckReady(ceed, op); CeedChk(ierr); 1058cae8b89aSjeremylt 1059250756a7Sjeremylt if (op->numelements) { 1060250756a7Sjeremylt // Standard Operator 1061250756a7Sjeremylt ierr = op->ApplyAdd(op, in, out, request); CeedChk(ierr); 1062250756a7Sjeremylt } else if (op->composite) { 1063250756a7Sjeremylt // Composite Operator 1064250756a7Sjeremylt if (op->ApplyAddComposite) { 1065250756a7Sjeremylt ierr = op->ApplyAddComposite(op, in, out, request); CeedChk(ierr); 1066cae8b89aSjeremylt } else { 1067250756a7Sjeremylt CeedInt numsub; 1068250756a7Sjeremylt ierr = CeedOperatorGetNumSub(op, &numsub); CeedChk(ierr); 1069250756a7Sjeremylt CeedOperator *suboperators; 1070250756a7Sjeremylt ierr = CeedOperatorGetSubList(op, &suboperators); CeedChk(ierr); 1071250756a7Sjeremylt 1072250756a7Sjeremylt for (CeedInt i=0; i<numsub; i++) { 1073250756a7Sjeremylt ierr = CeedOperatorApplyAdd(suboperators[i], in, out, request); 1074cae8b89aSjeremylt CeedChk(ierr); 10751d7d2407SJeremy L Thompson } 1076250756a7Sjeremylt } 1077250756a7Sjeremylt } 1078250756a7Sjeremylt 1079d7b241e6Sjeremylt return 0; 1080d7b241e6Sjeremylt } 1081d7b241e6Sjeremylt 1082d7b241e6Sjeremylt /** 1083b11c1e72Sjeremylt @brief Destroy a CeedOperator 1084d7b241e6Sjeremylt 1085d7b241e6Sjeremylt @param op CeedOperator to destroy 1086b11c1e72Sjeremylt 1087b11c1e72Sjeremylt @return An error code: 0 - success, otherwise - failure 1088dfdf5a53Sjeremylt 10897a982d89SJeremy L. Thompson @ref User 1090b11c1e72Sjeremylt **/ 1091d7b241e6Sjeremylt int CeedOperatorDestroy(CeedOperator *op) { 1092d7b241e6Sjeremylt int ierr; 1093d7b241e6Sjeremylt 1094d7b241e6Sjeremylt if (!*op || --(*op)->refcount > 0) return 0; 1095d7b241e6Sjeremylt if ((*op)->Destroy) { 1096d7b241e6Sjeremylt ierr = (*op)->Destroy(*op); CeedChk(ierr); 1097d7b241e6Sjeremylt } 1098fe2413ffSjeremylt ierr = CeedDestroy(&(*op)->ceed); CeedChk(ierr); 1099fe2413ffSjeremylt // Free fields 11001d102b48SJeremy L Thompson for (int i=0; i<(*op)->nfields; i++) 110152d6035fSJeremy L Thompson if ((*op)->inputfields[i]) { 110215910d16Sjeremylt if ((*op)->inputfields[i]->Erestrict != CEED_ELEMRESTRICTION_NONE) { 110371352a93Sjeremylt ierr = CeedElemRestrictionDestroy(&(*op)->inputfields[i]->Erestrict); 110471352a93Sjeremylt CeedChk(ierr); 110515910d16Sjeremylt } 110671352a93Sjeremylt if ((*op)->inputfields[i]->basis != CEED_BASIS_COLLOCATED) { 110771352a93Sjeremylt ierr = CeedBasisDestroy(&(*op)->inputfields[i]->basis); CeedChk(ierr); 110871352a93Sjeremylt } 110971352a93Sjeremylt if ((*op)->inputfields[i]->vec != CEED_VECTOR_ACTIVE && 111071352a93Sjeremylt (*op)->inputfields[i]->vec != CEED_VECTOR_NONE ) { 111171352a93Sjeremylt ierr = CeedVectorDestroy(&(*op)->inputfields[i]->vec); CeedChk(ierr); 111271352a93Sjeremylt } 1113fe2413ffSjeremylt ierr = CeedFree(&(*op)->inputfields[i]); CeedChk(ierr); 1114fe2413ffSjeremylt } 11151d102b48SJeremy L Thompson for (int i=0; i<(*op)->nfields; i++) 1116d0eb30a5Sjeremylt if ((*op)->outputfields[i]) { 111771352a93Sjeremylt ierr = CeedElemRestrictionDestroy(&(*op)->outputfields[i]->Erestrict); 111871352a93Sjeremylt CeedChk(ierr); 111971352a93Sjeremylt if ((*op)->outputfields[i]->basis != CEED_BASIS_COLLOCATED) { 112071352a93Sjeremylt ierr = CeedBasisDestroy(&(*op)->outputfields[i]->basis); CeedChk(ierr); 112171352a93Sjeremylt } 112271352a93Sjeremylt if ((*op)->outputfields[i]->vec != CEED_VECTOR_ACTIVE && 112371352a93Sjeremylt (*op)->outputfields[i]->vec != CEED_VECTOR_NONE ) { 112471352a93Sjeremylt ierr = CeedVectorDestroy(&(*op)->outputfields[i]->vec); CeedChk(ierr); 112571352a93Sjeremylt } 1126fe2413ffSjeremylt ierr = CeedFree(&(*op)->outputfields[i]); CeedChk(ierr); 1127fe2413ffSjeremylt } 112852d6035fSJeremy L Thompson // Destroy suboperators 11291d102b48SJeremy L Thompson for (int i=0; i<(*op)->numsub; i++) 113052d6035fSJeremy L Thompson if ((*op)->suboperators[i]) { 113152d6035fSJeremy L Thompson ierr = CeedOperatorDestroy(&(*op)->suboperators[i]); CeedChk(ierr); 113252d6035fSJeremy L Thompson } 1133d7b241e6Sjeremylt ierr = CeedQFunctionDestroy(&(*op)->qf); CeedChk(ierr); 1134d7b241e6Sjeremylt ierr = CeedQFunctionDestroy(&(*op)->dqf); CeedChk(ierr); 1135d7b241e6Sjeremylt ierr = CeedQFunctionDestroy(&(*op)->dqfT); CeedChk(ierr); 1136fe2413ffSjeremylt 11375107b09fSJeremy L Thompson // Destroy fallback 11385107b09fSJeremy L Thompson if ((*op)->opfallback) { 11395107b09fSJeremy L Thompson ierr = (*op)->qffallback->Destroy((*op)->qffallback); CeedChk(ierr); 11405107b09fSJeremy L Thompson ierr = CeedFree(&(*op)->qffallback); CeedChk(ierr); 11415107b09fSJeremy L Thompson ierr = (*op)->opfallback->Destroy((*op)->opfallback); CeedChk(ierr); 11425107b09fSJeremy L Thompson ierr = CeedFree(&(*op)->opfallback); CeedChk(ierr); 11435107b09fSJeremy L Thompson } 11445107b09fSJeremy L Thompson 1145fe2413ffSjeremylt ierr = CeedFree(&(*op)->inputfields); CeedChk(ierr); 1146fe2413ffSjeremylt ierr = CeedFree(&(*op)->outputfields); CeedChk(ierr); 114752d6035fSJeremy L Thompson ierr = CeedFree(&(*op)->suboperators); CeedChk(ierr); 1148d7b241e6Sjeremylt ierr = CeedFree(op); CeedChk(ierr); 1149d7b241e6Sjeremylt return 0; 1150d7b241e6Sjeremylt } 1151d7b241e6Sjeremylt 1152d7b241e6Sjeremylt /// @} 1153