From 7f621ea751e208a9ee5cedf168b4aeb7eec66ee8 Mon Sep 17 00:00:00 2001 From: Amishga Alphonius <107414376+AmishgaAlphonius@users.noreply.github.com> Date: Mon, 22 Jan 2024 09:15:35 -0500 Subject: [PATCH] Fix heat flux postprocessor material_id processing (#988) Description of the problem In the previous implementation of the heat flux postprocessor, there was a confusion related to the material_id of a cell. As a cell's material_id is defined in the mesh, fluids moving from a cell to another did not have the right material_id when outputted. Futhermore, some dealii meshs have an intrinc way of numbering material_ids (e.g. colorized subdivided_hyper_rectangle) that were leading to inadequate comparison of the material_ids. Description of the solution As of this PR, in the heat flux postprocessor, the material_id is only used to indentify cell material when solids are present. When simulating with fluids and solids, cells of the mesh in which the fluids lie should have a material_id of 0. For fluids-only simulations, at the moment, their heat fluxes will be output on the entire domain, as the current DataPostprocessor can only hold 1 dof_handler. During postprocessing, the user could clip the domain of interest of each fluid using the phase or filtered_phase fields outputted. How Has This Been Tested? Visual inspection of outputted results for a: single fluid simulation, 1 fluid + 1 solid simulation, and 2 fluids simulation. Future changes For multiple-fluids simulations, a DataPostprocessor handling multiple dof_handlers to allow postprocessing of quantities of interest of fluids in their respective subdomains will have to be added to the code. Former-commit-id: ea4aeae4dd5babf8108d5f994b855d95ea077de2 --- include/solvers/postprocessors.h | 108 +++++++++++++++++++++---------- source/solvers/heat_transfer.cc | 21 ++++-- 2 files changed, 92 insertions(+), 37 deletions(-) diff --git a/include/solvers/postprocessors.h b/include/solvers/postprocessors.h index dde4e2b84c..f21eed3a98 100644 --- a/include/solvers/postprocessors.h +++ b/include/solvers/postprocessors.h @@ -511,16 +511,20 @@ class DensityPostprocessor : public DataPostprocessorScalar /** - * @class Calculates heat flux (q = -k * grad T) of a material at every - * quadrature point. + * @class Postprocessor that calculates local heat fluxes \f$(q = -k * \nabla + * T)\f$ of a material at every evaluation point of the domain. + * * @param p_thermal_conductivity_model Thermal conductivity model of the * material * @param p_material_string Type of material: "f" (fluid) or "s" (solid) * @param p_id Either fluid ID or solid ID (as initialized in the physical * properties) used to name the postprocessed variable vector. - * @param p_material_id ID corresponding to the material. Material IDs the - * cumulative IDs of fluids and solids by convention, fluids are counted before - * solids. This is used to identify which cell lies in which material. + * @param p_mesh_material_id ID corresponding to the material as identified in + * the mesh. By convention, fluids are given the material_id 0 and solids follow + * with ids 1, 2, etc. This is used to identify which cell lies in which + * material. + * @param p_solid_phase_present Boolean indicating if the simulation has solids. + * If there is no solid, the @p material_id of cells is not checked. */ template class HeatFluxPostprocessor : public DataPostprocessorVector @@ -528,56 +532,94 @@ class HeatFluxPostprocessor : public DataPostprocessorVector public: HeatFluxPostprocessor(const std::shared_ptr &p_thermal_conductivity_model, - const std::string &p_material_string = "f", - const unsigned int &p_id = 0, - const unsigned int &p_material_id = 0) + const std::string &p_material_string = "f", + const unsigned int &p_id = 0, + const unsigned int &p_mesh_material_id = 0, + const bool &p_solid_phase_present = false) : DataPostprocessorVector("heat_flux_" + p_material_string + Utilities::to_string(p_id, 2), update_values | update_gradients) , thermal_conductivity_model(p_thermal_conductivity_model) - , material_id(p_material_id) + , mesh_material_id(p_mesh_material_id) + , solid_phase_present(p_solid_phase_present) {} + /** + * @brief Computes local heat flux quantities from temperature field and + * material properties. + * + * @param inputs Structure that contains all temperature information needed to + * evaluate local heat fluxes. + * @param computed_quantities Vector field containing local heat flux vectors + * evaluated at evaluation points. + */ virtual void evaluate_scalar_field( const DataPostprocessorInputs::Scalar &inputs, std::vector> &computed_quantities) const override { - const unsigned int n_quadrature_points = computed_quantities.size(); - AssertDimension(inputs.solution_gradients.size(), n_quadrature_points); - - // Get Cell - const auto cell = inputs.template get_cell(); + const unsigned int n_evaluation_points = computed_quantities.size(); + AssertDimension(inputs.solution_gradients.size(), n_evaluation_points); - // Check if the cell lies in the material of interest - if (cell->material_id() == material_id) + // If solid materials are present, identify each material's subdomain. + if (solid_phase_present) { - std::map field_values; - for (unsigned int q = 0; q < n_quadrature_points; ++q) - { - AssertDimension(computed_quantities[q].size(), dim); - field_values[field::temperature] = inputs.solution_values[q]; + // Get Cell + const auto cell = inputs.template get_cell(); - for (unsigned int i = 0; i < dim; ++i) - { - computed_quantities[q][i] = - -thermal_conductivity_model->value(field_values) * - inputs.solution_gradients[q][i]; - } + // Check if the cell lies in the material of interest + if (cell->material_id() == mesh_material_id) + { + compute_quantities(inputs, computed_quantities); + } + // Apply null values to irrelevant subdomain + else + { + Vector null_vector(dim); + std::vector> null_computed_quantities( + n_evaluation_points, null_vector); + computed_quantities = null_computed_quantities; } } - // Apply null values to irrelevant subdomain + // Only fluids are present else { - Vector null_vector(dim); - std::vector> null_computed_quantities( - n_quadrature_points, null_vector); - computed_quantities = null_computed_quantities; + compute_quantities(inputs, computed_quantities); } } private: + /** + * @brief Computes local heat flux quantities. This function is used to + * improve readability of the HeatFluxPostprocessor::evaluate_scalar_field + * function and avoid repeated lines. + * + * @param inputs Structure that contains all temperature information needed to + * evaluate local heat fluxes. + * @param computed_quantities Vector field containing local heat flux vectors + * evaluated at evaluation points. + */ + void + compute_quantities(const DataPostprocessorInputs::Scalar &inputs, + std::vector> &computed_quantities) const + { + std::map field_values; + for (unsigned int p = 0; p < computed_quantities.size(); ++p) + { + AssertDimension(computed_quantities[p].size(), dim); + field_values[field::temperature] = inputs.solution_values[p]; + + for (unsigned int i = 0; i < dim; ++i) + { + computed_quantities[p][i] = + -thermal_conductivity_model->value(field_values) * + inputs.solution_gradients[p][i]; + } + } + } + std::shared_ptr thermal_conductivity_model; - unsigned int material_id; + const unsigned int mesh_material_id; + const bool solid_phase_present; }; #endif diff --git a/source/solvers/heat_transfer.cc b/source/solvers/heat_transfer.cc index a41cbdbf7a..51ad829c78 100644 --- a/source/solvers/heat_transfer.cc +++ b/source/solvers/heat_transfer.cc @@ -696,6 +696,10 @@ HeatTransfer::attach_solution_to_output(DataOut &data_out) const unsigned int n_solids = this->simulation_parameters.physical_properties_manager .get_number_of_solids(); + unsigned int mesh_m_id = 0; + + // Check if there are solids + const bool solid_phase_present = (n_solids > 0); // Postprocess heat fluxes heat_flux_postprocessors.clear(); @@ -703,8 +707,12 @@ HeatTransfer::attach_solution_to_output(DataOut &data_out) // Heat fluxes in fluids for (unsigned int f_id = 0; f_id < n_fluids; ++f_id) { - heat_flux_postprocessors.push_back(HeatFluxPostprocessor( - thermal_conductivity_models[f_id], "f", f_id, f_id)); + heat_flux_postprocessors.push_back( + HeatFluxPostprocessor(thermal_conductivity_models[f_id], + "f", + f_id, + mesh_m_id, + solid_phase_present)); data_out.add_data_vector(this->dof_handler, this->present_solution, heat_flux_postprocessors[f_id]); @@ -712,8 +720,13 @@ HeatTransfer::attach_solution_to_output(DataOut &data_out) // Heat fluxes in solids for (unsigned int m_id = n_fluids; m_id < n_fluids + n_solids; ++m_id) { - heat_flux_postprocessors.push_back(HeatFluxPostprocessor( - thermal_conductivity_models[m_id], "s", m_id - n_fluids, m_id)); + mesh_m_id += 1; + heat_flux_postprocessors.push_back( + HeatFluxPostprocessor(thermal_conductivity_models[m_id], + "s", + m_id - n_fluids, + mesh_m_id, + solid_phase_present)); data_out.add_data_vector(this->dof_handler, this->present_solution, heat_flux_postprocessors[m_id]);