13d8e8822SJeremy L Thompson // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors. 23d8e8822SJeremy L Thompson // All Rights Reserved. See the top-level LICENSE and NOTICE files for details. 377841947SLeila Ghaffari // 43d8e8822SJeremy L Thompson // SPDX-License-Identifier: BSD-2-Clause 577841947SLeila Ghaffari // 63d8e8822SJeremy L Thompson // This file is part of CEED: http://github.com/ceed 777841947SLeila Ghaffari 877841947SLeila Ghaffari /// @file 977841947SLeila Ghaffari /// Setup libCEED for Navier-Stokes example using PETSc 1077841947SLeila Ghaffari 1177841947SLeila Ghaffari #include "../navierstokes.h" 1277841947SLeila Ghaffari 1377841947SLeila Ghaffari // Utility function - essential BC dofs are encoded in closure indices as -(i+1). 1477841947SLeila Ghaffari PetscInt Involute(PetscInt i) { 1577841947SLeila Ghaffari return i >= 0 ? i : -(i+1); 1677841947SLeila Ghaffari } 1777841947SLeila Ghaffari 1877841947SLeila Ghaffari // Utility function to create local CEED restriction 197ed3e4cdSJeremy L Thompson PetscErrorCode CreateRestrictionFromPlex(Ceed ceed, DM dm, CeedInt height, 207ed3e4cdSJeremy L Thompson DMLabel domain_label, CeedInt value, CeedElemRestriction *elem_restr) { 217ed3e4cdSJeremy L Thompson PetscInt num_elem, elem_size, num_dof, num_comp, *elem_restr_offsets; 2277841947SLeila Ghaffari PetscErrorCode ierr; 237ed3e4cdSJeremy L Thompson 2477841947SLeila Ghaffari PetscFunctionBeginUser; 257ed3e4cdSJeremy L Thompson ierr = DMPlexGetLocalOffsets(dm, domain_label, value, height, 0, &num_elem, 267ed3e4cdSJeremy L Thompson &elem_size, &num_comp, &num_dof, &elem_restr_offsets); 277ed3e4cdSJeremy L Thompson CHKERRQ(ierr); 2877841947SLeila Ghaffari 297ed3e4cdSJeremy L Thompson CeedElemRestrictionCreate(ceed, num_elem, elem_size, num_comp, 307ed3e4cdSJeremy L Thompson 1, num_dof, CEED_MEM_HOST, CEED_COPY_VALUES, 317ed3e4cdSJeremy L Thompson elem_restr_offsets, elem_restr); 327ed3e4cdSJeremy L Thompson ierr = PetscFree(elem_restr_offsets); CHKERRQ(ierr); 3377841947SLeila Ghaffari 3477841947SLeila Ghaffari PetscFunctionReturn(0); 3577841947SLeila Ghaffari } 3677841947SLeila Ghaffari 3777841947SLeila Ghaffari // Utility function to get Ceed Restriction for each domain 3877841947SLeila Ghaffari PetscErrorCode GetRestrictionForDomain(Ceed ceed, DM dm, CeedInt height, 3977841947SLeila Ghaffari DMLabel domain_label, PetscInt value, 407ed3e4cdSJeremy L Thompson CeedInt Q, CeedInt q_data_size, 4177841947SLeila Ghaffari CeedElemRestriction *elem_restr_q, 4277841947SLeila Ghaffari CeedElemRestriction *elem_restr_x, 4377841947SLeila Ghaffari CeedElemRestriction *elem_restr_qd_i) { 4477841947SLeila Ghaffari DM dm_coord; 4577841947SLeila Ghaffari CeedInt dim, loc_num_elem; 4677841947SLeila Ghaffari CeedInt Q_dim; 472534dcc8SJed Brown CeedElemRestriction elem_restr_tmp; 4877841947SLeila Ghaffari PetscErrorCode ierr; 4977841947SLeila Ghaffari PetscFunctionBeginUser; 5077841947SLeila Ghaffari 5177841947SLeila Ghaffari ierr = DMGetDimension(dm, &dim); CHKERRQ(ierr); 5277841947SLeila Ghaffari dim -= height; 5377841947SLeila Ghaffari Q_dim = CeedIntPow(Q, dim); 542534dcc8SJed Brown ierr = CreateRestrictionFromPlex(ceed, dm, height, domain_label, value, 552534dcc8SJed Brown &elem_restr_tmp); 562534dcc8SJed Brown CHKERRQ(ierr); 572534dcc8SJed Brown if (elem_restr_q) *elem_restr_q = elem_restr_tmp; 582534dcc8SJed Brown if (elem_restr_x) { 5977841947SLeila Ghaffari ierr = DMGetCoordinateDM(dm, &dm_coord); CHKERRQ(ierr); 6077841947SLeila Ghaffari ierr = DMPlexSetClosurePermutationTensor(dm_coord, PETSC_DETERMINE, NULL); 6177841947SLeila Ghaffari CHKERRQ(ierr); 627ed3e4cdSJeremy L Thompson ierr = CreateRestrictionFromPlex(ceed, dm_coord, height, domain_label, value, 637ed3e4cdSJeremy L Thompson elem_restr_x); 647ed3e4cdSJeremy L Thompson CHKERRQ(ierr); 652534dcc8SJed Brown } 662534dcc8SJed Brown if (elem_restr_qd_i) { 672534dcc8SJed Brown CeedElemRestrictionGetNumElements(elem_restr_tmp, &loc_num_elem); 6877841947SLeila Ghaffari CeedElemRestrictionCreateStrided(ceed, loc_num_elem, Q_dim, 6977841947SLeila Ghaffari q_data_size, q_data_size*loc_num_elem*Q_dim, 7077841947SLeila Ghaffari CEED_STRIDES_BACKEND, elem_restr_qd_i); 712534dcc8SJed Brown } 722534dcc8SJed Brown if (!elem_restr_q) CeedElemRestrictionDestroy(&elem_restr_tmp); 7377841947SLeila Ghaffari PetscFunctionReturn(0); 7477841947SLeila Ghaffari } 7577841947SLeila Ghaffari 7677841947SLeila Ghaffari // Utility function to create CEED Composite Operator for the entire domain 7777841947SLeila Ghaffari PetscErrorCode CreateOperatorForDomain(Ceed ceed, DM dm, SimpleBC bc, 7877841947SLeila Ghaffari CeedData ceed_data, Physics phys, 79e334ad8fSJed Brown CeedOperator op_apply_vol, 80e334ad8fSJed Brown CeedOperator op_apply_ijacobian_vol, 81e334ad8fSJed Brown CeedInt height, 82e334ad8fSJed Brown CeedInt P_sur, CeedInt Q_sur, 83e334ad8fSJed Brown CeedInt q_data_size_sur, CeedInt jac_data_size_sur, 84e334ad8fSJed Brown CeedOperator *op_apply, CeedOperator *op_apply_ijacobian) { 852fe7aee7SLeila Ghaffari //CeedInt dim; 8677841947SLeila Ghaffari DMLabel domain_label; 8777841947SLeila Ghaffari PetscErrorCode ierr; 8877841947SLeila Ghaffari PetscFunctionBeginUser; 8977841947SLeila Ghaffari 9077841947SLeila Ghaffari // Create Composite Operaters 9177841947SLeila Ghaffari CeedCompositeOperatorCreate(ceed, op_apply); 92e334ad8fSJed Brown if (op_apply_ijacobian) 93e334ad8fSJed Brown CeedCompositeOperatorCreate(ceed, op_apply_ijacobian); 9477841947SLeila Ghaffari 9577841947SLeila Ghaffari // --Apply Sub-Operator for the volume 9677841947SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_vol); 97e334ad8fSJed Brown if (op_apply_ijacobian) 98e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_ijacobian_vol); 9977841947SLeila Ghaffari 10077841947SLeila Ghaffari // -- Create Sub-Operator for in/outflow BCs 10188626eedSJames Wright if (phys->has_neumann || 1) { 10277841947SLeila Ghaffari // --- Setup 10377841947SLeila Ghaffari ierr = DMGetLabel(dm, "Face Sets", &domain_label); CHKERRQ(ierr); 1042fe7aee7SLeila Ghaffari //ierr = DMGetDimension(dm, &dim); CHKERRQ(ierr); 10577841947SLeila Ghaffari 10677841947SLeila Ghaffari // --- Get number of quadrature points for the boundaries 10777841947SLeila Ghaffari CeedInt num_qpts_sur; 10877841947SLeila Ghaffari CeedBasisGetNumQuadraturePoints(ceed_data->basis_q_sur, &num_qpts_sur); 10977841947SLeila Ghaffari 1102fe7aee7SLeila Ghaffari // --- Create Sub-Operator for inflow boundaries 1112fe7aee7SLeila Ghaffari for (CeedInt i=0; i < bc->num_inflow; i++) { 11277841947SLeila Ghaffari CeedVector q_data_sur; 113e334ad8fSJed Brown CeedOperator op_setup_sur, op_apply_inflow, 114e334ad8fSJed Brown op_apply_inflow_jacobian = NULL; 1152fe7aee7SLeila Ghaffari CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur; 11677841947SLeila Ghaffari 11777841947SLeila Ghaffari // ---- CEED Restriction 1182fe7aee7SLeila Ghaffari ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->inflows[i], 1192fe7aee7SLeila Ghaffari Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 1202fe7aee7SLeila Ghaffari &elem_restr_qd_i_sur); 12177841947SLeila Ghaffari CHKERRQ(ierr); 12277841947SLeila Ghaffari 12377841947SLeila Ghaffari // ---- CEED Vector 12477841947SLeila Ghaffari PetscInt loc_num_elem_sur; 12577841947SLeila Ghaffari CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 12677841947SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 12777841947SLeila Ghaffari &q_data_sur); 12877841947SLeila Ghaffari 12977841947SLeila Ghaffari // ---- CEED Operator 13077841947SLeila Ghaffari // ----- CEED Operator for Setup (geometric factors) 13177841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 13277841947SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 133e6225c47SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 13477841947SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 13577841947SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_NONE); 136a61c78d6SJeremy L Thompson CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 13777841947SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 13877841947SLeila Ghaffari 13977841947SLeila Ghaffari // ----- CEED Operator for Physics 1402fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_apply_inflow, NULL, NULL, 1412fe7aee7SLeila Ghaffari &op_apply_inflow); 1422fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "q", elem_restr_q_sur, 14377841947SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 144*e8b03feeSJames Wright CeedOperatorSetField(op_apply_inflow, "Grad_q", elem_restr_q_sur, 145*e8b03feeSJames Wright ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 1462fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "surface qdata", elem_restr_qd_i_sur, 14777841947SLeila Ghaffari CEED_BASIS_COLLOCATED, q_data_sur); 1482fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "x", elem_restr_x_sur, 14977841947SLeila Ghaffari ceed_data->basis_x_sur, ceed_data->x_coord); 1502fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "v", elem_restr_q_sur, 15177841947SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 15277841947SLeila Ghaffari 153e334ad8fSJed Brown if (ceed_data->qf_apply_inflow_jacobian) { 154e334ad8fSJed Brown CeedOperatorCreate(ceed, ceed_data->qf_apply_inflow_jacobian, NULL, NULL, 155e334ad8fSJed Brown &op_apply_inflow_jacobian); 156e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "dq", elem_restr_q_sur, 157e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 158e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "surface qdata", 159e334ad8fSJed Brown elem_restr_qd_i_sur, 160e334ad8fSJed Brown CEED_BASIS_COLLOCATED, q_data_sur); 161e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "x", elem_restr_x_sur, 162e334ad8fSJed Brown ceed_data->basis_x_sur, ceed_data->x_coord); 163e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "v", elem_restr_q_sur, 164e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 165e334ad8fSJed Brown } 166e334ad8fSJed Brown 16777841947SLeila Ghaffari // ----- Apply CEED operator for Setup 16877841947SLeila Ghaffari CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 16977841947SLeila Ghaffari CEED_REQUEST_IMMEDIATE); 17077841947SLeila Ghaffari 1712fe7aee7SLeila Ghaffari // ----- Apply Sub-Operator for Physics 1722fe7aee7SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_inflow); 173e334ad8fSJed Brown if (op_apply_ijacobian) 174e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_inflow_jacobian); 17577841947SLeila Ghaffari 17677841947SLeila Ghaffari // ----- Cleanup 17777841947SLeila Ghaffari CeedVectorDestroy(&q_data_sur); 17877841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_q_sur); 17977841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_x_sur); 18077841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 18177841947SLeila Ghaffari CeedOperatorDestroy(&op_setup_sur); 1822fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_apply_inflow); 183e334ad8fSJed Brown CeedOperatorDestroy(&op_apply_inflow_jacobian); 18477841947SLeila Ghaffari } 18577841947SLeila Ghaffari 1862fe7aee7SLeila Ghaffari // --- Create Sub-Operator for outflow boundaries 1872fe7aee7SLeila Ghaffari for (CeedInt i=0; i < bc->num_outflow; i++) { 188e334ad8fSJed Brown CeedVector q_data_sur, jac_data_sur; 189e334ad8fSJed Brown CeedOperator op_setup_sur, op_apply_outflow, 190e334ad8fSJed Brown op_apply_outflow_jacobian = NULL; 191e334ad8fSJed Brown CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur, 192e334ad8fSJed Brown elem_restr_jd_i_sur; 1932fe7aee7SLeila Ghaffari 1942fe7aee7SLeila Ghaffari // ---- CEED Restriction 1952fe7aee7SLeila Ghaffari ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->outflows[i], 1962fe7aee7SLeila Ghaffari Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 1972fe7aee7SLeila Ghaffari &elem_restr_qd_i_sur); 1982fe7aee7SLeila Ghaffari CHKERRQ(ierr); 19939c69132SJed Brown if (jac_data_size_sur > 0) { 20039c69132SJed Brown // State-dependent data will be passed from residual to Jacobian. This will be collocated. 201e334ad8fSJed Brown ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->outflows[i], 202e334ad8fSJed Brown Q_sur, jac_data_size_sur, NULL, NULL, 203e334ad8fSJed Brown &elem_restr_jd_i_sur); 204e334ad8fSJed Brown CHKERRQ(ierr); 20539c69132SJed Brown CeedElemRestrictionCreateVector(elem_restr_jd_i_sur, &jac_data_sur, NULL); 20639c69132SJed Brown } else { 20739c69132SJed Brown elem_restr_jd_i_sur = NULL; 20839c69132SJed Brown jac_data_sur = NULL; 20939c69132SJed Brown } 2102fe7aee7SLeila Ghaffari 2112fe7aee7SLeila Ghaffari // ---- CEED Vector 2122fe7aee7SLeila Ghaffari PetscInt loc_num_elem_sur; 2132fe7aee7SLeila Ghaffari CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 2142fe7aee7SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 2152fe7aee7SLeila Ghaffari &q_data_sur); 2162fe7aee7SLeila Ghaffari 2172fe7aee7SLeila Ghaffari // ---- CEED Operator 2182fe7aee7SLeila Ghaffari // ----- CEED Operator for Setup (geometric factors) 2192fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 2202fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 2212fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 2222fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 2232fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_NONE); 2242fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 2252fe7aee7SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 2262fe7aee7SLeila Ghaffari 2272fe7aee7SLeila Ghaffari // ----- CEED Operator for Physics 2282fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_apply_outflow, NULL, NULL, 2292fe7aee7SLeila Ghaffari &op_apply_outflow); 2302fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "q", elem_restr_q_sur, 2312fe7aee7SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 232*e8b03feeSJames Wright CeedOperatorSetField(op_apply_outflow, "Grad_q", elem_restr_q_sur, 233*e8b03feeSJames Wright ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 2342fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "surface qdata", elem_restr_qd_i_sur, 2352fe7aee7SLeila Ghaffari CEED_BASIS_COLLOCATED, q_data_sur); 2362fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "x", elem_restr_x_sur, 2372fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, ceed_data->x_coord); 2382fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "v", elem_restr_q_sur, 2392fe7aee7SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 24039c69132SJed Brown if (elem_restr_jd_i_sur) 241e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow, "surface jacobian data", 242e334ad8fSJed Brown elem_restr_jd_i_sur, 243e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data_sur); 244e334ad8fSJed Brown 245e334ad8fSJed Brown if (ceed_data->qf_apply_outflow_jacobian) { 246e334ad8fSJed Brown CeedOperatorCreate(ceed, ceed_data->qf_apply_outflow_jacobian, NULL, NULL, 247e334ad8fSJed Brown &op_apply_outflow_jacobian); 248e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "dq", elem_restr_q_sur, 249e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 250e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "surface qdata", 251e334ad8fSJed Brown elem_restr_qd_i_sur, 252e334ad8fSJed Brown CEED_BASIS_COLLOCATED, q_data_sur); 253e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "surface jacobian data", 254e334ad8fSJed Brown elem_restr_jd_i_sur, 255e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data_sur); 256e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "v", elem_restr_q_sur, 257e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 258e334ad8fSJed Brown } 2592fe7aee7SLeila Ghaffari 2602fe7aee7SLeila Ghaffari // ----- Apply CEED operator for Setup 2612fe7aee7SLeila Ghaffari CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 2622fe7aee7SLeila Ghaffari CEED_REQUEST_IMMEDIATE); 2632fe7aee7SLeila Ghaffari 2642fe7aee7SLeila Ghaffari // ----- Apply Sub-Operator for Physics 2652fe7aee7SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_outflow); 266e334ad8fSJed Brown if (op_apply_ijacobian) 267e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_outflow_jacobian); 2682fe7aee7SLeila Ghaffari 2692fe7aee7SLeila Ghaffari // ----- Cleanup 2702fe7aee7SLeila Ghaffari CeedVectorDestroy(&q_data_sur); 271e334ad8fSJed Brown CeedVectorDestroy(&jac_data_sur); 2722fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_q_sur); 2732fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_x_sur); 2742fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 275e334ad8fSJed Brown CeedElemRestrictionDestroy(&elem_restr_jd_i_sur); 2762fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_setup_sur); 2772fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_apply_outflow); 278e334ad8fSJed Brown CeedOperatorDestroy(&op_apply_outflow_jacobian); 2792fe7aee7SLeila Ghaffari } 2802fe7aee7SLeila Ghaffari } 28111436a05SJames Wright 28211436a05SJames Wright // ----- Get Context Labels for Operator 28311436a05SJames Wright CeedOperatorContextGetFieldLabel(*op_apply, "solution time", 28411436a05SJames Wright &phys->solution_time_label); 28588626eedSJames Wright CeedOperatorContextGetFieldLabel(*op_apply, "timestep size", 28688626eedSJames Wright &phys->timestep_size_label); 28711436a05SJames Wright 28877841947SLeila Ghaffari PetscFunctionReturn(0); 28977841947SLeila Ghaffari } 29077841947SLeila Ghaffari 29177841947SLeila Ghaffari PetscErrorCode SetupLibceed(Ceed ceed, CeedData ceed_data, DM dm, User user, 292c95f9967SJed Brown AppCtx app_ctx, ProblemData *problem, SimpleBC bc) { 29377841947SLeila Ghaffari PetscErrorCode ierr; 29477841947SLeila Ghaffari PetscFunctionBeginUser; 29577841947SLeila Ghaffari 29677841947SLeila Ghaffari // ***************************************************************************** 29777841947SLeila Ghaffari // Set up CEED objects for the interior domain (volume) 29877841947SLeila Ghaffari // ***************************************************************************** 29977841947SLeila Ghaffari const PetscInt num_comp_q = 5; 30077841947SLeila Ghaffari const CeedInt dim = problem->dim, 30177841947SLeila Ghaffari num_comp_x = problem->dim, 30277841947SLeila Ghaffari q_data_size_vol = problem->q_data_size_vol, 303e334ad8fSJed Brown jac_data_size_vol = num_comp_q + 6 + 3, 30477841947SLeila Ghaffari P = app_ctx->degree + 1, 30577841947SLeila Ghaffari Q = P + app_ctx->q_extra; 306a3ae0734SJed Brown CeedElemRestriction elem_restr_jd_i; 307a3ae0734SJed Brown CeedVector jac_data; 30877841947SLeila Ghaffari 30977841947SLeila Ghaffari // ----------------------------------------------------------------------------- 31077841947SLeila Ghaffari // CEED Bases 31177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 31277841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_q, P, Q, CEED_GAUSS, 31377841947SLeila Ghaffari &ceed_data->basis_q); 31477841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, Q, CEED_GAUSS, 31577841947SLeila Ghaffari &ceed_data->basis_x); 31677841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, P, 31777841947SLeila Ghaffari CEED_GAUSS_LOBATTO, &ceed_data->basis_xc); 31877841947SLeila Ghaffari 31977841947SLeila Ghaffari // ----------------------------------------------------------------------------- 32077841947SLeila Ghaffari // CEED Restrictions 32177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 32277841947SLeila Ghaffari // -- Create restriction 3237ed3e4cdSJeremy L Thompson ierr = GetRestrictionForDomain(ceed, dm, 0, 0, 0, Q, q_data_size_vol, 3247ed3e4cdSJeremy L Thompson &ceed_data->elem_restr_q, &ceed_data->elem_restr_x, 32577841947SLeila Ghaffari &ceed_data->elem_restr_qd_i); CHKERRQ(ierr); 326a3ae0734SJed Brown 327a3ae0734SJed Brown ierr = GetRestrictionForDomain(ceed, dm, 0, 0, 0, Q, jac_data_size_vol, 328a3ae0734SJed Brown NULL, NULL, 329a3ae0734SJed Brown &elem_restr_jd_i); CHKERRQ(ierr); 33077841947SLeila Ghaffari // -- Create E vectors 33177841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_ceed, NULL); 33277841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_dot_ceed, 33377841947SLeila Ghaffari NULL); 33477841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->g_ceed, NULL); 33577841947SLeila Ghaffari 33677841947SLeila Ghaffari // ----------------------------------------------------------------------------- 33777841947SLeila Ghaffari // CEED QFunctions 33877841947SLeila Ghaffari // ----------------------------------------------------------------------------- 33977841947SLeila Ghaffari // -- Create QFunction for quadrature data 34091e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->setup_vol.qfunction, 34191e5af17SJed Brown problem->setup_vol.qfunction_loc, 34277841947SLeila Ghaffari &ceed_data->qf_setup_vol); 343841e4c73SJed Brown if (problem->setup_vol.qfunction_context) { 344841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_setup_vol, 345841e4c73SJed Brown problem->setup_vol.qfunction_context); 346841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->setup_vol.qfunction_context); 347841e4c73SJed Brown } 34877841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_vol, "dx", num_comp_x*dim, 34977841947SLeila Ghaffari CEED_EVAL_GRAD); 35077841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_vol, "weight", 1, CEED_EVAL_WEIGHT); 351a61c78d6SJeremy L Thompson CeedQFunctionAddOutput(ceed_data->qf_setup_vol, "qdata", q_data_size_vol, 35277841947SLeila Ghaffari CEED_EVAL_NONE); 35377841947SLeila Ghaffari 35477841947SLeila Ghaffari // -- Create QFunction for ICs 35591e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->ics.qfunction, 35691e5af17SJed Brown problem->ics.qfunction_loc, 35777841947SLeila Ghaffari &ceed_data->qf_ics); 358841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_ics, problem->ics.qfunction_context); 359841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->ics.qfunction_context); 36077841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ics, "x", num_comp_x, CEED_EVAL_INTERP); 36177841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_ics, "q0", num_comp_q, CEED_EVAL_NONE); 36277841947SLeila Ghaffari 36377841947SLeila Ghaffari // -- Create QFunction for RHS 36491e5af17SJed Brown if (problem->apply_vol_rhs.qfunction) { 36591e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_rhs.qfunction, 36691e5af17SJed Brown problem->apply_vol_rhs.qfunction_loc, &ceed_data->qf_rhs_vol); 367841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_rhs_vol, 368841e4c73SJed Brown problem->apply_vol_rhs.qfunction_context); 369841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_rhs.qfunction_context); 37077841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "q", num_comp_q, CEED_EVAL_INTERP); 371a3ae0734SJed Brown CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "Grad_q", num_comp_q*dim, 37277841947SLeila Ghaffari CEED_EVAL_GRAD); 373a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "qdata", q_data_size_vol, 37477841947SLeila Ghaffari CEED_EVAL_NONE); 37577841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "x", num_comp_x, CEED_EVAL_INTERP); 37677841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "v", num_comp_q, 37777841947SLeila Ghaffari CEED_EVAL_INTERP); 378a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "Grad_v", num_comp_q*dim, 37977841947SLeila Ghaffari CEED_EVAL_GRAD); 38077841947SLeila Ghaffari } 38177841947SLeila Ghaffari 38277841947SLeila Ghaffari // -- Create QFunction for IFunction 38391e5af17SJed Brown if (problem->apply_vol_ifunction.qfunction) { 38491e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ifunction.qfunction, 38591e5af17SJed Brown problem->apply_vol_ifunction.qfunction_loc, &ceed_data->qf_ifunction_vol); 386841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_ifunction_vol, 387841e4c73SJed Brown problem->apply_vol_ifunction.qfunction_context); 388841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_ifunction.qfunction_context); 38977841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q", num_comp_q, 39077841947SLeila Ghaffari CEED_EVAL_INTERP); 391a3ae0734SJed Brown CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "Grad_q", num_comp_q*dim, 39277841947SLeila Ghaffari CEED_EVAL_GRAD); 393a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q dot", num_comp_q, 39477841947SLeila Ghaffari CEED_EVAL_INTERP); 395a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "qdata", q_data_size_vol, 39677841947SLeila Ghaffari CEED_EVAL_NONE); 39777841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "x", num_comp_x, 39877841947SLeila Ghaffari CEED_EVAL_INTERP); 39977841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "v", num_comp_q, 40077841947SLeila Ghaffari CEED_EVAL_INTERP); 401a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "Grad_v", num_comp_q*dim, 40277841947SLeila Ghaffari CEED_EVAL_GRAD); 403a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "jac_data", 404a3ae0734SJed Brown jac_data_size_vol, CEED_EVAL_NONE); 40577841947SLeila Ghaffari } 40677841947SLeila Ghaffari 407e334ad8fSJed Brown CeedQFunction qf_ijacobian_vol = NULL; 408e334ad8fSJed Brown if (problem->apply_vol_ijacobian.qfunction) { 409e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ijacobian.qfunction, 410e334ad8fSJed Brown problem->apply_vol_ijacobian.qfunction_loc, &qf_ijacobian_vol); 411e334ad8fSJed Brown CeedQFunctionSetContext(qf_ijacobian_vol, 412e334ad8fSJed Brown problem->apply_vol_ijacobian.qfunction_context); 413e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_ijacobian.qfunction_context); 414e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "dq", num_comp_q, 415e334ad8fSJed Brown CEED_EVAL_INTERP); 416e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "Grad_dq", num_comp_q*dim, 417e334ad8fSJed Brown CEED_EVAL_GRAD); 418e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "qdata", q_data_size_vol, 419e334ad8fSJed Brown CEED_EVAL_NONE); 420e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "x", num_comp_x, 421e334ad8fSJed Brown CEED_EVAL_INTERP); 422e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "jac_data", 423e334ad8fSJed Brown jac_data_size_vol, CEED_EVAL_NONE); 424e334ad8fSJed Brown CeedQFunctionAddOutput(qf_ijacobian_vol, "v", num_comp_q, 425e334ad8fSJed Brown CEED_EVAL_INTERP); 426e334ad8fSJed Brown CeedQFunctionAddOutput(qf_ijacobian_vol, "Grad_v", num_comp_q*dim, 427e334ad8fSJed Brown CEED_EVAL_GRAD); 428e334ad8fSJed Brown } 429e334ad8fSJed Brown 43077841947SLeila Ghaffari // --------------------------------------------------------------------------- 43177841947SLeila Ghaffari // Element coordinates 43277841947SLeila Ghaffari // --------------------------------------------------------------------------- 43377841947SLeila Ghaffari // -- Create CEED vector 43477841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_x, &ceed_data->x_coord, 43577841947SLeila Ghaffari NULL); 43677841947SLeila Ghaffari 43777841947SLeila Ghaffari // -- Copy PETSc vector in CEED vector 43877841947SLeila Ghaffari Vec X_loc; 43977841947SLeila Ghaffari const PetscScalar *X_loc_array; 44077841947SLeila Ghaffari ierr = DMGetCoordinatesLocal(dm, &X_loc); CHKERRQ(ierr); 4411864f1c2SLeila Ghaffari ierr = VecScale(X_loc, problem->dm_scale); CHKERRQ(ierr); 44277841947SLeila Ghaffari ierr = VecGetArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 44377841947SLeila Ghaffari CeedVectorSetArray(ceed_data->x_coord, CEED_MEM_HOST, CEED_COPY_VALUES, 44477841947SLeila Ghaffari (PetscScalar *)X_loc_array); 44577841947SLeila Ghaffari ierr = VecRestoreArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 44677841947SLeila Ghaffari 44777841947SLeila Ghaffari // ----------------------------------------------------------------------------- 44877841947SLeila Ghaffari // CEED vectors 44977841947SLeila Ghaffari // ----------------------------------------------------------------------------- 45077841947SLeila Ghaffari // -- Create CEED vector for geometric data 45177841947SLeila Ghaffari CeedInt num_qpts_vol; 45277841947SLeila Ghaffari PetscInt loc_num_elem_vol; 45377841947SLeila Ghaffari CeedBasisGetNumQuadraturePoints(ceed_data->basis_q, &num_qpts_vol); 45477841947SLeila Ghaffari CeedElemRestrictionGetNumElements(ceed_data->elem_restr_q, &loc_num_elem_vol); 45577841947SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_vol*loc_num_elem_vol*num_qpts_vol, 45677841947SLeila Ghaffari &ceed_data->q_data); 45777841947SLeila Ghaffari 458a3ae0734SJed Brown CeedElemRestrictionCreateVector(elem_restr_jd_i, &jac_data, NULL); 45977841947SLeila Ghaffari // ----------------------------------------------------------------------------- 46077841947SLeila Ghaffari // CEED Operators 46177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 46277841947SLeila Ghaffari // -- Create CEED operator for quadrature data 46377841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_vol, NULL, NULL, 46477841947SLeila Ghaffari &ceed_data->op_setup_vol); 46577841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_setup_vol, "dx", ceed_data->elem_restr_x, 46677841947SLeila Ghaffari ceed_data->basis_x, CEED_VECTOR_ACTIVE); 46777841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_setup_vol, "weight", 468e6225c47SLeila Ghaffari CEED_ELEMRESTRICTION_NONE, ceed_data->basis_x, CEED_VECTOR_NONE); 469a61c78d6SJeremy L Thompson CeedOperatorSetField(ceed_data->op_setup_vol, "qdata", 470e6225c47SLeila Ghaffari ceed_data->elem_restr_qd_i, CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 47177841947SLeila Ghaffari 47277841947SLeila Ghaffari // -- Create CEED operator for ICs 47377841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_ics, NULL, NULL, &ceed_data->op_ics); 47477841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_ics, "x", ceed_data->elem_restr_x, 47577841947SLeila Ghaffari ceed_data->basis_xc, CEED_VECTOR_ACTIVE); 47677841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_ics, "q0", ceed_data->elem_restr_q, 47777841947SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 478841e4c73SJed Brown CeedOperatorContextGetFieldLabel(ceed_data->op_ics, "evaluation time", 479841e4c73SJed Brown &user->phys->ics_time_label); 48077841947SLeila Ghaffari 48177841947SLeila Ghaffari // Create CEED operator for RHS 48277841947SLeila Ghaffari if (ceed_data->qf_rhs_vol) { 48377841947SLeila Ghaffari CeedOperator op; 48477841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_rhs_vol, NULL, NULL, &op); 48577841947SLeila Ghaffari CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 48677841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 487a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_q", ceed_data->elem_restr_q, ceed_data->basis_q, 48877841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 489a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 490e6225c47SLeila Ghaffari CEED_BASIS_COLLOCATED, ceed_data->q_data); 49177841947SLeila Ghaffari CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 49277841947SLeila Ghaffari ceed_data->x_coord); 49377841947SLeila Ghaffari CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 49477841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 495a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 49677841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 49777841947SLeila Ghaffari user->op_rhs_vol = op; 49877841947SLeila Ghaffari } 49977841947SLeila Ghaffari 50077841947SLeila Ghaffari // -- CEED operator for IFunction 50177841947SLeila Ghaffari if (ceed_data->qf_ifunction_vol) { 50277841947SLeila Ghaffari CeedOperator op; 50377841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_ifunction_vol, NULL, NULL, &op); 50477841947SLeila Ghaffari CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 50577841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 506a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_q", ceed_data->elem_restr_q, ceed_data->basis_q, 50777841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 508a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "q dot", ceed_data->elem_restr_q, ceed_data->basis_q, 50977841947SLeila Ghaffari user->q_dot_ceed); 510a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 511e6225c47SLeila Ghaffari CEED_BASIS_COLLOCATED, ceed_data->q_data); 51277841947SLeila Ghaffari CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 51377841947SLeila Ghaffari ceed_data->x_coord); 51477841947SLeila Ghaffari CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 51577841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 516a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 51777841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 518a3ae0734SJed Brown CeedOperatorSetField(op, "jac_data", elem_restr_jd_i, 519a3ae0734SJed Brown CEED_BASIS_COLLOCATED, jac_data); 520a3ae0734SJed Brown 52177841947SLeila Ghaffari user->op_ifunction_vol = op; 52277841947SLeila Ghaffari } 52377841947SLeila Ghaffari 524e334ad8fSJed Brown CeedOperator op_ijacobian_vol = NULL; 525e334ad8fSJed Brown if (qf_ijacobian_vol) { 526e334ad8fSJed Brown CeedOperator op; 527e334ad8fSJed Brown CeedOperatorCreate(ceed, qf_ijacobian_vol, NULL, NULL, &op); 528e334ad8fSJed Brown CeedOperatorSetField(op, "dq", ceed_data->elem_restr_q, ceed_data->basis_q, 529e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 530e334ad8fSJed Brown CeedOperatorSetField(op, "Grad_dq", ceed_data->elem_restr_q, ceed_data->basis_q, 531e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 532e334ad8fSJed Brown CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 533e334ad8fSJed Brown CEED_BASIS_COLLOCATED, ceed_data->q_data); 534e334ad8fSJed Brown CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 535e334ad8fSJed Brown ceed_data->x_coord); 536e334ad8fSJed Brown CeedOperatorSetField(op, "jac_data", elem_restr_jd_i, 537e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data); 538e334ad8fSJed Brown CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 539e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 540e334ad8fSJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 541e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 542e334ad8fSJed Brown op_ijacobian_vol = op; 543e334ad8fSJed Brown CeedQFunctionDestroy(&qf_ijacobian_vol); 544e334ad8fSJed Brown } 545e334ad8fSJed Brown 54677841947SLeila Ghaffari // ***************************************************************************** 54777841947SLeila Ghaffari // Set up CEED objects for the exterior domain (surface) 54877841947SLeila Ghaffari // ***************************************************************************** 54977841947SLeila Ghaffari CeedInt height = 1, 55077841947SLeila Ghaffari dim_sur = dim - height, 55177841947SLeila Ghaffari P_sur = app_ctx->degree + 1, 55277841947SLeila Ghaffari Q_sur = P_sur + app_ctx->q_extra; 553e334ad8fSJed Brown const CeedInt q_data_size_sur = problem->q_data_size_sur, 554e334ad8fSJed Brown jac_data_size_sur = problem->jac_data_size_sur; 55577841947SLeila Ghaffari 55677841947SLeila Ghaffari // ----------------------------------------------------------------------------- 55777841947SLeila Ghaffari // CEED Bases 55877841947SLeila Ghaffari // ----------------------------------------------------------------------------- 55977841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_q, P_sur, Q_sur, 56077841947SLeila Ghaffari CEED_GAUSS, &ceed_data->basis_q_sur); 56177841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_x, 2, Q_sur, CEED_GAUSS, 56277841947SLeila Ghaffari &ceed_data->basis_x_sur); 56377841947SLeila Ghaffari 56477841947SLeila Ghaffari // ----------------------------------------------------------------------------- 56577841947SLeila Ghaffari // CEED QFunctions 56677841947SLeila Ghaffari // ----------------------------------------------------------------------------- 56777841947SLeila Ghaffari // -- Create QFunction for quadrature data 56891e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->setup_sur.qfunction, 56991e5af17SJed Brown problem->setup_sur.qfunction_loc, 57077841947SLeila Ghaffari &ceed_data->qf_setup_sur); 571841e4c73SJed Brown if (problem->setup_sur.qfunction_context) { 572841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_setup_sur, 573841e4c73SJed Brown problem->setup_sur.qfunction_context); 574841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->setup_sur.qfunction_context); 575841e4c73SJed Brown } 57677841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_sur, "dx", num_comp_x*dim_sur, 57777841947SLeila Ghaffari CEED_EVAL_GRAD); 57877841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_sur, "weight", 1, CEED_EVAL_WEIGHT); 579a61c78d6SJeremy L Thompson CeedQFunctionAddOutput(ceed_data->qf_setup_sur, "surface qdata", 5802fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 58177841947SLeila Ghaffari 5822fe7aee7SLeila Ghaffari // -- Creat QFunction for inflow boundaries 58391e5af17SJed Brown if (problem->apply_inflow.qfunction) { 58491e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_inflow.qfunction, 58591e5af17SJed Brown problem->apply_inflow.qfunction_loc, &ceed_data->qf_apply_inflow); 586841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_inflow, 587841e4c73SJed Brown problem->apply_inflow.qfunction_context); 588841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_inflow.qfunction_context); 5892fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "q", num_comp_q, 59077841947SLeila Ghaffari CEED_EVAL_INTERP); 591*e8b03feeSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "Grad_q", num_comp_q*(dim-1), 592*e8b03feeSJames Wright CEED_EVAL_GRAD); 5932fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "surface qdata", 5942fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 5952fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "x", num_comp_x, 59677841947SLeila Ghaffari CEED_EVAL_INTERP); 5972fe7aee7SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_apply_inflow, "v", num_comp_q, 5982fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 5992fe7aee7SLeila Ghaffari } 600e334ad8fSJed Brown if (problem->apply_inflow_jacobian.qfunction) { 601e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_inflow_jacobian.qfunction, 602e334ad8fSJed Brown problem->apply_inflow_jacobian.qfunction_loc, 603e334ad8fSJed Brown &ceed_data->qf_apply_inflow_jacobian); 604e334ad8fSJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_inflow_jacobian, 605e334ad8fSJed Brown problem->apply_inflow_jacobian.qfunction_context); 606e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_inflow_jacobian.qfunction_context); 607e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "dq", num_comp_q, 608e334ad8fSJed Brown CEED_EVAL_INTERP); 609e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "surface qdata", 610e334ad8fSJed Brown q_data_size_sur, CEED_EVAL_NONE); 611e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "x", num_comp_x, 612e334ad8fSJed Brown CEED_EVAL_INTERP); 613e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_inflow_jacobian, "v", num_comp_q, 614e334ad8fSJed Brown CEED_EVAL_INTERP); 615e334ad8fSJed Brown } 6162fe7aee7SLeila Ghaffari 6172fe7aee7SLeila Ghaffari // -- Creat QFunction for outflow boundaries 61891e5af17SJed Brown if (problem->apply_outflow.qfunction) { 61991e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_outflow.qfunction, 62091e5af17SJed Brown problem->apply_outflow.qfunction_loc, &ceed_data->qf_apply_outflow); 621841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_outflow, 622841e4c73SJed Brown problem->apply_outflow.qfunction_context); 623841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_outflow.qfunction_context); 6242fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "q", num_comp_q, 6252fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 626*e8b03feeSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "Grad_q", num_comp_q*(dim-1), 627*e8b03feeSJames Wright CEED_EVAL_GRAD); 6282fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "surface qdata", 6292fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 6302fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "x", num_comp_x, 6312fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 6322fe7aee7SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_apply_outflow, "v", num_comp_q, 63377841947SLeila Ghaffari CEED_EVAL_INTERP); 63439c69132SJed Brown if (jac_data_size_sur) 635e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_outflow, "surface jacobian data", 636e334ad8fSJed Brown jac_data_size_sur, 637e334ad8fSJed Brown CEED_EVAL_NONE); 638e334ad8fSJed Brown } 639e334ad8fSJed Brown if (problem->apply_outflow_jacobian.qfunction) { 640e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_outflow_jacobian.qfunction, 641e334ad8fSJed Brown problem->apply_outflow_jacobian.qfunction_loc, 642e334ad8fSJed Brown &ceed_data->qf_apply_outflow_jacobian); 643e334ad8fSJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_outflow_jacobian, 644e334ad8fSJed Brown problem->apply_outflow_jacobian.qfunction_context); 645e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_outflow_jacobian.qfunction_context); 646e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "dq", num_comp_q, 647e334ad8fSJed Brown CEED_EVAL_INTERP); 648e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "surface qdata", 649e334ad8fSJed Brown q_data_size_sur, CEED_EVAL_NONE); 650e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, 651e334ad8fSJed Brown "surface jacobian data", 652e334ad8fSJed Brown jac_data_size_sur, CEED_EVAL_NONE); 653e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_outflow_jacobian, "v", num_comp_q, 654e334ad8fSJed Brown CEED_EVAL_INTERP); 65577841947SLeila Ghaffari } 65677841947SLeila Ghaffari 65777841947SLeila Ghaffari // ***************************************************************************** 65877841947SLeila Ghaffari // CEED Operator Apply 65977841947SLeila Ghaffari // ***************************************************************************** 66077841947SLeila Ghaffari // -- Apply CEED Operator for the geometric data 66177841947SLeila Ghaffari CeedOperatorApply(ceed_data->op_setup_vol, ceed_data->x_coord, 66277841947SLeila Ghaffari ceed_data->q_data, CEED_REQUEST_IMMEDIATE); 66377841947SLeila Ghaffari 66477841947SLeila Ghaffari // -- Create and apply CEED Composite Operator for the entire domain 66577841947SLeila Ghaffari if (!user->phys->implicit) { // RHS 66677841947SLeila Ghaffari ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 667e334ad8fSJed Brown user->op_rhs_vol, NULL, height, P_sur, Q_sur, 668e334ad8fSJed Brown q_data_size_sur, 0, 669e334ad8fSJed Brown &user->op_rhs, NULL); CHKERRQ(ierr); 67077841947SLeila Ghaffari } else { // IFunction 67177841947SLeila Ghaffari ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 672e334ad8fSJed Brown user->op_ifunction_vol, op_ijacobian_vol, 673e334ad8fSJed Brown height, P_sur, Q_sur, 674e334ad8fSJed Brown q_data_size_sur, jac_data_size_sur, 675e334ad8fSJed Brown &user->op_ifunction, 676e334ad8fSJed Brown op_ijacobian_vol ? &user->op_ijacobian : NULL); CHKERRQ(ierr); 677e334ad8fSJed Brown if (user->op_ijacobian) { 678e334ad8fSJed Brown CeedOperatorContextGetFieldLabel(user->op_ijacobian, "ijacobian time shift", 679e334ad8fSJed Brown &user->phys->ijacobian_time_shift_label); 680e334ad8fSJed Brown } 681e334ad8fSJed Brown 68277841947SLeila Ghaffari } 683a3ae0734SJed Brown CeedElemRestrictionDestroy(&elem_restr_jd_i); 68422d320f5SLeila Ghaffari CeedOperatorDestroy(&op_ijacobian_vol); 685a3ae0734SJed Brown CeedVectorDestroy(&jac_data); 68677841947SLeila Ghaffari PetscFunctionReturn(0); 68777841947SLeila Ghaffari } 688