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) { 8577841947SLeila Ghaffari DMLabel domain_label; 8677841947SLeila Ghaffari PetscErrorCode ierr; 8777841947SLeila Ghaffari PetscFunctionBeginUser; 8877841947SLeila Ghaffari 8977841947SLeila Ghaffari // Create Composite Operaters 9077841947SLeila Ghaffari CeedCompositeOperatorCreate(ceed, op_apply); 91e334ad8fSJed Brown if (op_apply_ijacobian) 92e334ad8fSJed Brown CeedCompositeOperatorCreate(ceed, op_apply_ijacobian); 9377841947SLeila Ghaffari 9477841947SLeila Ghaffari // --Apply Sub-Operator for the volume 9577841947SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_vol); 96e334ad8fSJed Brown if (op_apply_ijacobian) 97e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_ijacobian_vol); 9877841947SLeila Ghaffari 9977841947SLeila Ghaffari // -- Create Sub-Operator for in/outflow BCs 10088626eedSJames Wright if (phys->has_neumann || 1) { 10177841947SLeila Ghaffari // --- Setup 10277841947SLeila Ghaffari ierr = DMGetLabel(dm, "Face Sets", &domain_label); CHKERRQ(ierr); 10377841947SLeila Ghaffari 10477841947SLeila Ghaffari // --- Get number of quadrature points for the boundaries 10577841947SLeila Ghaffari CeedInt num_qpts_sur; 10677841947SLeila Ghaffari CeedBasisGetNumQuadraturePoints(ceed_data->basis_q_sur, &num_qpts_sur); 10777841947SLeila Ghaffari 1082fe7aee7SLeila Ghaffari // --- Create Sub-Operator for inflow boundaries 1092fe7aee7SLeila Ghaffari for (CeedInt i=0; i < bc->num_inflow; i++) { 11077841947SLeila Ghaffari CeedVector q_data_sur; 111e334ad8fSJed Brown CeedOperator op_setup_sur, op_apply_inflow, 112e334ad8fSJed Brown op_apply_inflow_jacobian = NULL; 1132fe7aee7SLeila Ghaffari CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur; 11477841947SLeila Ghaffari 11577841947SLeila Ghaffari // ---- CEED Restriction 1162fe7aee7SLeila Ghaffari ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->inflows[i], 1172fe7aee7SLeila Ghaffari Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 1182fe7aee7SLeila Ghaffari &elem_restr_qd_i_sur); 11977841947SLeila Ghaffari CHKERRQ(ierr); 12077841947SLeila Ghaffari 12177841947SLeila Ghaffari // ---- CEED Vector 12277841947SLeila Ghaffari PetscInt loc_num_elem_sur; 12377841947SLeila Ghaffari CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 12477841947SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 12577841947SLeila Ghaffari &q_data_sur); 12677841947SLeila Ghaffari 12777841947SLeila Ghaffari // ---- CEED Operator 12877841947SLeila Ghaffari // ----- CEED Operator for Setup (geometric factors) 12977841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 13077841947SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 131e6225c47SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 13277841947SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 13377841947SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_NONE); 134a61c78d6SJeremy L Thompson CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 13577841947SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 13677841947SLeila Ghaffari 13777841947SLeila Ghaffari // ----- CEED Operator for Physics 1382fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_apply_inflow, NULL, NULL, 1392fe7aee7SLeila Ghaffari &op_apply_inflow); 1402fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "q", elem_restr_q_sur, 14177841947SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 142e8b03feeSJames Wright CeedOperatorSetField(op_apply_inflow, "Grad_q", elem_restr_q_sur, 143e8b03feeSJames Wright ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 1442fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "surface qdata", elem_restr_qd_i_sur, 14577841947SLeila Ghaffari CEED_BASIS_COLLOCATED, q_data_sur); 1462fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "x", elem_restr_x_sur, 14777841947SLeila Ghaffari ceed_data->basis_x_sur, ceed_data->x_coord); 1482fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_inflow, "v", elem_restr_q_sur, 14977841947SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 15077841947SLeila Ghaffari 151e334ad8fSJed Brown if (ceed_data->qf_apply_inflow_jacobian) { 152e334ad8fSJed Brown CeedOperatorCreate(ceed, ceed_data->qf_apply_inflow_jacobian, NULL, NULL, 153e334ad8fSJed Brown &op_apply_inflow_jacobian); 154e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "dq", elem_restr_q_sur, 155e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 156e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "surface qdata", 157e334ad8fSJed Brown elem_restr_qd_i_sur, 158e334ad8fSJed Brown CEED_BASIS_COLLOCATED, q_data_sur); 159e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "x", elem_restr_x_sur, 160e334ad8fSJed Brown ceed_data->basis_x_sur, ceed_data->x_coord); 161e334ad8fSJed Brown CeedOperatorSetField(op_apply_inflow_jacobian, "v", elem_restr_q_sur, 162e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 163e334ad8fSJed Brown } 164e334ad8fSJed Brown 16577841947SLeila Ghaffari // ----- Apply CEED operator for Setup 16677841947SLeila Ghaffari CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 16777841947SLeila Ghaffari CEED_REQUEST_IMMEDIATE); 16877841947SLeila Ghaffari 1692fe7aee7SLeila Ghaffari // ----- Apply Sub-Operator for Physics 1702fe7aee7SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_inflow); 171e334ad8fSJed Brown if (op_apply_ijacobian) 172e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_inflow_jacobian); 17377841947SLeila Ghaffari 17477841947SLeila Ghaffari // ----- Cleanup 17577841947SLeila Ghaffari CeedVectorDestroy(&q_data_sur); 17677841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_q_sur); 17777841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_x_sur); 17877841947SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 17977841947SLeila Ghaffari CeedOperatorDestroy(&op_setup_sur); 1802fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_apply_inflow); 181e334ad8fSJed Brown CeedOperatorDestroy(&op_apply_inflow_jacobian); 18277841947SLeila Ghaffari } 18377841947SLeila Ghaffari 1842fe7aee7SLeila Ghaffari // --- Create Sub-Operator for outflow boundaries 1852fe7aee7SLeila Ghaffari for (CeedInt i=0; i < bc->num_outflow; i++) { 186e334ad8fSJed Brown CeedVector q_data_sur, jac_data_sur; 187e334ad8fSJed Brown CeedOperator op_setup_sur, op_apply_outflow, 188e334ad8fSJed Brown op_apply_outflow_jacobian = NULL; 189e334ad8fSJed Brown CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur, 190e334ad8fSJed Brown elem_restr_jd_i_sur; 1912fe7aee7SLeila Ghaffari 1922fe7aee7SLeila Ghaffari // ---- CEED Restriction 1932fe7aee7SLeila Ghaffari ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->outflows[i], 1942fe7aee7SLeila Ghaffari Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 1952fe7aee7SLeila Ghaffari &elem_restr_qd_i_sur); 1962fe7aee7SLeila Ghaffari CHKERRQ(ierr); 19739c69132SJed Brown if (jac_data_size_sur > 0) { 19839c69132SJed Brown // State-dependent data will be passed from residual to Jacobian. This will be collocated. 199e334ad8fSJed Brown ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->outflows[i], 200e334ad8fSJed Brown Q_sur, jac_data_size_sur, NULL, NULL, 201e334ad8fSJed Brown &elem_restr_jd_i_sur); 202e334ad8fSJed Brown CHKERRQ(ierr); 20339c69132SJed Brown CeedElemRestrictionCreateVector(elem_restr_jd_i_sur, &jac_data_sur, NULL); 20439c69132SJed Brown } else { 20539c69132SJed Brown elem_restr_jd_i_sur = NULL; 20639c69132SJed Brown jac_data_sur = NULL; 20739c69132SJed Brown } 2082fe7aee7SLeila Ghaffari 2092fe7aee7SLeila Ghaffari // ---- CEED Vector 2102fe7aee7SLeila Ghaffari PetscInt loc_num_elem_sur; 2112fe7aee7SLeila Ghaffari CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 2122fe7aee7SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 2132fe7aee7SLeila Ghaffari &q_data_sur); 2142fe7aee7SLeila Ghaffari 2152fe7aee7SLeila Ghaffari // ---- CEED Operator 2162fe7aee7SLeila Ghaffari // ----- CEED Operator for Setup (geometric factors) 2172fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 2182fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 2192fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 2202fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 2212fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, CEED_VECTOR_NONE); 2222fe7aee7SLeila Ghaffari CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 2232fe7aee7SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 2242fe7aee7SLeila Ghaffari 2252fe7aee7SLeila Ghaffari // ----- CEED Operator for Physics 2262fe7aee7SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_apply_outflow, NULL, NULL, 2272fe7aee7SLeila Ghaffari &op_apply_outflow); 2282fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "q", elem_restr_q_sur, 2292fe7aee7SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 230e8b03feeSJames Wright CeedOperatorSetField(op_apply_outflow, "Grad_q", elem_restr_q_sur, 231e8b03feeSJames Wright ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 2322fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "surface qdata", elem_restr_qd_i_sur, 2332fe7aee7SLeila Ghaffari CEED_BASIS_COLLOCATED, q_data_sur); 2342fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "x", elem_restr_x_sur, 2352fe7aee7SLeila Ghaffari ceed_data->basis_x_sur, ceed_data->x_coord); 2362fe7aee7SLeila Ghaffari CeedOperatorSetField(op_apply_outflow, "v", elem_restr_q_sur, 2372fe7aee7SLeila Ghaffari ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 23839c69132SJed Brown if (elem_restr_jd_i_sur) 239e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow, "surface jacobian data", 240e334ad8fSJed Brown elem_restr_jd_i_sur, 241e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data_sur); 242e334ad8fSJed Brown 243e334ad8fSJed Brown if (ceed_data->qf_apply_outflow_jacobian) { 244e334ad8fSJed Brown CeedOperatorCreate(ceed, ceed_data->qf_apply_outflow_jacobian, NULL, NULL, 245e334ad8fSJed Brown &op_apply_outflow_jacobian); 246e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "dq", elem_restr_q_sur, 247e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 248*0ec2498eSJames Wright CeedOperatorSetField(op_apply_outflow_jacobian, "Grad_dq", elem_restr_q_sur, 249*0ec2498eSJames Wright 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); 253*0ec2498eSJames Wright CeedOperatorSetField(op_apply_outflow_jacobian, "x", elem_restr_x_sur, 254*0ec2498eSJames Wright ceed_data->basis_x_sur, ceed_data->x_coord); 255e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "surface jacobian data", 256e334ad8fSJed Brown elem_restr_jd_i_sur, 257e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data_sur); 258e334ad8fSJed Brown CeedOperatorSetField(op_apply_outflow_jacobian, "v", elem_restr_q_sur, 259e334ad8fSJed Brown ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 260e334ad8fSJed Brown } 2612fe7aee7SLeila Ghaffari 2622fe7aee7SLeila Ghaffari // ----- Apply CEED operator for Setup 2632fe7aee7SLeila Ghaffari CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 2642fe7aee7SLeila Ghaffari CEED_REQUEST_IMMEDIATE); 2652fe7aee7SLeila Ghaffari 2662fe7aee7SLeila Ghaffari // ----- Apply Sub-Operator for Physics 2672fe7aee7SLeila Ghaffari CeedCompositeOperatorAddSub(*op_apply, op_apply_outflow); 268e334ad8fSJed Brown if (op_apply_ijacobian) 269e334ad8fSJed Brown CeedCompositeOperatorAddSub(*op_apply_ijacobian, op_apply_outflow_jacobian); 2702fe7aee7SLeila Ghaffari 2712fe7aee7SLeila Ghaffari // ----- Cleanup 2722fe7aee7SLeila Ghaffari CeedVectorDestroy(&q_data_sur); 273e334ad8fSJed Brown CeedVectorDestroy(&jac_data_sur); 2742fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_q_sur); 2752fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_x_sur); 2762fe7aee7SLeila Ghaffari CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 277e334ad8fSJed Brown CeedElemRestrictionDestroy(&elem_restr_jd_i_sur); 2782fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_setup_sur); 2792fe7aee7SLeila Ghaffari CeedOperatorDestroy(&op_apply_outflow); 280e334ad8fSJed Brown CeedOperatorDestroy(&op_apply_outflow_jacobian); 2812fe7aee7SLeila Ghaffari } 2822fe7aee7SLeila Ghaffari } 28311436a05SJames Wright 28411436a05SJames Wright // ----- Get Context Labels for Operator 28511436a05SJames Wright CeedOperatorContextGetFieldLabel(*op_apply, "solution time", 28611436a05SJames Wright &phys->solution_time_label); 28788626eedSJames Wright CeedOperatorContextGetFieldLabel(*op_apply, "timestep size", 28888626eedSJames Wright &phys->timestep_size_label); 28911436a05SJames Wright 29077841947SLeila Ghaffari PetscFunctionReturn(0); 29177841947SLeila Ghaffari } 29277841947SLeila Ghaffari 29377841947SLeila Ghaffari PetscErrorCode SetupLibceed(Ceed ceed, CeedData ceed_data, DM dm, User user, 294c95f9967SJed Brown AppCtx app_ctx, ProblemData *problem, SimpleBC bc) { 29577841947SLeila Ghaffari PetscErrorCode ierr; 29677841947SLeila Ghaffari PetscFunctionBeginUser; 29777841947SLeila Ghaffari 29877841947SLeila Ghaffari // ***************************************************************************** 29977841947SLeila Ghaffari // Set up CEED objects for the interior domain (volume) 30077841947SLeila Ghaffari // ***************************************************************************** 30177841947SLeila Ghaffari const PetscInt num_comp_q = 5; 30277841947SLeila Ghaffari const CeedInt dim = problem->dim, 30377841947SLeila Ghaffari num_comp_x = problem->dim, 30477841947SLeila Ghaffari q_data_size_vol = problem->q_data_size_vol, 305e334ad8fSJed Brown jac_data_size_vol = num_comp_q + 6 + 3, 30677841947SLeila Ghaffari P = app_ctx->degree + 1, 30777841947SLeila Ghaffari Q = P + app_ctx->q_extra; 308a3ae0734SJed Brown CeedElemRestriction elem_restr_jd_i; 309a3ae0734SJed Brown CeedVector jac_data; 31077841947SLeila Ghaffari 31177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 31277841947SLeila Ghaffari // CEED Bases 31377841947SLeila Ghaffari // ----------------------------------------------------------------------------- 31477841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_q, P, Q, CEED_GAUSS, 31577841947SLeila Ghaffari &ceed_data->basis_q); 31677841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, Q, CEED_GAUSS, 31777841947SLeila Ghaffari &ceed_data->basis_x); 31877841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, P, 31977841947SLeila Ghaffari CEED_GAUSS_LOBATTO, &ceed_data->basis_xc); 32077841947SLeila Ghaffari 32177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 32277841947SLeila Ghaffari // CEED Restrictions 32377841947SLeila Ghaffari // ----------------------------------------------------------------------------- 32477841947SLeila Ghaffari // -- Create restriction 3257ed3e4cdSJeremy L Thompson ierr = GetRestrictionForDomain(ceed, dm, 0, 0, 0, Q, q_data_size_vol, 3267ed3e4cdSJeremy L Thompson &ceed_data->elem_restr_q, &ceed_data->elem_restr_x, 32777841947SLeila Ghaffari &ceed_data->elem_restr_qd_i); CHKERRQ(ierr); 328a3ae0734SJed Brown 329a3ae0734SJed Brown ierr = GetRestrictionForDomain(ceed, dm, 0, 0, 0, Q, jac_data_size_vol, 330a3ae0734SJed Brown NULL, NULL, 331a3ae0734SJed Brown &elem_restr_jd_i); CHKERRQ(ierr); 33277841947SLeila Ghaffari // -- Create E vectors 33377841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_ceed, NULL); 33477841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_dot_ceed, 33577841947SLeila Ghaffari NULL); 33677841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->g_ceed, NULL); 33777841947SLeila Ghaffari 33877841947SLeila Ghaffari // ----------------------------------------------------------------------------- 33977841947SLeila Ghaffari // CEED QFunctions 34077841947SLeila Ghaffari // ----------------------------------------------------------------------------- 34177841947SLeila Ghaffari // -- Create QFunction for quadrature data 34291e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->setup_vol.qfunction, 34391e5af17SJed Brown problem->setup_vol.qfunction_loc, 34477841947SLeila Ghaffari &ceed_data->qf_setup_vol); 345841e4c73SJed Brown if (problem->setup_vol.qfunction_context) { 346841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_setup_vol, 347841e4c73SJed Brown problem->setup_vol.qfunction_context); 348841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->setup_vol.qfunction_context); 349841e4c73SJed Brown } 35077841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_vol, "dx", num_comp_x*dim, 35177841947SLeila Ghaffari CEED_EVAL_GRAD); 35277841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_vol, "weight", 1, CEED_EVAL_WEIGHT); 353a61c78d6SJeremy L Thompson CeedQFunctionAddOutput(ceed_data->qf_setup_vol, "qdata", q_data_size_vol, 35477841947SLeila Ghaffari CEED_EVAL_NONE); 35577841947SLeila Ghaffari 35677841947SLeila Ghaffari // -- Create QFunction for ICs 35791e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->ics.qfunction, 35891e5af17SJed Brown problem->ics.qfunction_loc, 35977841947SLeila Ghaffari &ceed_data->qf_ics); 360841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_ics, problem->ics.qfunction_context); 361841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->ics.qfunction_context); 36277841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ics, "x", num_comp_x, CEED_EVAL_INTERP); 36377841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_ics, "q0", num_comp_q, CEED_EVAL_NONE); 36477841947SLeila Ghaffari 36577841947SLeila Ghaffari // -- Create QFunction for RHS 36691e5af17SJed Brown if (problem->apply_vol_rhs.qfunction) { 36791e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_rhs.qfunction, 36891e5af17SJed Brown problem->apply_vol_rhs.qfunction_loc, &ceed_data->qf_rhs_vol); 369841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_rhs_vol, 370841e4c73SJed Brown problem->apply_vol_rhs.qfunction_context); 371841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_rhs.qfunction_context); 37277841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "q", num_comp_q, CEED_EVAL_INTERP); 373a3ae0734SJed Brown CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "Grad_q", num_comp_q*dim, 37477841947SLeila Ghaffari CEED_EVAL_GRAD); 375a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "qdata", q_data_size_vol, 37677841947SLeila Ghaffari CEED_EVAL_NONE); 37777841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "x", num_comp_x, CEED_EVAL_INTERP); 37877841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "v", num_comp_q, 37977841947SLeila Ghaffari CEED_EVAL_INTERP); 380a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "Grad_v", num_comp_q*dim, 38177841947SLeila Ghaffari CEED_EVAL_GRAD); 38277841947SLeila Ghaffari } 38377841947SLeila Ghaffari 38477841947SLeila Ghaffari // -- Create QFunction for IFunction 38591e5af17SJed Brown if (problem->apply_vol_ifunction.qfunction) { 38691e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ifunction.qfunction, 38791e5af17SJed Brown problem->apply_vol_ifunction.qfunction_loc, &ceed_data->qf_ifunction_vol); 388841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_ifunction_vol, 389841e4c73SJed Brown problem->apply_vol_ifunction.qfunction_context); 390841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_ifunction.qfunction_context); 39177841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q", num_comp_q, 39277841947SLeila Ghaffari CEED_EVAL_INTERP); 393a3ae0734SJed Brown CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "Grad_q", num_comp_q*dim, 39477841947SLeila Ghaffari CEED_EVAL_GRAD); 395a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q dot", num_comp_q, 39677841947SLeila Ghaffari CEED_EVAL_INTERP); 397a61c78d6SJeremy L Thompson CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "qdata", q_data_size_vol, 39877841947SLeila Ghaffari CEED_EVAL_NONE); 39977841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "x", num_comp_x, 40077841947SLeila Ghaffari CEED_EVAL_INTERP); 40177841947SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "v", num_comp_q, 40277841947SLeila Ghaffari CEED_EVAL_INTERP); 403a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "Grad_v", num_comp_q*dim, 40477841947SLeila Ghaffari CEED_EVAL_GRAD); 405a3ae0734SJed Brown CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "jac_data", 406a3ae0734SJed Brown jac_data_size_vol, CEED_EVAL_NONE); 40777841947SLeila Ghaffari } 40877841947SLeila Ghaffari 409e334ad8fSJed Brown CeedQFunction qf_ijacobian_vol = NULL; 410e334ad8fSJed Brown if (problem->apply_vol_ijacobian.qfunction) { 411e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ijacobian.qfunction, 412e334ad8fSJed Brown problem->apply_vol_ijacobian.qfunction_loc, &qf_ijacobian_vol); 413e334ad8fSJed Brown CeedQFunctionSetContext(qf_ijacobian_vol, 414e334ad8fSJed Brown problem->apply_vol_ijacobian.qfunction_context); 415e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_vol_ijacobian.qfunction_context); 416e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "dq", num_comp_q, 417e334ad8fSJed Brown CEED_EVAL_INTERP); 418e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "Grad_dq", num_comp_q*dim, 419e334ad8fSJed Brown CEED_EVAL_GRAD); 420e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "qdata", q_data_size_vol, 421e334ad8fSJed Brown CEED_EVAL_NONE); 422e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "x", num_comp_x, 423e334ad8fSJed Brown CEED_EVAL_INTERP); 424e334ad8fSJed Brown CeedQFunctionAddInput(qf_ijacobian_vol, "jac_data", 425e334ad8fSJed Brown jac_data_size_vol, CEED_EVAL_NONE); 426e334ad8fSJed Brown CeedQFunctionAddOutput(qf_ijacobian_vol, "v", num_comp_q, 427e334ad8fSJed Brown CEED_EVAL_INTERP); 428e334ad8fSJed Brown CeedQFunctionAddOutput(qf_ijacobian_vol, "Grad_v", num_comp_q*dim, 429e334ad8fSJed Brown CEED_EVAL_GRAD); 430e334ad8fSJed Brown } 431e334ad8fSJed Brown 43277841947SLeila Ghaffari // --------------------------------------------------------------------------- 43377841947SLeila Ghaffari // Element coordinates 43477841947SLeila Ghaffari // --------------------------------------------------------------------------- 43577841947SLeila Ghaffari // -- Create CEED vector 43677841947SLeila Ghaffari CeedElemRestrictionCreateVector(ceed_data->elem_restr_x, &ceed_data->x_coord, 43777841947SLeila Ghaffari NULL); 43877841947SLeila Ghaffari 43977841947SLeila Ghaffari // -- Copy PETSc vector in CEED vector 44077841947SLeila Ghaffari Vec X_loc; 44177841947SLeila Ghaffari const PetscScalar *X_loc_array; 44277841947SLeila Ghaffari ierr = DMGetCoordinatesLocal(dm, &X_loc); CHKERRQ(ierr); 4431864f1c2SLeila Ghaffari ierr = VecScale(X_loc, problem->dm_scale); CHKERRQ(ierr); 44477841947SLeila Ghaffari ierr = VecGetArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 44577841947SLeila Ghaffari CeedVectorSetArray(ceed_data->x_coord, CEED_MEM_HOST, CEED_COPY_VALUES, 44677841947SLeila Ghaffari (PetscScalar *)X_loc_array); 44777841947SLeila Ghaffari ierr = VecRestoreArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 44877841947SLeila Ghaffari 44977841947SLeila Ghaffari // ----------------------------------------------------------------------------- 45077841947SLeila Ghaffari // CEED vectors 45177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 45277841947SLeila Ghaffari // -- Create CEED vector for geometric data 45377841947SLeila Ghaffari CeedInt num_qpts_vol; 45477841947SLeila Ghaffari PetscInt loc_num_elem_vol; 45577841947SLeila Ghaffari CeedBasisGetNumQuadraturePoints(ceed_data->basis_q, &num_qpts_vol); 45677841947SLeila Ghaffari CeedElemRestrictionGetNumElements(ceed_data->elem_restr_q, &loc_num_elem_vol); 45777841947SLeila Ghaffari CeedVectorCreate(ceed, q_data_size_vol*loc_num_elem_vol*num_qpts_vol, 45877841947SLeila Ghaffari &ceed_data->q_data); 45977841947SLeila Ghaffari 460a3ae0734SJed Brown CeedElemRestrictionCreateVector(elem_restr_jd_i, &jac_data, NULL); 46177841947SLeila Ghaffari // ----------------------------------------------------------------------------- 46277841947SLeila Ghaffari // CEED Operators 46377841947SLeila Ghaffari // ----------------------------------------------------------------------------- 46477841947SLeila Ghaffari // -- Create CEED operator for quadrature data 46577841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_setup_vol, NULL, NULL, 46677841947SLeila Ghaffari &ceed_data->op_setup_vol); 46777841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_setup_vol, "dx", ceed_data->elem_restr_x, 46877841947SLeila Ghaffari ceed_data->basis_x, CEED_VECTOR_ACTIVE); 46977841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_setup_vol, "weight", 470e6225c47SLeila Ghaffari CEED_ELEMRESTRICTION_NONE, ceed_data->basis_x, CEED_VECTOR_NONE); 471a61c78d6SJeremy L Thompson CeedOperatorSetField(ceed_data->op_setup_vol, "qdata", 472e6225c47SLeila Ghaffari ceed_data->elem_restr_qd_i, CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 47377841947SLeila Ghaffari 47477841947SLeila Ghaffari // -- Create CEED operator for ICs 47577841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_ics, NULL, NULL, &ceed_data->op_ics); 47677841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_ics, "x", ceed_data->elem_restr_x, 47777841947SLeila Ghaffari ceed_data->basis_xc, CEED_VECTOR_ACTIVE); 47877841947SLeila Ghaffari CeedOperatorSetField(ceed_data->op_ics, "q0", ceed_data->elem_restr_q, 47977841947SLeila Ghaffari CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 480841e4c73SJed Brown CeedOperatorContextGetFieldLabel(ceed_data->op_ics, "evaluation time", 481841e4c73SJed Brown &user->phys->ics_time_label); 48277841947SLeila Ghaffari 48377841947SLeila Ghaffari // Create CEED operator for RHS 48477841947SLeila Ghaffari if (ceed_data->qf_rhs_vol) { 48577841947SLeila Ghaffari CeedOperator op; 48677841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_rhs_vol, NULL, NULL, &op); 48777841947SLeila Ghaffari CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 48877841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 489a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_q", ceed_data->elem_restr_q, ceed_data->basis_q, 49077841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 491a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 492e6225c47SLeila Ghaffari CEED_BASIS_COLLOCATED, ceed_data->q_data); 49377841947SLeila Ghaffari CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 49477841947SLeila Ghaffari ceed_data->x_coord); 49577841947SLeila Ghaffari CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 49677841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 497a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 49877841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 49977841947SLeila Ghaffari user->op_rhs_vol = op; 50077841947SLeila Ghaffari } 50177841947SLeila Ghaffari 50277841947SLeila Ghaffari // -- CEED operator for IFunction 50377841947SLeila Ghaffari if (ceed_data->qf_ifunction_vol) { 50477841947SLeila Ghaffari CeedOperator op; 50577841947SLeila Ghaffari CeedOperatorCreate(ceed, ceed_data->qf_ifunction_vol, NULL, NULL, &op); 50677841947SLeila Ghaffari CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 50777841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 508a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_q", ceed_data->elem_restr_q, ceed_data->basis_q, 50977841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 510a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "q dot", ceed_data->elem_restr_q, ceed_data->basis_q, 51177841947SLeila Ghaffari user->q_dot_ceed); 512a61c78d6SJeremy L Thompson CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 513e6225c47SLeila Ghaffari CEED_BASIS_COLLOCATED, ceed_data->q_data); 51477841947SLeila Ghaffari CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 51577841947SLeila Ghaffari ceed_data->x_coord); 51677841947SLeila Ghaffari CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 51777841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 518a3ae0734SJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 51977841947SLeila Ghaffari CEED_VECTOR_ACTIVE); 520a3ae0734SJed Brown CeedOperatorSetField(op, "jac_data", elem_restr_jd_i, 521a3ae0734SJed Brown CEED_BASIS_COLLOCATED, jac_data); 522a3ae0734SJed Brown 52377841947SLeila Ghaffari user->op_ifunction_vol = op; 52477841947SLeila Ghaffari } 52577841947SLeila Ghaffari 526e334ad8fSJed Brown CeedOperator op_ijacobian_vol = NULL; 527e334ad8fSJed Brown if (qf_ijacobian_vol) { 528e334ad8fSJed Brown CeedOperator op; 529e334ad8fSJed Brown CeedOperatorCreate(ceed, qf_ijacobian_vol, NULL, NULL, &op); 530e334ad8fSJed Brown CeedOperatorSetField(op, "dq", ceed_data->elem_restr_q, ceed_data->basis_q, 531e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 532e334ad8fSJed Brown CeedOperatorSetField(op, "Grad_dq", ceed_data->elem_restr_q, ceed_data->basis_q, 533e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 534e334ad8fSJed Brown CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 535e334ad8fSJed Brown CEED_BASIS_COLLOCATED, ceed_data->q_data); 536e334ad8fSJed Brown CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 537e334ad8fSJed Brown ceed_data->x_coord); 538e334ad8fSJed Brown CeedOperatorSetField(op, "jac_data", elem_restr_jd_i, 539e334ad8fSJed Brown CEED_BASIS_COLLOCATED, jac_data); 540e334ad8fSJed Brown CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 541e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 542e334ad8fSJed Brown CeedOperatorSetField(op, "Grad_v", ceed_data->elem_restr_q, ceed_data->basis_q, 543e334ad8fSJed Brown CEED_VECTOR_ACTIVE); 544e334ad8fSJed Brown op_ijacobian_vol = op; 545e334ad8fSJed Brown CeedQFunctionDestroy(&qf_ijacobian_vol); 546e334ad8fSJed Brown } 547e334ad8fSJed Brown 54877841947SLeila Ghaffari // ***************************************************************************** 54977841947SLeila Ghaffari // Set up CEED objects for the exterior domain (surface) 55077841947SLeila Ghaffari // ***************************************************************************** 55177841947SLeila Ghaffari CeedInt height = 1, 55277841947SLeila Ghaffari dim_sur = dim - height, 55377841947SLeila Ghaffari P_sur = app_ctx->degree + 1, 55477841947SLeila Ghaffari Q_sur = P_sur + app_ctx->q_extra; 555e334ad8fSJed Brown const CeedInt q_data_size_sur = problem->q_data_size_sur, 556e334ad8fSJed Brown jac_data_size_sur = problem->jac_data_size_sur; 55777841947SLeila Ghaffari 55877841947SLeila Ghaffari // ----------------------------------------------------------------------------- 55977841947SLeila Ghaffari // CEED Bases 56077841947SLeila Ghaffari // ----------------------------------------------------------------------------- 56177841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_q, P_sur, Q_sur, 56277841947SLeila Ghaffari CEED_GAUSS, &ceed_data->basis_q_sur); 56377841947SLeila Ghaffari CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_x, 2, Q_sur, CEED_GAUSS, 56477841947SLeila Ghaffari &ceed_data->basis_x_sur); 56577841947SLeila Ghaffari 56677841947SLeila Ghaffari // ----------------------------------------------------------------------------- 56777841947SLeila Ghaffari // CEED QFunctions 56877841947SLeila Ghaffari // ----------------------------------------------------------------------------- 56977841947SLeila Ghaffari // -- Create QFunction for quadrature data 57091e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->setup_sur.qfunction, 57191e5af17SJed Brown problem->setup_sur.qfunction_loc, 57277841947SLeila Ghaffari &ceed_data->qf_setup_sur); 573841e4c73SJed Brown if (problem->setup_sur.qfunction_context) { 574841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_setup_sur, 575841e4c73SJed Brown problem->setup_sur.qfunction_context); 576841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->setup_sur.qfunction_context); 577841e4c73SJed Brown } 57877841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_sur, "dx", num_comp_x*dim_sur, 57977841947SLeila Ghaffari CEED_EVAL_GRAD); 58077841947SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_setup_sur, "weight", 1, CEED_EVAL_WEIGHT); 581a61c78d6SJeremy L Thompson CeedQFunctionAddOutput(ceed_data->qf_setup_sur, "surface qdata", 5822fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 58377841947SLeila Ghaffari 5842fe7aee7SLeila Ghaffari // -- Creat QFunction for inflow boundaries 58591e5af17SJed Brown if (problem->apply_inflow.qfunction) { 58691e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_inflow.qfunction, 58791e5af17SJed Brown problem->apply_inflow.qfunction_loc, &ceed_data->qf_apply_inflow); 588841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_inflow, 589841e4c73SJed Brown problem->apply_inflow.qfunction_context); 590841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_inflow.qfunction_context); 5912fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "q", num_comp_q, 59277841947SLeila Ghaffari CEED_EVAL_INTERP); 593e8b03feeSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "Grad_q", num_comp_q*(dim-1), 594e8b03feeSJames Wright CEED_EVAL_GRAD); 5952fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "surface qdata", 5962fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 5972fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "x", num_comp_x, 59877841947SLeila Ghaffari CEED_EVAL_INTERP); 5992fe7aee7SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_apply_inflow, "v", num_comp_q, 6002fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 6012fe7aee7SLeila Ghaffari } 602e334ad8fSJed Brown if (problem->apply_inflow_jacobian.qfunction) { 603e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_inflow_jacobian.qfunction, 604e334ad8fSJed Brown problem->apply_inflow_jacobian.qfunction_loc, 605e334ad8fSJed Brown &ceed_data->qf_apply_inflow_jacobian); 606e334ad8fSJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_inflow_jacobian, 607e334ad8fSJed Brown problem->apply_inflow_jacobian.qfunction_context); 608e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_inflow_jacobian.qfunction_context); 609e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "dq", num_comp_q, 610e334ad8fSJed Brown CEED_EVAL_INTERP); 611e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "surface qdata", 612e334ad8fSJed Brown q_data_size_sur, CEED_EVAL_NONE); 613e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_inflow_jacobian, "x", num_comp_x, 614e334ad8fSJed Brown CEED_EVAL_INTERP); 615e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_inflow_jacobian, "v", num_comp_q, 616e334ad8fSJed Brown CEED_EVAL_INTERP); 617e334ad8fSJed Brown } 6182fe7aee7SLeila Ghaffari 6192fe7aee7SLeila Ghaffari // -- Creat QFunction for outflow boundaries 62091e5af17SJed Brown if (problem->apply_outflow.qfunction) { 62191e5af17SJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_outflow.qfunction, 62291e5af17SJed Brown problem->apply_outflow.qfunction_loc, &ceed_data->qf_apply_outflow); 623841e4c73SJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_outflow, 624841e4c73SJed Brown problem->apply_outflow.qfunction_context); 625841e4c73SJed Brown CeedQFunctionContextDestroy(&problem->apply_outflow.qfunction_context); 6262fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "q", num_comp_q, 6272fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 628e8b03feeSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "Grad_q", num_comp_q*(dim-1), 629e8b03feeSJames Wright CEED_EVAL_GRAD); 6302fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "surface qdata", 6312fe7aee7SLeila Ghaffari q_data_size_sur, CEED_EVAL_NONE); 6322fe7aee7SLeila Ghaffari CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "x", num_comp_x, 6332fe7aee7SLeila Ghaffari CEED_EVAL_INTERP); 6342fe7aee7SLeila Ghaffari CeedQFunctionAddOutput(ceed_data->qf_apply_outflow, "v", num_comp_q, 63577841947SLeila Ghaffari CEED_EVAL_INTERP); 63639c69132SJed Brown if (jac_data_size_sur) 637e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_outflow, "surface jacobian data", 638e334ad8fSJed Brown jac_data_size_sur, 639e334ad8fSJed Brown CEED_EVAL_NONE); 640e334ad8fSJed Brown } 641e334ad8fSJed Brown if (problem->apply_outflow_jacobian.qfunction) { 642e334ad8fSJed Brown CeedQFunctionCreateInterior(ceed, 1, problem->apply_outflow_jacobian.qfunction, 643e334ad8fSJed Brown problem->apply_outflow_jacobian.qfunction_loc, 644e334ad8fSJed Brown &ceed_data->qf_apply_outflow_jacobian); 645e334ad8fSJed Brown CeedQFunctionSetContext(ceed_data->qf_apply_outflow_jacobian, 646e334ad8fSJed Brown problem->apply_outflow_jacobian.qfunction_context); 647e334ad8fSJed Brown CeedQFunctionContextDestroy(&problem->apply_outflow_jacobian.qfunction_context); 648e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "dq", num_comp_q, 649e334ad8fSJed Brown CEED_EVAL_INTERP); 650*0ec2498eSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "Grad_dq", num_comp_q*dim_sur, 651*0ec2498eSJames Wright CEED_EVAL_GRAD); 652e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "surface qdata", 653e334ad8fSJed Brown q_data_size_sur, CEED_EVAL_NONE); 654*0ec2498eSJames Wright CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, "x", num_comp_x, 655*0ec2498eSJames Wright CEED_EVAL_INTERP); 656e334ad8fSJed Brown CeedQFunctionAddInput(ceed_data->qf_apply_outflow_jacobian, 657e334ad8fSJed Brown "surface jacobian data", 658e334ad8fSJed Brown jac_data_size_sur, CEED_EVAL_NONE); 659e334ad8fSJed Brown CeedQFunctionAddOutput(ceed_data->qf_apply_outflow_jacobian, "v", num_comp_q, 660e334ad8fSJed Brown CEED_EVAL_INTERP); 66177841947SLeila Ghaffari } 66277841947SLeila Ghaffari 66377841947SLeila Ghaffari // ***************************************************************************** 66477841947SLeila Ghaffari // CEED Operator Apply 66577841947SLeila Ghaffari // ***************************************************************************** 66677841947SLeila Ghaffari // -- Apply CEED Operator for the geometric data 66777841947SLeila Ghaffari CeedOperatorApply(ceed_data->op_setup_vol, ceed_data->x_coord, 66877841947SLeila Ghaffari ceed_data->q_data, CEED_REQUEST_IMMEDIATE); 66977841947SLeila Ghaffari 67077841947SLeila Ghaffari // -- Create and apply CEED Composite Operator for the entire domain 67177841947SLeila Ghaffari if (!user->phys->implicit) { // RHS 67277841947SLeila Ghaffari ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 673e334ad8fSJed Brown user->op_rhs_vol, NULL, height, P_sur, Q_sur, 674e334ad8fSJed Brown q_data_size_sur, 0, 675e334ad8fSJed Brown &user->op_rhs, NULL); CHKERRQ(ierr); 67677841947SLeila Ghaffari } else { // IFunction 67777841947SLeila Ghaffari ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 678e334ad8fSJed Brown user->op_ifunction_vol, op_ijacobian_vol, 679e334ad8fSJed Brown height, P_sur, Q_sur, 680e334ad8fSJed Brown q_data_size_sur, jac_data_size_sur, 681e334ad8fSJed Brown &user->op_ifunction, 682e334ad8fSJed Brown op_ijacobian_vol ? &user->op_ijacobian : NULL); CHKERRQ(ierr); 683e334ad8fSJed Brown if (user->op_ijacobian) { 684e334ad8fSJed Brown CeedOperatorContextGetFieldLabel(user->op_ijacobian, "ijacobian time shift", 685e334ad8fSJed Brown &user->phys->ijacobian_time_shift_label); 686e334ad8fSJed Brown } 687e334ad8fSJed Brown 68877841947SLeila Ghaffari } 689a3ae0734SJed Brown CeedElemRestrictionDestroy(&elem_restr_jd_i); 69022d320f5SLeila Ghaffari CeedOperatorDestroy(&op_ijacobian_vol); 691a3ae0734SJed Brown CeedVectorDestroy(&jac_data); 69277841947SLeila Ghaffari PetscFunctionReturn(0); 69377841947SLeila Ghaffari } 694