Lines Matching +full:- +full:- +full:ceed
1 // Copyright (c) 2017-2026, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
4 // SPDX-License-Identifier: BSD-2-Clause
6 // This file is part of CEED: http://github.com/ceed
9 /// Setup libCEED for Navier-Stokes example using PETSc
11 #include <ceed.h>
18 Ceed ceed = user->ceed; in CreateKSPMassOperator_Unstabilized() local
31 PetscCallCeed(ceed, CeedOperatorCompositeGetSubList(user->op_rhs_ctx->op, &sub_ops)); in CreateKSPMassOperator_Unstabilized()
32 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "q", &field)); in CreateKSPMassOperator_Unstabilized()
33 PetscCallCeed(ceed, CeedOperatorFieldGetData(field, NULL, &elem_restr_q, &basis_q, NULL)); in CreateKSPMassOperator_Unstabilized()
35 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "qdata", &field)); in CreateKSPMassOperator_Unstabilized()
36 PetscCallCeed(ceed, CeedOperatorFieldGetData(field, NULL, &elem_restr_qd_i, NULL, &q_data)); in CreateKSPMassOperator_Unstabilized()
39 PetscCallCeed(ceed, CeedElemRestrictionGetNumComponents(elem_restr_q, &num_comp_q)); in CreateKSPMassOperator_Unstabilized()
40 PetscCallCeed(ceed, CeedElemRestrictionGetNumComponents(elem_restr_qd_i, &q_data_size)); in CreateKSPMassOperator_Unstabilized()
42 PetscCall(CreateMassQFunction(ceed, num_comp_q, q_data_size, &qf_mass)); in CreateKSPMassOperator_Unstabilized()
43 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_mass, NULL, NULL, op_mass)); in CreateKSPMassOperator_Unstabilized()
44 …PetscCallCeed(ceed, CeedOperatorSetField(*op_mass, "u", elem_restr_q, basis_q, CEED_VECTOR_ACTIVE)… in CreateKSPMassOperator_Unstabilized()
45 …PetscCallCeed(ceed, CeedOperatorSetField(*op_mass, "qdata", elem_restr_qd_i, CEED_BASIS_NONE, q_da… in CreateKSPMassOperator_Unstabilized()
46 …PetscCallCeed(ceed, CeedOperatorSetField(*op_mass, "v", elem_restr_q, basis_q, CEED_VECTOR_ACTIVE)… in CreateKSPMassOperator_Unstabilized()
48 PetscCallCeed(ceed, CeedVectorDestroy(&q_data)); in CreateKSPMassOperator_Unstabilized()
49 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_q)); in CreateKSPMassOperator_Unstabilized()
50 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_qd_i)); in CreateKSPMassOperator_Unstabilized()
51 PetscCallCeed(ceed, CeedBasisDestroy(&basis_q)); in CreateKSPMassOperator_Unstabilized()
52 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_mass)); in CreateKSPMassOperator_Unstabilized()
58 Ceed ceed = user->ceed; in CreateKSPMass() local
59 DM dm = user->dm; in CreateKSPMass()
63 if (problem->create_mass_operator) PetscCall(problem->create_mass_operator(user, &op_mass)); in CreateKSPMass()
66 { // -- Setup KSP for mass operator in CreateKSPMass()
76 PetscCall(KSPCreate(comm, &user->mass_ksp)); in CreateKSPMass()
77 PetscCall(KSPSetOptionsPrefix(user->mass_ksp, "mass_")); in CreateKSPMass()
80 PetscCall(KSPGetPC(user->mass_ksp, &pc)); in CreateKSPMass()
83 PetscCall(KSPSetType(user->mass_ksp, KSPPREONLY)); in CreateKSPMass()
85 PetscCall(KSPSetFromOptions_WithMatCeed(user->mass_ksp, mat_mass)); in CreateKSPMass()
90 PetscCallCeed(ceed, CeedOperatorDestroy(&op_mass)); in CreateKSPMass()
94 static PetscErrorCode AddBCSubOperator(Ceed ceed, DM dm, CeedData ceed_data, DMLabel domain_label, … in AddBCSubOperator() argument
104 …PetscCall(DMPlexCeedElemRestrictionCreate(ceed, dm, domain_label, label_value, height, dm_field, &… in AddBCSubOperator()
105 …PetscCall(DMPlexCeedElemRestrictionCoordinateCreate(ceed, dm, domain_label, label_value, height, &… in AddBCSubOperator()
107 // State-dependent data will be passed from residual to Jacobian. This will be collocated. in AddBCSubOperator()
108 …PetscCall(DMPlexCeedElemRestrictionQDataCreate(ceed, dm, domain_label, label_value, height, jac_da… in AddBCSubOperator()
109 PetscCallCeed(ceed, CeedElemRestrictionCreateVector(elem_restr_jd_i_sur, &jac_data_sur, NULL)); in AddBCSubOperator()
112 …PetscCall(QDataBoundaryGet(ceed, dm, domain_label, label_value, elem_restr_x_sur, basis_x_sur, cee… in AddBCSubOperator()
115 // CEED Operator for Physics in AddBCSubOperator()
116 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_apply_bc, NULL, NULL, &op_apply_bc)); in AddBCSubOperator()
117 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "q", elem_restr_q_sur, basis_q_sur, CEED_VEC… in AddBCSubOperator()
118 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "Grad_q", elem_restr_q_sur, basis_q_sur, CEE… in AddBCSubOperator()
119 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "surface qdata", elem_restr_qd_i_sur, CEED_B… in AddBCSubOperator()
120 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "x", elem_restr_x_sur, basis_x_sur, ceed_dat… in AddBCSubOperator()
121 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "v", elem_restr_q_sur, basis_q_sur, CEED_VEC… in AddBCSubOperator()
123 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc, "surface jacobian data", elem_restr_jd_i_sur… in AddBCSubOperator()
127 …PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_apply_bc_jacobian, NULL, NULL, &op_apply_bc_jacobi… in AddBCSubOperator()
128 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "dq", elem_restr_q_sur, basis_q_sur… in AddBCSubOperator()
129 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "Grad_dq", elem_restr_q_sur, basis_… in AddBCSubOperator()
130 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "surface qdata", elem_restr_qd_i_su… in AddBCSubOperator()
131 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "x", elem_restr_x_sur, basis_x_sur,… in AddBCSubOperator()
132 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "surface jacobian data", elem_restr… in AddBCSubOperator()
133 …PetscCallCeed(ceed, CeedOperatorSetField(op_apply_bc_jacobian, "v", elem_restr_q_sur, basis_q_sur,… in AddBCSubOperator()
136 // Apply Sub-Operator for Physics in AddBCSubOperator()
137 PetscCallCeed(ceed, CeedOperatorCompositeAddSub(op_apply, op_apply_bc)); in AddBCSubOperator()
138 …if (op_apply_bc_jacobian) PetscCallCeed(ceed, CeedOperatorCompositeAddSub(op_apply_ijacobian, op_a… in AddBCSubOperator()
140 PetscCallCeed(ceed, CeedVectorDestroy(&q_data_sur)); in AddBCSubOperator()
141 PetscCallCeed(ceed, CeedVectorDestroy(&jac_data_sur)); in AddBCSubOperator()
142 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_q_sur)); in AddBCSubOperator()
143 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_x_sur)); in AddBCSubOperator()
144 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_qd_i_sur)); in AddBCSubOperator()
145 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_jd_i_sur)); in AddBCSubOperator()
146 PetscCallCeed(ceed, CeedOperatorDestroy(&op_apply_bc)); in AddBCSubOperator()
147 PetscCallCeed(ceed, CeedOperatorDestroy(&op_apply_bc_jacobian)); in AddBCSubOperator()
151 static PetscErrorCode SetupBCQFunctions(Ceed ceed, PetscInt dim_sur, PetscInt num_comp_x, PetscInt … in SetupBCQFunctions() argument
156 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, apply_bc.qfunction, apply_bc.qfunction_lo… in SetupBCQFunctions()
157 PetscCallCeed(ceed, CeedQFunctionSetContext(*qf_apply_bc, apply_bc.qfunction_context)); in SetupBCQFunctions()
158 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(*qf_apply_bc, 0)); in SetupBCQFunctions()
159 PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc, "q", num_comp_q, CEED_EVAL_INTERP)); in SetupBCQFunctions()
160 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc, "Grad_q", num_comp_q * dim_sur, CEED_EVAL_… in SetupBCQFunctions()
161 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc, "surface qdata", q_data_size_sur, CEED_EVA… in SetupBCQFunctions()
162 PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc, "x", num_comp_x, CEED_EVAL_INTERP)); in SetupBCQFunctions()
163 PetscCallCeed(ceed, CeedQFunctionAddOutput(*qf_apply_bc, "v", num_comp_q, CEED_EVAL_INTERP)); in SetupBCQFunctions()
164 …if (jac_data_size_sur) PetscCallCeed(ceed, CeedQFunctionAddOutput(*qf_apply_bc, "surface jacobian … in SetupBCQFunctions()
167 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, apply_bc_jacobian.qfunction, apply_bc_jac… in SetupBCQFunctions()
168 …PetscCallCeed(ceed, CeedQFunctionSetContext(*qf_apply_bc_jacobian, apply_bc_jacobian.qfunction_con… in SetupBCQFunctions()
169 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(*qf_apply_bc_jacobian, 0)); in SetupBCQFunctions()
170 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc_jacobian, "dq", num_comp_q, CEED_EVAL_INTER… in SetupBCQFunctions()
171 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc_jacobian, "Grad_dq", num_comp_q * dim_sur, … in SetupBCQFunctions()
172 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc_jacobian, "surface qdata", q_data_size_sur,… in SetupBCQFunctions()
173 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc_jacobian, "x", num_comp_x, CEED_EVAL_INTERP… in SetupBCQFunctions()
174 …PetscCallCeed(ceed, CeedQFunctionAddInput(*qf_apply_bc_jacobian, "surface jacobian data", jac_data… in SetupBCQFunctions()
175 …PetscCallCeed(ceed, CeedQFunctionAddOutput(*qf_apply_bc_jacobian, "v", num_comp_q, CEED_EVAL_INTER… in SetupBCQFunctions()
181 static PetscErrorCode AddBCSubOperators(User user, Ceed ceed, DM dm, SimpleBC bc, ProblemData probl… in AddBCSubOperators() argument
184 …CeedInt P_sur = user->app_ctx->degree + 1, Q_sur = P_sur + user->app_ctx->q_extra, dim_sur, … in AddBCSubOperators()
185 const CeedInt jac_data_size_sur = user->phys->implicit ? problem->jac_data_size_sur : 0; in AddBCSubOperators()
193 dim_sur = dim - height; in AddBCSubOperators()
200 PetscCallCeed(ceed, CeedOperatorCompositeGetSubList(op_apply, &sub_ops)); in AddBCSubOperators()
201 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "q", &field)); in AddBCSubOperators()
202 PetscCallCeed(ceed, CeedOperatorFieldGetElemRestriction(field, &elem_restr_q)); in AddBCSubOperators()
203 PetscCallCeed(ceed, CeedElemRestrictionGetNumComponents(elem_restr_q, &num_comp_q)); in AddBCSubOperators()
204 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_q)); in AddBCSubOperators()
206 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "x", &field)); in AddBCSubOperators()
207 PetscCallCeed(ceed, CeedOperatorFieldGetElemRestriction(field, &elem_restr_x)); in AddBCSubOperators()
208 PetscCallCeed(ceed, CeedElemRestrictionGetNumComponents(elem_restr_x, &num_comp_x)); in AddBCSubOperators()
209 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_x)); in AddBCSubOperators()
219 PetscCall(CreateBasisFromPlex(ceed, dm, label, label_value, height, field, &basis_q_sur)); in AddBCSubOperators()
220 PetscCall(CreateBasisFromPlex(ceed, dm_coord, label, label_value, height, field, &basis_x_sur)); in AddBCSubOperators()
225 { // --- Create Sub-Operator for inflow boundaries in AddBCSubOperators()
228 …PetscCall(SetupBCQFunctions(ceed, dim_sur, num_comp_x, num_comp_q, q_data_size_sur, jac_data_size_… in AddBCSubOperators()
229 … problem->apply_inflow_jacobian, &qf_apply_inflow, &qf_apply_inflow_jacobian)); in AddBCSubOperators()
230 for (CeedInt i = 0; i < bc->num_inflow; i++) { in AddBCSubOperators()
231 …PetscCall(AddBCSubOperator(ceed, dm, ceed_data, face_sets_label, bc->inflows[i], height, Q_sur, q_… in AddBCSubOperators()
234 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_inflow)); in AddBCSubOperators()
235 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_inflow_jacobian)); in AddBCSubOperators()
238 { // --- Create Sub-Operator for outflow boundaries in AddBCSubOperators()
241 …PetscCall(SetupBCQFunctions(ceed, dim_sur, num_comp_x, num_comp_q, q_data_size_sur, jac_data_size_… in AddBCSubOperators()
242 … problem->apply_outflow_jacobian, &qf_apply_outflow, &qf_apply_outflow_jacobian)); in AddBCSubOperators()
243 for (CeedInt i = 0; i < bc->num_outflow; i++) { in AddBCSubOperators()
244 …PetscCall(AddBCSubOperator(ceed, dm, ceed_data, face_sets_label, bc->outflows[i], height, Q_sur, q… in AddBCSubOperators()
247 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_outflow)); in AddBCSubOperators()
248 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_outflow_jacobian)); in AddBCSubOperators()
251 { // --- Create Sub-Operator for freestream boundaries in AddBCSubOperators()
254 …PetscCall(SetupBCQFunctions(ceed, dim_sur, num_comp_x, num_comp_q, q_data_size_sur, jac_data_size_… in AddBCSubOperators()
255 … problem->apply_freestream_jacobian, &qf_apply_freestream, &qf_apply_freestream_jacobian)); in AddBCSubOperators()
256 for (CeedInt i = 0; i < bc->num_freestream; i++) { in AddBCSubOperators()
257 …PetscCall(AddBCSubOperator(ceed, dm, ceed_data, face_sets_label, bc->freestreams[i], height, Q_sur… in AddBCSubOperators()
260 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_freestream)); in AddBCSubOperators()
261 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_freestream_jacobian)); in AddBCSubOperators()
264 { // --- Create Sub-Operator for slip boundaries in AddBCSubOperators()
267 …PetscCall(SetupBCQFunctions(ceed, dim_sur, num_comp_x, num_comp_q, q_data_size_sur, jac_data_size_… in AddBCSubOperators()
268 … problem->apply_slip_jacobian, &qf_apply_slip, &qf_apply_slip_jacobian)); in AddBCSubOperators()
269 for (CeedInt i = 0; i < bc->num_slip; i++) { in AddBCSubOperators()
270 …PetscCall(AddBCSubOperator(ceed, dm, ceed_data, face_sets_label, bc->slips[i], height, Q_sur, q_da… in AddBCSubOperators()
273 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_slip)); in AddBCSubOperators()
274 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_apply_slip_jacobian)); in AddBCSubOperators()
277 PetscCallCeed(ceed, CeedBasisDestroy(&basis_q_sur)); in AddBCSubOperators()
278 PetscCallCeed(ceed, CeedBasisDestroy(&basis_x_sur)); in AddBCSubOperators()
282 PetscErrorCode SetupLibceed(Ceed ceed, CeedData ceed_data, DM dm, User user, AppCtx app_ctx, Proble… in SetupLibceed() argument
284 const CeedInt dim = problem->dim, num_comp_x = problem->dim; in SetupLibceed()
292 if (problem->apply_vol_ifunction.qfunction && problem->uses_newtonian) { in SetupLibceed()
294 …PetscCallCeed(ceed, CeedQFunctionContextGetDataRead(problem->apply_vol_ifunction.qfunction_context… in SetupLibceed()
295 jac_data_size_vol += (gas->idl_enable ? 1 : 0); in SetupLibceed()
296 …PetscCallCeed(ceed, CeedQFunctionContextRestoreDataRead(problem->apply_vol_ifunction.qfunction_con… in SetupLibceed()
305 …PetscCall(CreateBasisFromPlex(ceed, dm, domain_label, label_value, height, dm_field, &ceed_data->b… in SetupLibceed()
306 …PetscCall(CreateBasisFromPlex(ceed, dm_coord, domain_label, label_value, height, dm_field, &ceed_d… in SetupLibceed()
308 …PetscCall(DMPlexCeedElemRestrictionCreate(ceed, dm, domain_label, label_value, height, 0, &ceed_da… in SetupLibceed()
309 …cCall(DMPlexCeedElemRestrictionCoordinateCreate(ceed, dm, domain_label, label_value, height, &ceed… in SetupLibceed()
310 …PetscCall(DMPlexCeedElemRestrictionQDataCreate(ceed, dm, domain_label, label_value, height, jac_da… in SetupLibceed()
312 …PetscCallCeed(ceed, CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_ceed, NULL)); in SetupLibceed()
313 …PetscCallCeed(ceed, CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->q_dot_ceed, NU… in SetupLibceed()
314 …PetscCallCeed(ceed, CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &user->g_ceed, NULL)); in SetupLibceed()
315 …PetscCallCeed(ceed, CeedElemRestrictionCreateVector(ceed_data->elem_restr_x, &ceed_data->x_coord, … in SetupLibceed()
316 PetscCallCeed(ceed, CeedElemRestrictionCreateVector(elem_restr_jd_i, &jac_data, NULL)); in SetupLibceed()
318 { // -- Copy PETSc coordinate vector into CEED vector in SetupLibceed()
328 PetscCall(VecScale(X_loc, problem->dm_scale)); in SetupLibceed()
329 PetscCall(VecCopyPetscToCeed(X_loc, ceed_data->x_coord)); in SetupLibceed()
332 …PetscCall(QDataGet(ceed, dm, domain_label, label_value, ceed_data->elem_restr_x, ceed_data->basis_… in SetupLibceed()
333 &ceed_data->elem_restr_qd_i, &ceed_data->q_data, &problem->q_data_size_vol)); in SetupLibceed()
336 { // -- Create QFunction for ICs in SetupLibceed()
341 … PetscCallCeed(ceed, CeedBasisCreateProjection(ceed_data->basis_x, ceed_data->basis_q, &basis_xc)); in SetupLibceed()
342 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, problem->ics.qfunction, problem->ics.qfun… in SetupLibceed()
343 PetscCallCeed(ceed, CeedQFunctionSetContext(qf_ics, problem->ics.qfunction_context)); in SetupLibceed()
344 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(qf_ics, 0)); in SetupLibceed()
345 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ics, "x", num_comp_x, CEED_EVAL_INTERP)); in SetupLibceed()
346 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ics, "dx", num_comp_x * dim, CEED_EVAL_GRAD)); in SetupLibceed()
347 PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ics, "q0", num_comp_q, CEED_EVAL_NONE)); in SetupLibceed()
349 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_ics, NULL, NULL, &op_ics)); in SetupLibceed()
350 …PetscCallCeed(ceed, CeedOperatorSetField(op_ics, "x", ceed_data->elem_restr_x, basis_xc, CEED_VECT… in SetupLibceed()
351 …PetscCallCeed(ceed, CeedOperatorSetField(op_ics, "dx", ceed_data->elem_restr_x, basis_xc, CEED_VEC… in SetupLibceed()
352 …PetscCallCeed(ceed, CeedOperatorSetField(op_ics, "q0", ceed_data->elem_restr_q, CEED_BASIS_NONE, C… in SetupLibceed()
353 …PetscCallCeed(ceed, CeedOperatorGetContextFieldLabel(op_ics, "evaluation time", &user->phys->ics_t… in SetupLibceed()
354 …peratorApplyContextCreate(NULL, dm, user->ceed, op_ics, ceed_data->x_coord, NULL, NULL, user->Q_lo… in SetupLibceed()
356 PetscCallCeed(ceed, CeedBasisDestroy(&basis_xc)); in SetupLibceed()
357 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_ics)); in SetupLibceed()
358 PetscCallCeed(ceed, CeedOperatorDestroy(&op_ics)); in SetupLibceed()
361 if (problem->apply_vol_rhs.qfunction) { in SetupLibceed()
364 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_rhs.qfunction, problem… in SetupLibceed()
365 …PetscCallCeed(ceed, CeedQFunctionSetContext(qf_rhs_vol, problem->apply_vol_rhs.qfunction_context)); in SetupLibceed()
366 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(qf_rhs_vol, 0)); in SetupLibceed()
367 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_rhs_vol, "q", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
368 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_rhs_vol, "Grad_q", num_comp_q * dim, CEED_EVAL_GRAD)); in SetupLibceed()
369 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_rhs_vol, "qdata", problem->q_data_size_vol, CEED_EVAL… in SetupLibceed()
370 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_rhs_vol, "x", num_comp_x, CEED_EVAL_INTERP)); in SetupLibceed()
371 PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_rhs_vol, "v", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
372 …PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_rhs_vol, "Grad_v", num_comp_q * dim, CEED_EVAL_GRAD)… in SetupLibceed()
374 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_rhs_vol, NULL, NULL, &op_rhs_vol)); in SetupLibceed()
375 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "q", ceed_data->elem_restr_q, ceed_data->basi… in SetupLibceed()
376 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "Grad_q", ceed_data->elem_restr_q, ceed_data-… in SetupLibceed()
377 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "qdata", ceed_data->elem_restr_qd_i, CEED_BAS… in SetupLibceed()
378 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "x", ceed_data->elem_restr_x, ceed_data->basi… in SetupLibceed()
379 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "v", ceed_data->elem_restr_q, ceed_data->basi… in SetupLibceed()
380 …PetscCallCeed(ceed, CeedOperatorSetField(op_rhs_vol, "Grad_v", ceed_data->elem_restr_q, ceed_data-… in SetupLibceed()
382 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_rhs_vol)); in SetupLibceed()
385 if (problem->apply_vol_ifunction.qfunction) { in SetupLibceed()
388 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ifunction.qfunction, p… in SetupLibceed()
390 …PetscCallCeed(ceed, CeedQFunctionSetContext(qf_ifunction_vol, problem->apply_vol_ifunction.qfuncti… in SetupLibceed()
391 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(qf_ifunction_vol, 0)); in SetupLibceed()
392 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ifunction_vol, "q", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
393 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ifunction_vol, "Grad_q", num_comp_q * dim, CEED_EVAL_… in SetupLibceed()
394 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ifunction_vol, "q dot", num_comp_q, CEED_EVAL_INTERP)… in SetupLibceed()
395 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ifunction_vol, "qdata", problem->q_data_size_vol, CEE… in SetupLibceed()
396 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ifunction_vol, "x", num_comp_x, CEED_EVAL_INTERP)); in SetupLibceed()
397 … PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ifunction_vol, "v", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
398 …PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ifunction_vol, "Grad_v", num_comp_q * dim, CEED_EVAL… in SetupLibceed()
399 …PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ifunction_vol, "jac_data", jac_data_size_vol, CEED_E… in SetupLibceed()
401 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_ifunction_vol, NULL, NULL, &op_ifunction_vol)); in SetupLibceed()
402 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "q", ceed_data->elem_restr_q, ceed_data… in SetupLibceed()
403 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "Grad_q", ceed_data->elem_restr_q, ceed… in SetupLibceed()
404 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "q dot", ceed_data->elem_restr_q, ceed_… in SetupLibceed()
405 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "qdata", ceed_data->elem_restr_qd_i, CE… in SetupLibceed()
406 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "x", ceed_data->elem_restr_x, ceed_data… in SetupLibceed()
407 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "v", ceed_data->elem_restr_q, ceed_data… in SetupLibceed()
408 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "Grad_v", ceed_data->elem_restr_q, ceed… in SetupLibceed()
409 …PetscCallCeed(ceed, CeedOperatorSetField(op_ifunction_vol, "jac_data", elem_restr_jd_i, CEED_BASIS… in SetupLibceed()
411 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_ifunction_vol)); in SetupLibceed()
414 if (problem->apply_vol_ijacobian.qfunction) { in SetupLibceed()
417 …PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, problem->apply_vol_ijacobian.qfunction, p… in SetupLibceed()
419 …PetscCallCeed(ceed, CeedQFunctionSetContext(qf_ijacobian_vol, problem->apply_vol_ijacobian.qfuncti… in SetupLibceed()
420 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(qf_ijacobian_vol, 0)); in SetupLibceed()
421 … PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ijacobian_vol, "dq", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
422 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ijacobian_vol, "Grad_dq", num_comp_q * dim, CEED_EVAL… in SetupLibceed()
423 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ijacobian_vol, "qdata", problem->q_data_size_vol, CEE… in SetupLibceed()
424 …PetscCallCeed(ceed, CeedQFunctionAddInput(qf_ijacobian_vol, "jac_data", jac_data_size_vol, CEED_EV… in SetupLibceed()
425 … PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ijacobian_vol, "v", num_comp_q, CEED_EVAL_INTERP)); in SetupLibceed()
426 …PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_ijacobian_vol, "Grad_v", num_comp_q * dim, CEED_EVAL… in SetupLibceed()
428 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_ijacobian_vol, NULL, NULL, &op_ijacobian_vol)); in SetupLibceed()
429 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "dq", ceed_data->elem_restr_q, ceed_dat… in SetupLibceed()
430 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "Grad_dq", ceed_data->elem_restr_q, cee… in SetupLibceed()
431 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "qdata", ceed_data->elem_restr_qd_i, CE… in SetupLibceed()
432 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "jac_data", elem_restr_jd_i, CEED_BASIS… in SetupLibceed()
433 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "v", ceed_data->elem_restr_q, ceed_data… in SetupLibceed()
434 …PetscCallCeed(ceed, CeedOperatorSetField(op_ijacobian_vol, "Grad_v", ceed_data->elem_restr_q, ceed… in SetupLibceed()
436 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_ijacobian_vol)); in SetupLibceed()
439 // -- Create and apply CEED Composite Operator for the entire domain in SetupLibceed()
440 if (!user->phys->implicit) { // RHS in SetupLibceed()
443 PetscCallCeed(ceed, CeedOperatorCreateComposite(ceed, &op_rhs)); in SetupLibceed()
444 PetscCallCeed(ceed, CeedOperatorCompositeAddSub(op_rhs, op_rhs_vol)); in SetupLibceed()
445 PetscCall(AddBCSubOperators(user, ceed, dm, bc, problem, ceed_data, op_rhs, NULL)); in SetupLibceed()
447 …tscCall(OperatorApplyContextCreate(dm, dm, ceed, op_rhs, user->q_ceed, user->g_ceed, user->Q_loc, … in SetupLibceed()
449 // ----- Get Context Labels for Operator in SetupLibceed()
450 …PetscCallCeed(ceed, CeedOperatorGetContextFieldLabel(op_rhs, "solution time", &user->phys->solutio… in SetupLibceed()
451 …PetscCallCeed(ceed, CeedOperatorGetContextFieldLabel(op_rhs, "timestep size", &user->phys->timeste… in SetupLibceed()
453 PetscCallCeed(ceed, CeedOperatorDestroy(&op_rhs)); in SetupLibceed()
459 PetscCallCeed(ceed, CeedOperatorCreateComposite(ceed, &user->op_ifunction)); in SetupLibceed()
460 PetscCallCeed(ceed, CeedOperatorCompositeAddSub(user->op_ifunction, op_ifunction_vol)); in SetupLibceed()
462 PetscCallCeed(ceed, CeedOperatorCreateComposite(ceed, &op_ijacobian)); in SetupLibceed()
463 PetscCallCeed(ceed, CeedOperatorCompositeAddSub(op_ijacobian, op_ijacobian_vol)); in SetupLibceed()
465 …PetscCall(AddBCSubOperators(user, ceed, dm, bc, problem, ceed_data, user->op_ifunction, op_ijacobi… in SetupLibceed()
467 // ----- Get Context Labels for Operator in SetupLibceed()
468 …PetscCallCeed(ceed, CeedOperatorGetContextFieldLabel(user->op_ifunction, "solution time", &user->p… in SetupLibceed()
469 …PetscCallCeed(ceed, CeedOperatorGetContextFieldLabel(user->op_ifunction, "timestep size", &user->p… in SetupLibceed()
472 PetscCall(MatCreateCeed(user->dm, user->dm, op_ijacobian, NULL, &user->mat_ijacobian)); in SetupLibceed()
473 PetscCall(MatCeedSetLocalVectors(user->mat_ijacobian, user->Q_dot_loc, NULL)); in SetupLibceed()
474 PetscCallCeed(ceed, CeedOperatorDestroy(&op_ijacobian)); in SetupLibceed()
478 …if (problem->use_strong_bc_ceed) PetscCall(SetupStrongBC_Ceed(ceed, ceed_data, dm, user, problem, … in SetupLibceed()
479 …if (app_ctx->turb_spanstats_enable) PetscCall(TurbulenceStatisticsSetup(ceed, user, ceed_data, pro… in SetupLibceed()
480 …if (app_ctx->diff_filter_monitor && !user->diff_filter) PetscCall(DifferentialFilterSetup(ceed, us… in SetupLibceed()
482 PetscCallCeed(ceed, CeedVectorDestroy(&jac_data)); in SetupLibceed()
483 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_jd_i)); in SetupLibceed()
484 PetscCallCeed(ceed, CeedOperatorDestroy(&op_ijacobian_vol)); in SetupLibceed()
485 PetscCallCeed(ceed, CeedOperatorDestroy(&op_ifunction_vol)); in SetupLibceed()
486 PetscCallCeed(ceed, CeedOperatorDestroy(&op_rhs_vol)); in SetupLibceed()