1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 #include <navierstokes.h> 5 6 PetscClassId HONEE_CLASSID; 7 8 // @brief Initalize `HONEE` class. 9 static PetscErrorCode HoneeInitClass_Private() { 10 static PetscBool registered = PETSC_FALSE; 11 12 PetscFunctionBeginUser; 13 if (registered) PetscFunctionReturn(PETSC_SUCCESS); 14 PetscCall(PetscClassIdRegister("Honee", &HONEE_CLASSID)); 15 registered = PETSC_TRUE; 16 PetscFunctionReturn(PETSC_SUCCESS); 17 } 18 19 /** 20 @brief Initialize `Honee` object 21 22 @param[in] comm `MPI_Comm` for the object 23 @param[out] honee The initialized `Honee` object 24 **/ 25 PetscErrorCode HoneeInit(MPI_Comm comm, Honee *honee) { 26 Honee honee_; 27 AppCtx app_ctx; 28 ProblemData problem; 29 Physics phys_ctx; 30 Units units; 31 32 PetscFunctionBeginUser; 33 PetscCall(HoneeInitClass_Private()); 34 PetscCall(PetscHeaderCreate(honee_, HONEE_CLASSID, "Honee", "HONEE", "Honee", comm, HoneeDestroy, NULL)); 35 PetscCall(PetscNew(&app_ctx)); 36 PetscCall(PetscNew(&problem)); 37 PetscCall(PetscNew(&phys_ctx)); 38 PetscCall(PetscNew(&units)); 39 honee_->app_ctx = app_ctx; 40 honee_->units = units; 41 honee_->phys = phys_ctx; 42 honee_->problem_data = problem; 43 problem->set_bc_from_ics = PETSC_TRUE; 44 honee_->comm = PETSC_COMM_WORLD; 45 honee_->start_time = time(NULL); 46 47 PetscCall(RegisterLogEvents()); 48 49 *honee = honee_; 50 PetscFunctionReturn(PETSC_SUCCESS); 51 } 52 53 /** 54 @brief Destroy `Honee` object 55 56 @param[in] honee Object to be destroyed 57 **/ 58 PetscErrorCode HoneeDestroy(Honee *honee) { 59 Honee honee_ = *honee; 60 Ceed ceed; 61 62 PetscFunctionBeginUser; 63 if (!honee_) PetscFunctionReturn(PETSC_SUCCESS); 64 PetscValidHeaderSpecific(honee_, HONEE_CLASSID, 1); 65 66 ceed = honee_->ceed; 67 MPI_Comm comm = honee_->comm; 68 AppCtx app_ctx = honee_->app_ctx; 69 ProblemData problem = honee_->problem_data; 70 71 PetscCall(NodalProjectionDataDestroy(honee_->grad_velo_proj)); 72 PetscCall(SgsDDDataDestroy(honee_->sgs_dd_data)); 73 PetscCall(DifferentialFilterDataDestroy(honee_->diff_filter)); 74 PetscCall(SGS_DD_TrainingDataDestroy(honee_->sgs_dd_train)); 75 PetscCall(SmartSimDataDestroy(honee_->smartsim)); 76 PetscCall(QDataClearStoredData()); 77 PetscCall(DivDiffFluxProjectionDataDestroy(honee_->diff_flux_proj)); 78 79 // -- Vectors 80 PetscCallCeed(ceed, CeedVectorDestroy(&honee_->x_coord)); 81 PetscCallCeed(ceed, CeedVectorDestroy(&honee_->q_ceed)); 82 PetscCallCeed(ceed, CeedVectorDestroy(&honee_->q_dot_ceed)); 83 PetscCallCeed(ceed, CeedVectorDestroy(&honee_->g_ceed)); 84 85 // -- Bases 86 PetscCallCeed(ceed, CeedBasisDestroy(&honee_->basis_q)); 87 PetscCallCeed(ceed, CeedBasisDestroy(&honee_->basis_x)); 88 89 // -- Restrictions 90 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&honee_->elem_restr_q)); 91 PetscCallCeed(ceed, CeedElemRestrictionDestroy(&honee_->elem_restr_x)); 92 93 // Destroy QFunction contexts after using 94 { 95 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->ics.qfctx)); 96 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_rhs.qfctx)); 97 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_ifunction.qfctx)); 98 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->apply_vol_ijacobian.qfctx)); 99 } 100 101 // -- Operators 102 PetscCall(OperatorApplyContextDestroy(honee_->op_ics_ctx)); 103 PetscCall(OperatorApplyContextDestroy(honee_->op_rhs_ctx)); 104 PetscCall(OperatorApplyContextDestroy(honee_->op_strong_bc_ctx)); 105 PetscCallCeed(ceed, CeedOperatorDestroy(&honee_->op_ifunction)); 106 107 // -- Ceed 108 PetscCheck(CeedDestroy(&ceed) == CEED_ERROR_SUCCESS, comm, PETSC_ERR_LIB, "Destroying Ceed object failed"); 109 110 // --------------------------------------------------------------------------- 111 // Clean up PETSc 112 // --------------------------------------------------------------------------- 113 // -- Vectors 114 PetscCall(VecDestroy(&honee_->Q_loc)); 115 PetscCall(VecDestroy(&honee_->Q_dot_loc)); 116 117 PetscCall(KSPDestroy(&honee_->mass_ksp)); 118 119 // -- Matrices 120 PetscCall(MatDestroy(&honee_->interp_viz)); 121 PetscCall(MatDestroy(&honee_->mat_ijacobian)); 122 123 // -- DM 124 PetscCall(DMDestroy(&honee_->dm_viz)); 125 126 // -- Function list 127 PetscCall(PetscFunctionListDestroy(&app_ctx->problems)); 128 129 PetscCall(PetscFree(app_ctx->amat_type)); 130 PetscCall(PetscFree(app_ctx->wall_forces.walls)); 131 PetscCall(PetscViewerDestroy(&app_ctx->wall_forces.viewer)); 132 133 // -- Structs 134 for (PetscInt i = 0; i < problem->num_bc_defs; i++) { 135 PetscCall(BCDefinitionDestroy(&problem->bc_defs[i])); 136 } 137 PetscCall(PetscFree(problem->bc_defs)); 138 PetscCall(PetscFree(honee_->units)); 139 PetscCall(PetscFree(problem)); 140 PetscCall(PetscFree(honee_->phys)); 141 PetscCall(PetscFree(app_ctx)); 142 PetscCall(PetscHeaderDestroy(&honee_)); 143 *honee = NULL; 144 PetscFunctionReturn(PETSC_SUCCESS); 145 } 146