1 // Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2 // the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3 // reserved. See files LICENSE and NOTICE for details. 4 // 5 // This file is part of CEED, a collection of benchmarks, miniapps, software 6 // libraries and APIs for efficient high-order finite element and spectral 7 // element discretizations for exascale applications. For more information and 8 // source code availability see http://github.com/ceed. 9 // 10 // The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11 // a collaborative effort of two U.S. Department of Energy organizations (Office 12 // of Science and the National Nuclear Security Administration) responsible for 13 // the planning and preparation of a capable exascale ecosystem, including 14 // software, applications, hardware, advanced system engineering and early 15 // testbed platforms, in support of the nation's exascale computing imperative. 16 17 /// @file 18 /// Setup libCEED for Navier-Stokes example using PETSc 19 20 #include "../navierstokes.h" 21 22 // Utility function - essential BC dofs are encoded in closure indices as -(i+1). 23 PetscInt Involute(PetscInt i) { 24 return i >= 0 ? i : -(i+1); 25 } 26 27 // Utility function to create local CEED restriction 28 PetscErrorCode CreateRestrictionFromPlex(Ceed ceed, DM dm, CeedInt height, 29 DMLabel domain_label, CeedInt value, CeedElemRestriction *elem_restr) { 30 PetscInt num_elem, elem_size, num_dof, num_comp, *elem_restr_offsets; 31 PetscErrorCode ierr; 32 33 PetscFunctionBeginUser; 34 ierr = DMPlexGetLocalOffsets(dm, domain_label, value, height, 0, &num_elem, 35 &elem_size, &num_comp, &num_dof, &elem_restr_offsets); 36 CHKERRQ(ierr); 37 38 CeedElemRestrictionCreate(ceed, num_elem, elem_size, num_comp, 39 1, num_dof, CEED_MEM_HOST, CEED_COPY_VALUES, 40 elem_restr_offsets, elem_restr); 41 ierr = PetscFree(elem_restr_offsets); CHKERRQ(ierr); 42 43 PetscFunctionReturn(0); 44 } 45 46 // Utility function to get Ceed Restriction for each domain 47 PetscErrorCode GetRestrictionForDomain(Ceed ceed, DM dm, CeedInt height, 48 DMLabel domain_label, PetscInt value, 49 CeedInt Q, CeedInt q_data_size, 50 CeedElemRestriction *elem_restr_q, 51 CeedElemRestriction *elem_restr_x, 52 CeedElemRestriction *elem_restr_qd_i) { 53 DM dm_coord; 54 CeedInt dim, loc_num_elem; 55 CeedInt Q_dim; 56 PetscErrorCode ierr; 57 PetscFunctionBeginUser; 58 59 ierr = DMGetDimension(dm, &dim); CHKERRQ(ierr); 60 dim -= height; 61 Q_dim = CeedIntPow(Q, dim); 62 ierr = DMGetCoordinateDM(dm, &dm_coord); CHKERRQ(ierr); 63 ierr = DMPlexSetClosurePermutationTensor(dm_coord, PETSC_DETERMINE, NULL); 64 CHKERRQ(ierr); 65 ierr = CreateRestrictionFromPlex(ceed, dm, height, domain_label, value, 66 elem_restr_q); 67 CHKERRQ(ierr); 68 ierr = CreateRestrictionFromPlex(ceed, dm_coord, height, domain_label, value, 69 elem_restr_x); 70 CHKERRQ(ierr); 71 CeedElemRestrictionGetNumElements(*elem_restr_q, &loc_num_elem); 72 CeedElemRestrictionCreateStrided(ceed, loc_num_elem, Q_dim, 73 q_data_size, q_data_size*loc_num_elem*Q_dim, 74 CEED_STRIDES_BACKEND, elem_restr_qd_i); 75 PetscFunctionReturn(0); 76 } 77 78 // Utility function to create CEED Composite Operator for the entire domain 79 PetscErrorCode CreateOperatorForDomain(Ceed ceed, DM dm, SimpleBC bc, 80 CeedData ceed_data, Physics phys, 81 CeedOperator op_apply_vol, CeedInt height, 82 CeedInt P_sur, CeedInt Q_sur, CeedInt q_data_size_sur, 83 CeedOperator *op_apply) { 84 //CeedInt dim; 85 DMLabel domain_label; 86 PetscErrorCode ierr; 87 PetscFunctionBeginUser; 88 89 // Create Composite Operaters 90 CeedCompositeOperatorCreate(ceed, op_apply); 91 92 // --Apply Sub-Operator for the volume 93 CeedCompositeOperatorAddSub(*op_apply, op_apply_vol); 94 95 // -- Create Sub-Operator for in/outflow BCs 96 if (phys->has_neumann) { 97 // --- Setup 98 ierr = DMGetLabel(dm, "Face Sets", &domain_label); CHKERRQ(ierr); 99 //ierr = DMGetDimension(dm, &dim); CHKERRQ(ierr); 100 101 // --- Get number of quadrature points for the boundaries 102 CeedInt num_qpts_sur; 103 CeedBasisGetNumQuadraturePoints(ceed_data->basis_q_sur, &num_qpts_sur); 104 105 // --- Create Sub-Operator for inflow boundaries 106 for (CeedInt i=0; i < bc->num_inflow; i++) { 107 CeedVector q_data_sur; 108 CeedOperator op_setup_sur, op_apply_inflow; 109 CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur; 110 111 // ---- CEED Restriction 112 ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->inflows[i], 113 Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 114 &elem_restr_qd_i_sur); 115 CHKERRQ(ierr); 116 117 // ---- CEED Vector 118 PetscInt loc_num_elem_sur; 119 CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 120 CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 121 &q_data_sur); 122 123 // ---- CEED Operator 124 // ----- CEED Operator for Setup (geometric factors) 125 CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 126 CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 127 ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 128 CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 129 ceed_data->basis_x_sur, CEED_VECTOR_NONE); 130 CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 131 CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 132 133 // ----- CEED Operator for Physics 134 CeedOperatorCreate(ceed, ceed_data->qf_apply_inflow, NULL, NULL, 135 &op_apply_inflow); 136 CeedOperatorSetField(op_apply_inflow, "q", elem_restr_q_sur, 137 ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 138 CeedOperatorSetField(op_apply_inflow, "surface qdata", elem_restr_qd_i_sur, 139 CEED_BASIS_COLLOCATED, q_data_sur); 140 CeedOperatorSetField(op_apply_inflow, "x", elem_restr_x_sur, 141 ceed_data->basis_x_sur, ceed_data->x_coord); 142 CeedOperatorSetField(op_apply_inflow, "v", elem_restr_q_sur, 143 ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 144 145 // ----- Apply CEED operator for Setup 146 CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 147 CEED_REQUEST_IMMEDIATE); 148 149 // ----- Apply Sub-Operator for Physics 150 CeedCompositeOperatorAddSub(*op_apply, op_apply_inflow); 151 152 // ----- Cleanup 153 CeedVectorDestroy(&q_data_sur); 154 CeedElemRestrictionDestroy(&elem_restr_q_sur); 155 CeedElemRestrictionDestroy(&elem_restr_x_sur); 156 CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 157 CeedOperatorDestroy(&op_setup_sur); 158 CeedOperatorDestroy(&op_apply_inflow); 159 } 160 161 // --- Create Sub-Operator for outflow boundaries 162 for (CeedInt i=0; i < bc->num_outflow; i++) { 163 CeedVector q_data_sur; 164 CeedOperator op_setup_sur, op_apply_outflow; 165 CeedElemRestriction elem_restr_x_sur, elem_restr_q_sur, elem_restr_qd_i_sur; 166 167 // ---- CEED Restriction 168 ierr = GetRestrictionForDomain(ceed, dm, height, domain_label, bc->outflows[i], 169 Q_sur, q_data_size_sur, &elem_restr_q_sur, &elem_restr_x_sur, 170 &elem_restr_qd_i_sur); 171 CHKERRQ(ierr); 172 173 // ---- CEED Vector 174 PetscInt loc_num_elem_sur; 175 CeedElemRestrictionGetNumElements(elem_restr_q_sur, &loc_num_elem_sur); 176 CeedVectorCreate(ceed, q_data_size_sur*loc_num_elem_sur*num_qpts_sur, 177 &q_data_sur); 178 179 // ---- CEED Operator 180 // ----- CEED Operator for Setup (geometric factors) 181 CeedOperatorCreate(ceed, ceed_data->qf_setup_sur, NULL, NULL, &op_setup_sur); 182 CeedOperatorSetField(op_setup_sur, "dx", elem_restr_x_sur, 183 ceed_data->basis_x_sur, CEED_VECTOR_ACTIVE); 184 CeedOperatorSetField(op_setup_sur, "weight", CEED_ELEMRESTRICTION_NONE, 185 ceed_data->basis_x_sur, CEED_VECTOR_NONE); 186 CeedOperatorSetField(op_setup_sur, "surface qdata", elem_restr_qd_i_sur, 187 CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 188 189 // ----- CEED Operator for Physics 190 CeedOperatorCreate(ceed, ceed_data->qf_apply_outflow, NULL, NULL, 191 &op_apply_outflow); 192 CeedOperatorSetField(op_apply_outflow, "q", elem_restr_q_sur, 193 ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 194 CeedOperatorSetField(op_apply_outflow, "surface qdata", elem_restr_qd_i_sur, 195 CEED_BASIS_COLLOCATED, q_data_sur); 196 CeedOperatorSetField(op_apply_outflow, "x", elem_restr_x_sur, 197 ceed_data->basis_x_sur, ceed_data->x_coord); 198 CeedOperatorSetField(op_apply_outflow, "v", elem_restr_q_sur, 199 ceed_data->basis_q_sur, CEED_VECTOR_ACTIVE); 200 201 // ----- Apply CEED operator for Setup 202 CeedOperatorApply(op_setup_sur, ceed_data->x_coord, q_data_sur, 203 CEED_REQUEST_IMMEDIATE); 204 205 // ----- Apply Sub-Operator for Physics 206 CeedCompositeOperatorAddSub(*op_apply, op_apply_outflow); 207 208 // ----- Cleanup 209 CeedVectorDestroy(&q_data_sur); 210 CeedElemRestrictionDestroy(&elem_restr_q_sur); 211 CeedElemRestrictionDestroy(&elem_restr_x_sur); 212 CeedElemRestrictionDestroy(&elem_restr_qd_i_sur); 213 CeedOperatorDestroy(&op_setup_sur); 214 CeedOperatorDestroy(&op_apply_outflow); 215 } 216 } 217 PetscFunctionReturn(0); 218 } 219 220 PetscErrorCode SetupLibceed(Ceed ceed, CeedData ceed_data, DM dm, User user, 221 AppCtx app_ctx, ProblemData *problem, SimpleBC bc) { 222 PetscErrorCode ierr; 223 PetscFunctionBeginUser; 224 225 // ***************************************************************************** 226 // Set up CEED objects for the interior domain (volume) 227 // ***************************************************************************** 228 const PetscInt num_comp_q = 5; 229 const CeedInt dim = problem->dim, 230 num_comp_x = problem->dim, 231 q_data_size_vol = problem->q_data_size_vol, 232 P = app_ctx->degree + 1, 233 Q = P + app_ctx->q_extra; 234 235 // ----------------------------------------------------------------------------- 236 // CEED Bases 237 // ----------------------------------------------------------------------------- 238 CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_q, P, Q, CEED_GAUSS, 239 &ceed_data->basis_q); 240 CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, Q, CEED_GAUSS, 241 &ceed_data->basis_x); 242 CeedBasisCreateTensorH1Lagrange(ceed, dim, num_comp_x, 2, P, 243 CEED_GAUSS_LOBATTO, &ceed_data->basis_xc); 244 245 // ----------------------------------------------------------------------------- 246 // CEED Restrictions 247 // ----------------------------------------------------------------------------- 248 // -- Create restriction 249 ierr = GetRestrictionForDomain(ceed, dm, 0, 0, 0, Q, q_data_size_vol, 250 &ceed_data->elem_restr_q, &ceed_data->elem_restr_x, 251 &ceed_data->elem_restr_qd_i); CHKERRQ(ierr); 252 // -- Create E vectors 253 CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_ceed, NULL); 254 CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_dot_ceed, 255 NULL); 256 CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->g_ceed, NULL); 257 258 // ----------------------------------------------------------------------------- 259 // CEED QFunctions 260 // ----------------------------------------------------------------------------- 261 // -- Create QFunction for quadrature data 262 CeedQFunctionCreateInterior(ceed, 1, problem->setup_vol, problem->setup_vol_loc, 263 &ceed_data->qf_setup_vol); 264 CeedQFunctionAddInput(ceed_data->qf_setup_vol, "dx", num_comp_x*dim, 265 CEED_EVAL_GRAD); 266 CeedQFunctionAddInput(ceed_data->qf_setup_vol, "weight", 1, CEED_EVAL_WEIGHT); 267 CeedQFunctionAddOutput(ceed_data->qf_setup_vol, "qdata", q_data_size_vol, 268 CEED_EVAL_NONE); 269 270 // -- Create QFunction for ICs 271 CeedQFunctionCreateInterior(ceed, 1, problem->ics, problem->ics_loc, 272 &ceed_data->qf_ics); 273 CeedQFunctionAddInput(ceed_data->qf_ics, "x", num_comp_x, CEED_EVAL_INTERP); 274 CeedQFunctionAddOutput(ceed_data->qf_ics, "q0", num_comp_q, CEED_EVAL_NONE); 275 276 // -- Create QFunction for RHS 277 if (problem->apply_vol_rhs) { 278 CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_rhs, 279 problem->apply_vol_rhs_loc, &ceed_data->qf_rhs_vol); 280 CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "q", num_comp_q, CEED_EVAL_INTERP); 281 CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "dq", num_comp_q*dim, 282 CEED_EVAL_GRAD); 283 CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "qdata", q_data_size_vol, 284 CEED_EVAL_NONE); 285 CeedQFunctionAddInput(ceed_data->qf_rhs_vol, "x", num_comp_x, CEED_EVAL_INTERP); 286 CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "v", num_comp_q, 287 CEED_EVAL_INTERP); 288 CeedQFunctionAddOutput(ceed_data->qf_rhs_vol, "dv", num_comp_q*dim, 289 CEED_EVAL_GRAD); 290 } 291 292 // -- Create QFunction for IFunction 293 if (problem->apply_vol_ifunction) { 294 CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ifunction, 295 problem->apply_vol_ifunction_loc, &ceed_data->qf_ifunction_vol); 296 CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q", num_comp_q, 297 CEED_EVAL_INTERP); 298 CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "dq", num_comp_q*dim, 299 CEED_EVAL_GRAD); 300 CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "q dot", num_comp_q, 301 CEED_EVAL_INTERP); 302 CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "qdata", q_data_size_vol, 303 CEED_EVAL_NONE); 304 CeedQFunctionAddInput(ceed_data->qf_ifunction_vol, "x", num_comp_x, 305 CEED_EVAL_INTERP); 306 CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "v", num_comp_q, 307 CEED_EVAL_INTERP); 308 CeedQFunctionAddOutput(ceed_data->qf_ifunction_vol, "dv", num_comp_q*dim, 309 CEED_EVAL_GRAD); 310 } 311 312 // --------------------------------------------------------------------------- 313 // Element coordinates 314 // --------------------------------------------------------------------------- 315 // -- Create CEED vector 316 CeedElemRestrictionCreateVector(ceed_data->elem_restr_x, &ceed_data->x_coord, 317 NULL); 318 319 // -- Copy PETSc vector in CEED vector 320 Vec X_loc; 321 const PetscScalar *X_loc_array; 322 ierr = DMGetCoordinatesLocal(dm, &X_loc); CHKERRQ(ierr); 323 ierr = VecScale(X_loc, problem->dm_scale); CHKERRQ(ierr); 324 ierr = VecGetArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 325 CeedVectorSetArray(ceed_data->x_coord, CEED_MEM_HOST, CEED_COPY_VALUES, 326 (PetscScalar *)X_loc_array); 327 ierr = VecRestoreArrayRead(X_loc, &X_loc_array); CHKERRQ(ierr); 328 329 // ----------------------------------------------------------------------------- 330 // CEED vectors 331 // ----------------------------------------------------------------------------- 332 // -- Create CEED vector for geometric data 333 CeedInt num_qpts_vol; 334 PetscInt loc_num_elem_vol; 335 CeedBasisGetNumQuadraturePoints(ceed_data->basis_q, &num_qpts_vol); 336 CeedElemRestrictionGetNumElements(ceed_data->elem_restr_q, &loc_num_elem_vol); 337 CeedVectorCreate(ceed, q_data_size_vol*loc_num_elem_vol*num_qpts_vol, 338 &ceed_data->q_data); 339 340 // ----------------------------------------------------------------------------- 341 // CEED Operators 342 // ----------------------------------------------------------------------------- 343 // -- Create CEED operator for quadrature data 344 CeedOperatorCreate(ceed, ceed_data->qf_setup_vol, NULL, NULL, 345 &ceed_data->op_setup_vol); 346 CeedOperatorSetField(ceed_data->op_setup_vol, "dx", ceed_data->elem_restr_x, 347 ceed_data->basis_x, CEED_VECTOR_ACTIVE); 348 CeedOperatorSetField(ceed_data->op_setup_vol, "weight", 349 CEED_ELEMRESTRICTION_NONE, ceed_data->basis_x, CEED_VECTOR_NONE); 350 CeedOperatorSetField(ceed_data->op_setup_vol, "qdata", 351 ceed_data->elem_restr_qd_i, CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 352 353 // -- Create CEED operator for ICs 354 CeedOperatorCreate(ceed, ceed_data->qf_ics, NULL, NULL, &ceed_data->op_ics); 355 CeedOperatorSetField(ceed_data->op_ics, "x", ceed_data->elem_restr_x, 356 ceed_data->basis_xc, CEED_VECTOR_ACTIVE); 357 CeedOperatorSetField(ceed_data->op_ics, "q0", ceed_data->elem_restr_q, 358 CEED_BASIS_COLLOCATED, CEED_VECTOR_ACTIVE); 359 360 // Create CEED operator for RHS 361 if (ceed_data->qf_rhs_vol) { 362 CeedOperator op; 363 CeedOperatorCreate(ceed, ceed_data->qf_rhs_vol, NULL, NULL, &op); 364 CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 365 CEED_VECTOR_ACTIVE); 366 CeedOperatorSetField(op, "dq", ceed_data->elem_restr_q, ceed_data->basis_q, 367 CEED_VECTOR_ACTIVE); 368 CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 369 CEED_BASIS_COLLOCATED, ceed_data->q_data); 370 CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 371 ceed_data->x_coord); 372 CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 373 CEED_VECTOR_ACTIVE); 374 CeedOperatorSetField(op, "dv", ceed_data->elem_restr_q, ceed_data->basis_q, 375 CEED_VECTOR_ACTIVE); 376 user->op_rhs_vol = op; 377 } 378 379 // -- CEED operator for IFunction 380 if (ceed_data->qf_ifunction_vol) { 381 CeedOperator op; 382 CeedOperatorCreate(ceed, ceed_data->qf_ifunction_vol, NULL, NULL, &op); 383 CeedOperatorSetField(op, "q", ceed_data->elem_restr_q, ceed_data->basis_q, 384 CEED_VECTOR_ACTIVE); 385 CeedOperatorSetField(op, "dq", ceed_data->elem_restr_q, ceed_data->basis_q, 386 CEED_VECTOR_ACTIVE); 387 CeedOperatorSetField(op, "q dot", ceed_data->elem_restr_q, ceed_data->basis_q, 388 user->q_dot_ceed); 389 CeedOperatorSetField(op, "qdata", ceed_data->elem_restr_qd_i, 390 CEED_BASIS_COLLOCATED, ceed_data->q_data); 391 CeedOperatorSetField(op, "x", ceed_data->elem_restr_x, ceed_data->basis_x, 392 ceed_data->x_coord); 393 CeedOperatorSetField(op, "v", ceed_data->elem_restr_q, ceed_data->basis_q, 394 CEED_VECTOR_ACTIVE); 395 CeedOperatorSetField(op, "dv", ceed_data->elem_restr_q, ceed_data->basis_q, 396 CEED_VECTOR_ACTIVE); 397 user->op_ifunction_vol = op; 398 } 399 400 // ***************************************************************************** 401 // Set up CEED objects for the exterior domain (surface) 402 // ***************************************************************************** 403 CeedInt height = 1, 404 dim_sur = dim - height, 405 P_sur = app_ctx->degree + 1, 406 Q_sur = P_sur + app_ctx->q_extra; 407 const CeedInt q_data_size_sur = problem->q_data_size_sur; 408 409 // ----------------------------------------------------------------------------- 410 // CEED Bases 411 // ----------------------------------------------------------------------------- 412 CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_q, P_sur, Q_sur, 413 CEED_GAUSS, &ceed_data->basis_q_sur); 414 CeedBasisCreateTensorH1Lagrange(ceed, dim_sur, num_comp_x, 2, Q_sur, CEED_GAUSS, 415 &ceed_data->basis_x_sur); 416 417 // ----------------------------------------------------------------------------- 418 // CEED QFunctions 419 // ----------------------------------------------------------------------------- 420 // -- Create QFunction for quadrature data 421 CeedQFunctionCreateInterior(ceed, 1, problem->setup_sur, problem->setup_sur_loc, 422 &ceed_data->qf_setup_sur); 423 CeedQFunctionAddInput(ceed_data->qf_setup_sur, "dx", num_comp_x*dim_sur, 424 CEED_EVAL_GRAD); 425 CeedQFunctionAddInput(ceed_data->qf_setup_sur, "weight", 1, CEED_EVAL_WEIGHT); 426 CeedQFunctionAddOutput(ceed_data->qf_setup_sur, "surface qdata", 427 q_data_size_sur, CEED_EVAL_NONE); 428 429 // -- Creat QFunction for inflow boundaries 430 if (problem->apply_inflow) { 431 CeedQFunctionCreateInterior(ceed, 1, problem->apply_inflow, 432 problem->apply_inflow_loc, &ceed_data->qf_apply_inflow); 433 CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "q", num_comp_q, 434 CEED_EVAL_INTERP); 435 CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "surface qdata", 436 q_data_size_sur, CEED_EVAL_NONE); 437 CeedQFunctionAddInput(ceed_data->qf_apply_inflow, "x", num_comp_x, 438 CEED_EVAL_INTERP); 439 CeedQFunctionAddOutput(ceed_data->qf_apply_inflow, "v", num_comp_q, 440 CEED_EVAL_INTERP); 441 } 442 443 // -- Creat QFunction for outflow boundaries 444 if (problem->apply_outflow) { 445 CeedQFunctionCreateInterior(ceed, 1, problem->apply_outflow, 446 problem->apply_outflow_loc, &ceed_data->qf_apply_outflow); 447 CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "q", num_comp_q, 448 CEED_EVAL_INTERP); 449 CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "surface qdata", 450 q_data_size_sur, CEED_EVAL_NONE); 451 CeedQFunctionAddInput(ceed_data->qf_apply_outflow, "x", num_comp_x, 452 CEED_EVAL_INTERP); 453 CeedQFunctionAddOutput(ceed_data->qf_apply_outflow, "v", num_comp_q, 454 CEED_EVAL_INTERP); 455 } 456 457 // ***************************************************************************** 458 // CEED Operator Apply 459 // ***************************************************************************** 460 // -- Apply CEED Operator for the geometric data 461 CeedOperatorApply(ceed_data->op_setup_vol, ceed_data->x_coord, 462 ceed_data->q_data, CEED_REQUEST_IMMEDIATE); 463 464 // -- Create and apply CEED Composite Operator for the entire domain 465 if (!user->phys->implicit) { // RHS 466 ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 467 user->op_rhs_vol, height, P_sur, Q_sur, 468 q_data_size_sur, &user->op_rhs); CHKERRQ(ierr); 469 } else { // IFunction 470 ierr = CreateOperatorForDomain(ceed, dm, bc, ceed_data, user->phys, 471 user->op_ifunction_vol, height, P_sur, Q_sur, 472 q_data_size_sur, &user->op_ifunction); CHKERRQ(ierr); 473 } 474 475 PetscFunctionReturn(0); 476 } 477