1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 /// @file 4 /// Functions for setting up and projecting the divergence of the diffusive flux 5 6 #include "../qfunctions/diff_flux_projection.h" 7 8 #include <petscdmplex.h> 9 10 #include <navierstokes.h> 11 12 /** 13 @brief Create `DivDiffFluxProjectionData` for solution DM in `user` 14 15 @param[in] user `User` context 16 @param[in] num_diff_flux_comps Number of components that makes up the diffusive flux (e.g. 1 for scalar advection-diffusion) 17 @param[out] pdiff_flux_proj The `DivDiffFluxProjectionData` object created, or `NULL` if not created 18 **/ 19 PetscErrorCode DivDiffFluxProjectionCreate(User user, PetscInt num_diff_flux_comps, DivDiffFluxProjectionData *pdiff_flux_proj) { 20 PetscInt label_value = 0, height = 0, dm_field = 0, dim, degree = user->app_ctx->degree, q_extra = user->app_ctx->q_extra; 21 DMLabel domain_label = NULL; 22 DivDiffFluxProjectionData diff_flux_proj; 23 NodalProjectionData projection; 24 25 PetscFunctionBeginUser; 26 if (user->app_ctx->divFdiffproj_method == DIV_DIFF_FLUX_PROJ_NONE) { 27 *pdiff_flux_proj = NULL; 28 PetscFunctionReturn(PETSC_SUCCESS); 29 } 30 PetscCall(PetscNew(pdiff_flux_proj)); 31 diff_flux_proj = *pdiff_flux_proj; 32 PetscCall(PetscNew(&user->diff_flux_proj->projection)); 33 projection = user->diff_flux_proj->projection; 34 diff_flux_proj->method = user->app_ctx->divFdiffproj_method; 35 diff_flux_proj->num_diff_flux_comps = num_diff_flux_comps; 36 37 PetscCall(DMClone(user->dm, &projection->dm)); 38 PetscCall(DMSetMatrixPreallocateSkip(projection->dm, PETSC_TRUE)); 39 PetscCall(DMGetDimension(projection->dm, &dim)); 40 switch (diff_flux_proj->method) { 41 case DIV_DIFF_FLUX_PROJ_DIRECT: { 42 projection->num_comp = diff_flux_proj->num_diff_flux_comps; 43 PetscCall(PetscObjectSetName((PetscObject)projection->dm, "DivDiffFluxProj")); 44 PetscCall(DMSetupByOrder_FEM(PETSC_TRUE, PETSC_TRUE, degree, 1, q_extra, 1, &projection->num_comp, projection->dm)); 45 46 PetscCall(DMPlexCeedElemRestrictionCreate(user->ceed, projection->dm, domain_label, label_value, height, dm_field, 47 &diff_flux_proj->elem_restr_div_diff_flux)); 48 PetscCallCeed(user->ceed, CeedElemRestrictionCreateVector(diff_flux_proj->elem_restr_div_diff_flux, &diff_flux_proj->div_diff_flux_ceed, NULL)); 49 PetscCall(CreateBasisFromPlex(user->ceed, projection->dm, domain_label, label_value, height, dm_field, &diff_flux_proj->basis_div_diff_flux)); 50 diff_flux_proj->eval_mode_div_diff_flux = CEED_EVAL_INTERP; 51 52 { // Create face labels on projection->dm for boundary integrals 53 DMLabel face_sets_label; 54 PetscInt num_face_set_values, *face_set_values; 55 56 PetscCall(DMGetLabel(user->dm, "Face Sets", &face_sets_label)); 57 PetscCall(DMLabelCreateGlobalValueArray(user->dm, face_sets_label, &num_face_set_values, &face_set_values)); 58 for (PetscInt f = 0; f < num_face_set_values; f++) { 59 DMLabel face_orientation_label; 60 char *face_orientation_label_name; 61 62 PetscCall(DMPlexCreateFaceLabel(user->dm, face_set_values[f], &face_orientation_label_name)); 63 PetscCall(DMGetLabel(user->dm, face_orientation_label_name, &face_orientation_label)); 64 PetscCall(DMAddLabel(projection->dm, face_orientation_label)); 65 PetscCall(PetscFree(face_orientation_label_name)); 66 } 67 PetscCall(PetscFree(face_set_values)); 68 } 69 } break; 70 case DIV_DIFF_FLUX_PROJ_INDIRECT: { 71 projection->num_comp = diff_flux_proj->num_diff_flux_comps * dim; 72 PetscCall(PetscObjectSetName((PetscObject)projection->dm, "DiffFluxProj")); 73 PetscCall(DMSetupByOrder_FEM(PETSC_TRUE, PETSC_TRUE, degree, 1, q_extra, 1, &projection->num_comp, projection->dm)); 74 75 PetscCall(DMPlexCeedElemRestrictionQDataCreate(user->ceed, projection->dm, domain_label, label_value, height, 76 diff_flux_proj->num_diff_flux_comps, &diff_flux_proj->elem_restr_div_diff_flux)); 77 PetscCallCeed(user->ceed, CeedElemRestrictionCreateVector(diff_flux_proj->elem_restr_div_diff_flux, &diff_flux_proj->div_diff_flux_ceed, NULL)); 78 diff_flux_proj->basis_div_diff_flux = CEED_BASIS_NONE; 79 diff_flux_proj->eval_mode_div_diff_flux = CEED_EVAL_NONE; 80 } break; 81 case DIV_DIFF_FLUX_PROJ_NONE: 82 SETERRQ(PetscObjectComm((PetscObject)user->dm), PETSC_ERR_ARG_WRONG, "Should not reach here with div_diff_flux_projection_method %s", 83 DivDiffFluxProjectionMethods[user->app_ctx->divFdiffproj_method]); 84 break; 85 } 86 PetscFunctionReturn(PETSC_SUCCESS); 87 }; 88 89 /** 90 @brief Setup direct projection of divergence of diffusive flux 91 92 @param[in] ceed `Ceed` context 93 @param[in] user `User` context 94 @param[in] ceed_data `CeedData` context 95 @param[in] problem `ProblemData` context 96 **/ 97 static PetscErrorCode DivDiffFluxProjectionSetup_Direct(Ceed ceed, User user, CeedData ceed_data, ProblemData problem) { 98 DivDiffFluxProjectionData diff_flux_proj = user->diff_flux_proj; 99 NodalProjectionData projection = diff_flux_proj->projection; 100 CeedOperator op_rhs; 101 CeedBasis basis_diff_flux; 102 CeedElemRestriction elem_restr_diff_flux_volume, elem_restr_qd; 103 CeedVector q_data; 104 CeedInt num_comp_q, q_data_size; 105 PetscInt dim, label_value = 0; 106 DMLabel domain_label = NULL; 107 MPI_Comm comm = PetscObjectComm((PetscObject)projection->dm); 108 109 PetscFunctionBeginUser; 110 // -- Get Pre-requisite things 111 PetscCall(DMGetDimension(projection->dm, &dim)); 112 PetscCallCeed(ceed, CeedBasisGetNumComponents(ceed_data->basis_q, &num_comp_q)); 113 114 { // Get elem_restr_diff_flux and basis_diff_flux 115 CeedOperator *sub_ops; 116 CeedOperatorField op_field; 117 PetscInt sub_op_index = 0; // will be 0 for the volume op 118 119 PetscCallCeed(ceed, CeedCompositeOperatorGetSubList(user->op_ifunction, &sub_ops)); 120 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "div F_diff", &op_field)); 121 PetscCallCeed(ceed, CeedOperatorFieldGetElemRestriction(op_field, &elem_restr_diff_flux_volume)); 122 PetscCallCeed(ceed, CeedOperatorFieldGetBasis(op_field, &basis_diff_flux)); 123 } 124 PetscCall(QDataGet(ceed, projection->dm, domain_label, label_value, ceed_data->elem_restr_x, ceed_data->basis_x, ceed_data->x_coord, &elem_restr_qd, 125 &q_data, &q_data_size)); 126 127 { // Create Projection RHS OperatorApplyContext 128 CeedOperator op_rhs; 129 130 PetscCheck(diff_flux_proj->CreateRHSOperator_Direct, comm, PETSC_ERR_ARG_WRONGSTATE, 131 "Must define CreateRHSOperator_Direct to use indirect div_diff_flux projection"); 132 PetscCall(diff_flux_proj->CreateRHSOperator_Direct(user, ceed_data, diff_flux_proj, &op_rhs)); 133 PetscCall(DMCreateLocalVector(projection->dm, &diff_flux_proj->DivDiffFlux_loc)); 134 diff_flux_proj->ceed_vec_has_array = PETSC_FALSE; 135 PetscCall(OperatorApplyContextCreate(user->dm, projection->dm, ceed, op_rhs, NULL, NULL, NULL, diff_flux_proj->DivDiffFlux_loc, 136 &projection->l2_rhs_ctx)); 137 PetscCallCeed(ceed, CeedOperatorDestroy(&op_rhs)); 138 } 139 140 { // -- Build Mass operator 141 CeedQFunction qf_mass; 142 CeedOperator op_mass; 143 144 PetscCall(CreateMassQFunction(ceed, projection->num_comp, q_data_size, &qf_mass)); 145 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_mass, NULL, NULL, &op_mass)); 146 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "u", elem_restr_diff_flux_volume, basis_diff_flux, CEED_VECTOR_ACTIVE)); 147 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "qdata", elem_restr_qd, CEED_BASIS_NONE, q_data)); 148 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "v", elem_restr_diff_flux_volume, basis_diff_flux, CEED_VECTOR_ACTIVE)); 149 150 { // -- Setup KSP for L^2 projection 151 Mat mat_mass; 152 153 PetscCall(MatCreateCeed(projection->dm, projection->dm, op_mass, NULL, &mat_mass)); 154 155 PetscCall(KSPCreate(comm, &projection->ksp)); 156 PetscCall(KSPSetOptionsPrefix(projection->ksp, "div_diff_flux_projection_")); 157 { // lumped by default 158 PC pc; 159 PetscCall(KSPGetPC(projection->ksp, &pc)); 160 PetscCall(PCSetType(pc, PCJACOBI)); 161 PetscCall(PCJacobiSetType(pc, PC_JACOBI_ROWSUM)); 162 PetscCall(KSPSetType(projection->ksp, KSPPREONLY)); 163 } 164 PetscCall(KSPSetFromOptions_WithMatCeed(projection->ksp, mat_mass)); 165 PetscCall(MatDestroy(&mat_mass)); 166 } 167 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_mass)); 168 PetscCallCeed(ceed, CeedOperatorDestroy(&op_mass)); 169 } 170 PetscCallCeed(ceed, CeedVectorDestroy(&q_data)); 171 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_qd)); 172 PetscCallCeed(ceed, CeedOperatorDestroy(&op_rhs)); 173 PetscFunctionReturn(PETSC_SUCCESS); 174 } 175 176 /** 177 @brief Setup indirect projection of divergence of diffusive flux 178 179 @param[in] ceed `Ceed` context 180 @param[in,out] user `User` context 181 @param[in] ceed_data `CeedData` context 182 @param[in] problem `ProblemData` context 183 **/ 184 static PetscErrorCode DivDiffFluxProjectionSetup_Indirect(Ceed ceed, User user, CeedData ceed_data, ProblemData problem) { 185 DivDiffFluxProjectionData diff_flux_proj = user->diff_flux_proj; 186 NodalProjectionData projection = diff_flux_proj->projection; 187 CeedBasis basis_diff_flux; 188 CeedElemRestriction elem_restr_diff_flux, elem_restr_qd; 189 CeedVector q_data; 190 CeedInt num_comp_q, q_data_size; 191 PetscInt dim; 192 PetscInt label_value = 0, height = 0, dm_field = 0; 193 DMLabel domain_label = NULL; 194 MPI_Comm comm = PetscObjectComm((PetscObject)projection->dm); 195 196 PetscFunctionBeginUser; 197 PetscCall(DMGetDimension(projection->dm, &dim)); 198 PetscCallCeed(ceed, CeedBasisGetNumComponents(ceed_data->basis_q, &num_comp_q)); 199 200 PetscCall(DMPlexCeedElemRestrictionCreate(ceed, projection->dm, domain_label, label_value, height, dm_field, &elem_restr_diff_flux)); 201 PetscCall(CreateBasisFromPlex(ceed, projection->dm, domain_label, label_value, height, dm_field, &basis_diff_flux)); 202 PetscCall(QDataGet(ceed, projection->dm, domain_label, label_value, ceed_data->elem_restr_x, ceed_data->basis_x, ceed_data->x_coord, &elem_restr_qd, 203 &q_data, &q_data_size)); 204 205 { 206 CeedOperator op_rhs; 207 208 PetscCheck(diff_flux_proj->CreateRHSOperator_Indirect, comm, PETSC_ERR_ARG_WRONGSTATE, 209 "Must define CreateRHSOperator_Indirect to use indirect div_diff_flux projection"); 210 PetscCall(diff_flux_proj->CreateRHSOperator_Indirect(user, ceed_data, diff_flux_proj, &op_rhs)); 211 PetscCall(OperatorApplyContextCreate(user->dm, projection->dm, ceed, op_rhs, NULL, NULL, NULL, NULL, &projection->l2_rhs_ctx)); 212 PetscCallCeed(ceed, CeedOperatorDestroy(&op_rhs)); 213 } 214 215 { // -- Build Mass operator 216 CeedQFunction qf_mass; 217 CeedOperator op_mass; 218 219 PetscCall(CreateMassQFunction(ceed, projection->num_comp, q_data_size, &qf_mass)); 220 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_mass, NULL, NULL, &op_mass)); 221 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "u", elem_restr_diff_flux, basis_diff_flux, CEED_VECTOR_ACTIVE)); 222 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "qdata", elem_restr_qd, CEED_BASIS_NONE, q_data)); 223 PetscCallCeed(ceed, CeedOperatorSetField(op_mass, "v", elem_restr_diff_flux, basis_diff_flux, CEED_VECTOR_ACTIVE)); 224 225 { // -- Setup KSP for L^2 projection 226 Mat mat_mass; 227 228 PetscCall(MatCreateCeed(projection->dm, projection->dm, op_mass, NULL, &mat_mass)); 229 230 PetscCall(KSPCreate(comm, &projection->ksp)); 231 PetscCall(KSPSetOptionsPrefix(projection->ksp, "div_diff_flux_projection_")); 232 { // lumped by default 233 PC pc; 234 PetscCall(KSPGetPC(projection->ksp, &pc)); 235 PetscCall(PCSetType(pc, PCJACOBI)); 236 PetscCall(PCJacobiSetType(pc, PC_JACOBI_ROWSUM)); 237 PetscCall(KSPSetType(projection->ksp, KSPPREONLY)); 238 } 239 PetscCall(KSPSetFromOptions_WithMatCeed(projection->ksp, mat_mass)); 240 PetscCall(MatDestroy(&mat_mass)); 241 } 242 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_mass)); 243 PetscCallCeed(ceed, CeedOperatorDestroy(&op_mass)); 244 } 245 246 { // Create OperatorApplyContext to calculate divergence at quadrature points 247 CeedQFunction qf_calc_divergence; 248 CeedOperator op_calc_divergence; 249 CeedElemRestriction elem_restr_div_diff_flux; 250 251 { // Get elem_restr_div_diff_flux 252 CeedOperator *sub_ops; 253 CeedOperatorField op_field; 254 PetscInt sub_op_index = 0; // will be 0 for the volume op 255 256 PetscCallCeed(ceed, CeedCompositeOperatorGetSubList(user->op_ifunction, &sub_ops)); 257 PetscCallCeed(ceed, CeedOperatorGetFieldByName(sub_ops[sub_op_index], "div F_diff", &op_field)); 258 PetscCallCeed(ceed, CeedOperatorFieldGetElemRestriction(op_field, &elem_restr_div_diff_flux)); 259 } 260 261 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, ComputeDivDiffusiveFlux3D_4, ComputeDivDiffusiveFlux3D_4_loc, &qf_calc_divergence)); 262 263 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_calc_divergence, "Grad F_diff", projection->num_comp * dim, CEED_EVAL_GRAD)); 264 PetscCallCeed(ceed, CeedQFunctionAddInput(qf_calc_divergence, "qdata", q_data_size, CEED_EVAL_NONE)); 265 PetscCallCeed(ceed, CeedQFunctionAddOutput(qf_calc_divergence, "Div F_diff", 4, CEED_EVAL_NONE)); 266 267 PetscCallCeed(ceed, CeedOperatorCreate(ceed, qf_calc_divergence, NULL, NULL, &op_calc_divergence)); 268 PetscCallCeed(ceed, CeedOperatorSetField(op_calc_divergence, "Grad F_diff", elem_restr_diff_flux, basis_diff_flux, CEED_VECTOR_ACTIVE)); 269 PetscCallCeed(ceed, CeedOperatorSetField(op_calc_divergence, "qdata", elem_restr_qd, CEED_BASIS_NONE, q_data)); 270 PetscCallCeed( 271 ceed, CeedOperatorSetField(op_calc_divergence, "Div F_diff", elem_restr_div_diff_flux, CEED_BASIS_NONE, diff_flux_proj->div_diff_flux_ceed)); 272 273 PetscCall( 274 OperatorApplyContextCreate(projection->dm, NULL, ceed, op_calc_divergence, NULL, NULL, NULL, NULL, &diff_flux_proj->calc_div_diff_flux)); 275 276 PetscCallCeed(ceed, CeedQFunctionDestroy(&qf_calc_divergence)); 277 PetscCallCeed(ceed, CeedOperatorDestroy(&op_calc_divergence)); 278 } 279 PetscCallCeed(ceed, CeedBasisDestroy(&basis_diff_flux)); 280 PetscCallCeed(ceed, CeedVectorDestroy(&q_data)); 281 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_qd)); 282 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&elem_restr_diff_flux)); 283 PetscFunctionReturn(PETSC_SUCCESS); 284 } 285 286 /** 287 @brief Setup projection of divergence of diffusive flux 288 289 @param[in] ceed `Ceed` context 290 @param[in,out] user `User` context 291 @param[in] ceed_data `CeedData` context 292 @param[in] problem `ProblemData` context 293 **/ 294 PetscErrorCode DivDiffFluxProjectionSetup(Ceed ceed, User user, CeedData ceed_data, ProblemData problem) { 295 PetscFunctionBeginUser; 296 switch (user->app_ctx->divFdiffproj_method) { 297 case DIV_DIFF_FLUX_PROJ_DIRECT: 298 PetscCall(DivDiffFluxProjectionSetup_Direct(ceed, user, ceed_data, problem)); 299 break; 300 case DIV_DIFF_FLUX_PROJ_INDIRECT: 301 PetscCall(DivDiffFluxProjectionSetup_Indirect(ceed, user, ceed_data, problem)); 302 break; 303 case DIV_DIFF_FLUX_PROJ_NONE: 304 SETERRQ(PetscObjectComm((PetscObject)user->dm), PETSC_ERR_ARG_WRONG, "Should not reach here with div_diff_flux_projection_method %s", 305 DivDiffFluxProjectionMethods[user->app_ctx->divFdiffproj_method]); 306 break; 307 } 308 PetscFunctionReturn(PETSC_SUCCESS); 309 } 310 311 /** 312 @brief Project the divergence of diffusive flux 313 314 This implicitly sets the `CeedVector` input (`div_diff_flux_ceed`) to the divergence of diffusive flux. 315 316 @param[in] diff_flux_proj `NodalProjectionData` for the projection 317 @param[in] Q_loc Localized solution vector 318 **/ 319 PetscErrorCode DivDiffFluxProjectionApply(DivDiffFluxProjectionData diff_flux_proj, Vec Q_loc) { 320 NodalProjectionData projection = diff_flux_proj->projection; 321 322 PetscFunctionBeginUser; 323 PetscCall(PetscLogEventBegin(FLUIDS_DivDiffFluxProjection, Q_loc, 0, 0, 0)); 324 switch (diff_flux_proj->method) { 325 case DIV_DIFF_FLUX_PROJ_DIRECT: { 326 Vec DivDiffFlux; 327 328 PetscCall(DMGetGlobalVector(projection->dm, &DivDiffFlux)); 329 if (diff_flux_proj->ceed_vec_has_array) { 330 PetscCall(VecReadCeedToPetsc(diff_flux_proj->div_diff_flux_ceed, diff_flux_proj->DivDiffFlux_memtype, diff_flux_proj->DivDiffFlux_loc)); 331 diff_flux_proj->ceed_vec_has_array = PETSC_FALSE; 332 } 333 PetscCall(ApplyCeedOperatorLocalToGlobal(Q_loc, DivDiffFlux, projection->l2_rhs_ctx)); 334 PetscCall(VecViewFromOptions(DivDiffFlux, NULL, "-div_diff_flux_projection_rhs_view")); 335 336 PetscCall(KSPSolve(projection->ksp, DivDiffFlux, DivDiffFlux)); 337 PetscCall(VecViewFromOptions(DivDiffFlux, NULL, "-div_diff_flux_projection_view")); 338 339 PetscCall(DMGlobalToLocal(projection->dm, DivDiffFlux, INSERT_VALUES, diff_flux_proj->DivDiffFlux_loc)); 340 PetscCall(VecReadPetscToCeed(diff_flux_proj->DivDiffFlux_loc, &diff_flux_proj->DivDiffFlux_memtype, diff_flux_proj->div_diff_flux_ceed)); 341 diff_flux_proj->ceed_vec_has_array = PETSC_TRUE; 342 343 PetscCall(DMRestoreGlobalVector(projection->dm, &DivDiffFlux)); 344 break; 345 } 346 case DIV_DIFF_FLUX_PROJ_INDIRECT: { 347 Vec DiffFlux; 348 349 PetscCall(DMGetGlobalVector(projection->dm, &DiffFlux)); 350 PetscCall(ApplyCeedOperatorLocalToGlobal(Q_loc, DiffFlux, projection->l2_rhs_ctx)); 351 PetscCall(VecViewFromOptions(DiffFlux, NULL, "-div_diff_flux_projection_rhs_view")); 352 353 PetscCall(KSPSolve(projection->ksp, DiffFlux, DiffFlux)); 354 PetscCall(VecViewFromOptions(DiffFlux, NULL, "-div_diff_flux_projection_view")); 355 356 PetscCall(ApplyCeedOperatorGlobalToLocal(DiffFlux, NULL, diff_flux_proj->calc_div_diff_flux)); 357 PetscCall(DMRestoreGlobalVector(projection->dm, &DiffFlux)); 358 } break; 359 case DIV_DIFF_FLUX_PROJ_NONE: 360 SETERRQ(PetscObjectComm((PetscObject)projection->dm), PETSC_ERR_ARG_WRONG, "Should not reach here with div_diff_flux_projection_method %s", 361 DivDiffFluxProjectionMethods[diff_flux_proj->method]); 362 break; 363 } 364 PetscCall(PetscLogEventEnd(FLUIDS_DivDiffFluxProjection, Q_loc, 0, 0, 0)); 365 PetscFunctionReturn(PETSC_SUCCESS); 366 } 367 368 /** 369 @brief Destroy `DivDiffFluxProjectionData` object 370 371 @param[in,out] diff_flux_proj Object to destroy 372 **/ 373 PetscErrorCode DivDiffFluxProjectionDataDestroy(DivDiffFluxProjectionData diff_flux_proj) { 374 PetscFunctionBeginUser; 375 if (diff_flux_proj == NULL) PetscFunctionReturn(PETSC_SUCCESS); 376 PetscCall(NodalProjectionDataDestroy(diff_flux_proj->projection)); 377 PetscCall(OperatorApplyContextDestroy(diff_flux_proj->calc_div_diff_flux)); 378 if (diff_flux_proj->ceed_vec_has_array) { 379 PetscCall(VecReadCeedToPetsc(diff_flux_proj->div_diff_flux_ceed, diff_flux_proj->DivDiffFlux_memtype, diff_flux_proj->DivDiffFlux_loc)); 380 diff_flux_proj->ceed_vec_has_array = PETSC_FALSE; 381 } 382 PetscCall(CeedVectorDestroy(&diff_flux_proj->div_diff_flux_ceed)); 383 PetscCall(VecDestroy(&diff_flux_proj->DivDiffFlux_loc)); 384 PetscCall(PetscFree(diff_flux_proj)); 385 PetscFunctionReturn(PETSC_SUCCESS); 386 } 387