24 #include "petscviewer.h" 25 #include "petscerror.h" 43 *
this += water_content.name(
"water_content")
46 .description(R
"(Water content. 47 It is a fraction of water volume to the whole volume.)"); 48 *this += conductivity_richards.name(
"conductivity_richards")
49 .units(
UnitSI().m().s(-1) )
51 .description(
"Computed isotropic scalar conductivity by the soil model.");
53 *
this += water_content_saturated.name(
"water_content_saturated")
54 .description(R
"(Saturated water content (($ \theta_s $)). 55 Relative volume of water in a reference volume of a saturated porous media.)") 59 *
this += water_content_residual.name(
"water_content_residual")
60 .description(R
"(Residual water content (($ \theta_r $)). 61 Relative volume of water in a reference volume of an ideally dry porous media.)") 65 *
this += genuchten_p_head_scale.name(
"genuchten_p_head_scale")
66 .description(R
"(The van Genuchten pressure head scaling parameter (($ \alpha $)). 67 It is related to the inverse of the air entry pressure, i.e. the pressure 68 where the relative water content starts to decrease below 1.)") 72 *
this += genuchten_n_exponent.name(
"genuchten_n_exponent")
73 .description(
"The van Genuchten exponent parameter (($ n $)).")
90 auto soil_rec =
it::Record(
"SoilModel",
"Soil model settings.")
93 "Selection of the globally applied soil model. In future we replace this key by a field for selection of the model." 94 "That will allow usage of different soil model in a single simulation.")
96 "Fraction of the water content where we cut and rescale the curve.")
101 return it::Record(
"Flow_Richards_LMH",
"Lumped Mixed-Hybrid solver for unsteady unsaturated Darcy flow.")
105 "Input data for Darcy flow model.")
107 IT::Default(
"{ \"fields\": [ \"pressure_p0\", \"velocity_p0\" ] }"),
108 "Specification of output fields and output times.")
110 "Soil model settings.")
116 Input::register_class< RichardsLMH, Mesh &, const Input::Record >(
"Flow_Richards_LMH") +
124 data_ = make_shared<EqData>();
135 double fraction= model_rec.val<
double>(
"cut_fraction");
137 data_->soil_model_ = std::make_shared<SoilModel_VanGenuchten>(fraction);
139 data_->soil_model_ = std::make_shared<SoilModel_Irmay>(fraction);
144 data_->water_content_previous_time =
data_->dh_cr_disc_->create_vector();
145 data_->capacity =
data_->dh_cr_disc_->create_vector();
151 data_->water_content_ptr = create_field_fe< 3, FieldValue<3>::Scalar >(
data_->dh_cr_disc_);
154 data_->conductivity_ptr = create_field_fe< 3, FieldValue<3>::Scalar >(
data_->dh_p_);
158 data_->multidim_assembler = AssemblyBase::create< AssemblyRichards >(
data_);
168 data_->multidim_assembler[dh_cell.elm().dim()-1]->update_water_content(dh_cell);
175 data_->p_edge_solution_previous_time.copy_from(
data_->p_edge_solution);
177 data_->water_content_previous_time.copy_from(water_content_vec);
179 data_->p_edge_solution_previous_time.local_to_ghost_begin();
180 data_->p_edge_solution_previous_time.local_to_ghost_end();
185 return (
data_->storativity.input_list_size() == 0)
186 && (
data_->water_content_saturated.input_list_size() == 0);
201 data_->p_edge_solution.local_to_ghost_begin();
202 data_->p_edge_solution.local_to_ghost_end();
void initialize_specific() override
Output class for darcy_flow_mh model.
RegionSet get_region_set(const std::string &set_name) const
void initial_condition_postprocess() override
static const Input::Type::Record & get_input_type()
virtual void start_add_assembly()
void assembly_mh_matrix(MultidimAssembly &assembler)
virtual PetscErrorCode mat_zero_entries()
bool zero_time_term(bool time_global=false) override
Cell accessor allow iterate over DOF handler cells.
virtual void finish_assembly()=0
static const Input::Type::Instance & get_input_type(FieldSet &eq_data, const std::string &equation_name)
const RegionDB & region_db() const
#define ASSERT(expr)
Allow use shorter versions of macro names if these names is not used with external library...
Mixed-hybrid model of linear Darcy flow, possibly unsteady.
static const std::string field_descriptor_record_description(const string &record_name)
Basic time management functionality for unsteady (and steady) solvers (class Equation).
Basic time management class.
Assembly explicit Schur complement for the given linear system. Provides method for resolution of the...
Global macros to enhance readability and debugging, general constants.
static constexpr Mask equation_result
Match result fields. These are never given by input or copy of input.
Edge lumped mixed-hybrid solution of unsteady Darcy flow.
#define START_TIMER(tag)
Starts a timer with specified tag.
static Input::Type::Abstract & get_input_type()
std::shared_ptr< EqData > data_
virtual PetscErrorCode rhs_zero_entries()
std::shared_ptr< EqData > data_
LinSys & lin_sys_schur()
Getter for the linear system of the 2. Schur complement.
#define ASSERT_PTR(ptr)
Definition of assert macro checking non-null pointer (PTR)
Input::Record input_record_
static const Input::Type::Record & get_input_type()
RichardsLMH(Mesh &mesh, const Input::Record in_rec, TimeGovernor *tm=nullptr)
void assembly_linear_system() override
void set_matrix_changed()
static const Input::Type::Record & type_field_descriptor()
void accept_time_step() override
postprocess velocity field (add sources)
static const int registrar
Registrar of class to factory.
Class for representation SI units of Fields.
static UnitSI & dimensionless()
Returns dimensionless unit.
#define FLOW123D_FORCE_LINK_IN_CHILD(x)