OGS
CreateLiquidFlowProcess.cpp
Go to the documentation of this file.
1
13
14#include <algorithm>
15#include <typeinfo>
16
17#include "LiquidFlowProcess.h"
24#include "ParameterLib/Utils.h"
27
28namespace ProcessLib
29{
30namespace LiquidFlow
31{
33 MeshLib::Mesh const& mesh,
35 bool const is_equation_type_volume)
36{
38
39 std::array const required_medium_properties = {
44
45 std::array const required_liquid_properties = {
48
49 // Check Constant-type density.
50 if (is_equation_type_volume)
51 {
52 for (auto const& medium : media_map.media())
53 {
54 // auto const& medium = *media_map.getMedium(element_id);
55 auto const& fluid_phase_density =
57 if (typeid(fluid_phase_density) !=
59 {
61 "Since `equation_balance_type` is set to `volume`,the "
62 "phase density type must be `Constant`. Note: by "
63 "default, `equation_balance_type` is set to `volume`.");
64 }
65 }
66 }
67
68 for (auto const& medium : media_map.media())
69 {
70 checkRequiredProperties(*medium, required_medium_properties);
71 checkRequiredProperties(fluidPhase(*medium),
72 required_liquid_properties);
73 }
74 DBUG("Media properties verified.");
75}
76
78 std::string_view const type_in_string)
79{
80 if (type_in_string == "volume")
81 {
83 }
84 if (type_in_string == "mass")
85 {
87 }
89 "The value of `equation_balance_type` must be either `volume` or "
90 "`mass`");
91};
92
93std::unique_ptr<Process> createLiquidFlowProcess(
94 std::string const& name,
95 MeshLib::Mesh& mesh,
96 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
97 std::vector<ProcessVariable> const& variables,
98 std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
99 unsigned const integration_order,
100 BaseLib::ConfigTree const& config,
101 std::vector<std::unique_ptr<MeshLib::Mesh>> const& meshes,
102 std::map<int, std::shared_ptr<MaterialPropertyLib::Medium>> const& media)
103{
105 config.checkConfigParameter("type", "LIQUID_FLOW");
106
107 DBUG("Create LiquidFlowProcess.");
108
110
112 auto const pv_config = config.getConfigSubtree("process_variables");
113
115 auto per_process_variables = findProcessVariables(
116 variables, pv_config,
117 {
118 "process_variable"});
119 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>
120 process_variables;
121 process_variables.push_back(std::move(per_process_variables));
122
123 SecondaryVariableCollection secondary_variables;
124
125 ProcessLib::createSecondaryVariables(config, secondary_variables);
126
128 std::vector<double> const b =
130 config.getConfigParameter<std::vector<double>>("specific_body_force");
131 int const mesh_space_dimension =
133 if (static_cast<int>(b.size()) != mesh_space_dimension)
134 {
135 OGS_FATAL(
136 "specific body force (gravity vector) has {:d} components, but the "
137 "space dimension is {:d}.",
138 b.size(), mesh_space_dimension);
139 }
140
141 Eigen::VectorXd specific_body_force(b.size());
142 bool const has_gravity = MathLib::toVector(b).norm() > 0;
143 if (has_gravity)
144 {
145 std::copy_n(b.data(), b.size(), specific_body_force.data());
146 }
147
148 std::unique_ptr<ProcessLib::SurfaceFluxData> surfaceflux;
149 auto calculatesurfaceflux_config =
151 config.getConfigSubtreeOptional("calculatesurfaceflux");
152 if (calculatesurfaceflux_config)
153 {
155 *calculatesurfaceflux_config, meshes);
156 }
157
158 auto media_map =
160
161 auto const is_linear =
163 config.getConfigParameter<bool>("linear", false);
164 if (is_linear)
165 {
166 INFO("LiquidFlow process is set to be linear.");
167 }
168
169 auto const equation_balance_type_str =
171 config.getConfigParameter<std::string>("equation_balance_type",
172 "volume");
173
174 DBUG("Check the media properties of LiquidFlow process ...");
175 checkMPLProperties(mesh, media_map, equation_balance_type_str == "volume");
176 DBUG("Media properties verified.");
177
178 auto const* aperture_size_parameter = &ParameterLib::findParameter<double>(
180
181 auto const aperture_config =
183 config.getConfigSubtreeOptional("aperture_size");
184 if (aperture_config)
185 {
186 aperture_size_parameter = &ParameterLib::findParameter<double>(
188 *aperture_config, "parameter", parameters, 1);
189 }
190
191 // The uniqueness of phase has already been checked in
192 // `checkMPLProperties`.
193 MaterialPropertyLib::Variable const phase_variable =
194 (*ranges::begin(media_map.media()))->hasPhase("Gas")
197
198 LiquidFlowData process_data{
199 covertEquationBalanceTypeFromString(equation_balance_type_str),
200 std::move(media_map),
202 mesh_space_dimension, mesh.getDimension(), mesh.getElements()),
203 mesh_space_dimension,
204 std::move(specific_body_force),
205 has_gravity,
206 *aperture_size_parameter,
207 NumLib::ShapeMatrixCache{integration_order},
208 phase_variable};
209
210 return std::make_unique<LiquidFlowProcess>(
211 std::move(name), mesh, std::move(jacobian_assembler), parameters,
212 integration_order, std::move(process_variables),
213 std::move(process_data), std::move(secondary_variables),
214 std::move(surfaceflux), is_linear);
215}
216
217} // namespace LiquidFlow
218} // namespace ProcessLib
#define OGS_FATAL(...)
Definition Error.h:26
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
void DBUG(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:30
std::optional< ConfigTree > getConfigSubtreeOptional(std::string const &root) const
T getConfigParameter(std::string const &param) const
ConfigTree getConfigSubtree(std::string const &root) const
void checkConfigParameter(std::string const &param, std::string_view const value) const
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:106
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:109
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
Definition Mesh.h:88
static PROCESSLIB_EXPORT const std::string constant_one_parameter_name
Definition Process.h:47
Handles configuration of several secondary variables from the project file.
void checkMPLPhasesForSinglePhaseFlow(MeshLib::Mesh const &mesh, MaterialPropertyLib::MaterialSpatialDistributionMap const &media_map)
MaterialSpatialDistributionMap createMaterialSpatialDistributionMap(std::map< int, std::shared_ptr< Medium > > const &media, MeshLib::Mesh const &mesh)
Eigen::Map< const Vector > toVector(std::vector< double > const &data, Eigen::VectorXd::Index size)
Creates an Eigen mapped vector from the given data vector.
std::vector< Eigen::MatrixXd > getElementRotationMatrices(int const space_dimension, int const mesh_dimension, std::vector< Element * > const &elements)
Element rotation matrix computation.
int getSpaceDimension(std::vector< Node * > const &nodes)
Computes dimension of the embedding space containing the set of given points.
void checkMPLProperties(MeshLib::Mesh const &mesh, MaterialPropertyLib::MaterialSpatialDistributionMap const &media_map, bool const is_equation_type_volume)
EquationBalanceType
Governing equation balance type.
std::unique_ptr< Process > createLiquidFlowProcess(std::string const &name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< ProcessVariable > const &variables, std::vector< std::unique_ptr< ParameterLib::ParameterBase > > const &parameters, unsigned const integration_order, BaseLib::ConfigTree const &config, std::vector< std::unique_ptr< MeshLib::Mesh > > const &meshes, std::map< int, std::shared_ptr< MaterialPropertyLib::Medium > > const &media)
EquationBalanceType covertEquationBalanceTypeFromString(std::string_view const type_in_string)
std::vector< std::reference_wrapper< ProcessVariable > > findProcessVariables(std::vector< ProcessVariable > const &variables, BaseLib::ConfigTree const &pv_config, std::initializer_list< std::string > tags)
void createSecondaryVariables(BaseLib::ConfigTree const &config, SecondaryVariableCollection &secondary_variables)
static std::unique_ptr< ProcessLib::SurfaceFluxData > createSurfaceFluxData(BaseLib::ConfigTree const &calculatesurfaceflux_config, std::vector< std::unique_ptr< MeshLib::Mesh > > const &meshes)