22 #include "petscviewer.h" 23 #include "petscerror.h" 44 *
this += water_content_saturated.name(
"water_content_saturated")
45 .description(R
"(Saturated water content (($ \theta_s $)). 46 Relative volume of water in a reference volume of a saturated porous media.)") 50 *
this += water_content_residual.name(
"water_content_residual")
51 .description(R
"(Residual water content (($ \theta_r $)). 52 Relative volume of water in a reference volume of an ideally dry porous media.)") 56 *
this += genuchten_p_head_scale.name(
"genuchten_p_head_scale")
57 .description(R
"(The van Genuchten pressure head scaling parameter (($ \alpha $)). 58 It is related to the inverse of the air entry pressure, i.e. the pressure 59 where the relative water content starts to decrease below 1.)") 63 *
this += genuchten_n_exponent.name(
"genuchten_n_exponent")
64 .description(
"The van Genuchten exponent parameter (($ n $)).")
81 auto soil_rec =
it::Record(
"SoilModel",
"Soil model settings.")
84 "Selection of the globally applied soil model. In future we replace this key by a field for selection of the model." 85 "That will allow usage of different soil model in a single simulation.")
87 "Fraction of the water content where we cut and rescale the curve.")
90 return it::Record(
"Flow_Richards_LMH",
"Lumped Mixed-Hybrid solver for unsteady unsaturated Darcy flow.")
94 "Input data for Darcy flow model.")
96 "Soil model settings.")
102 Input::register_class< RichardsLMH, Mesh &, const Input::Record >(
"Flow_Richards_LMH") +
110 data_ = make_shared<EqData>();
120 double fraction= model_rec.val<
double>(
"cut_fraction");
122 data_->soil_model_ = std::make_shared<SoilModel_VanGenuchten>(fraction);
124 data_->soil_model_ = std::make_shared<SoilModel_Irmay>(fraction);
131 data_->phead_edge_.resize( n_local_edges);
132 data_->water_content_previous_it.resize(n_local_sides);
133 data_->water_content_previous_time.resize(n_local_sides);
134 data_->capacity.resize(n_local_sides);
136 Distribution ds_split_edges(n_local_edges, PETSC_COMM_WORLD);
143 ISCreateGeneral(PETSC_COMM_SELF, local_edge_rows.size(),
144 &(local_edge_rows[0]), PETSC_COPY_VALUES, &(is_loc));
148 chkerr(ISDestroy(&is_loc));
165 for (
unsigned int i_loc_el = 0; i_loc_el <
mh_dh.
el_ds->
lsize(); i_loc_el++) {
168 init_value =
data_->init_pressure.value(ele_ac.centre(), ele_ac.element_accessor());
170 for (
unsigned int i=0; i<ele_ac.element_accessor()->n_sides(); i++) {
171 int edge_row = ele_ac.edge_row(i);
172 uint n_sides_of_edge = ele_ac.element_accessor().side(i)->edge().n_sides();
198 data_->water_content_previous_time.copy(
data_->water_content_previous_it);
204 return (
data_->storativity.input_list_size() == 0)
205 && (
data_->water_content_saturated.input_list_size() == 0);
233 auto multidim_assembler = AssemblyBase::create< AssemblyLMH >(
data_);
290 auto multidim_assembler = AssemblyBase::create< AssemblyLMH >(
data_);
295 multidim_assembler[ele_ac.dim()-1]->update_water_content(ele_ac);
297 double ele_scale = ele_ac.
measure() *
298 data_->cross_section.value(ele_ac.centre(), ele_ac.element_accessor()) / ele_ac.n_sides();
299 double ele_source =
data_->water_source_density.value(ele_ac.centre(), ele_ac.element_accessor())
300 +
data_->extra_source.value(ele_ac.centre(), ele_ac.element_accessor());
303 for (
unsigned int i=0; i<ele_ac.element_accessor()->n_sides(); i++) {
305 side_rows[i] = ele_ac.side_row(i);
306 double water_content =
data_->water_content_previous_it[ele_ac.side_local_idx(i)];
307 double water_content_previous_time =
data_->water_content_previous_time[ele_ac.side_local_idx(i)];
309 values[i] = ele_scale * ele_source - ele_scale * (water_content - water_content_previous_time) /
time_->
dt();
311 VecSetValues(
schur0->
get_solution(), ele_ac.n_sides(), side_rows, values, ADD_VALUES);
void initialize_specific() override
Output class for darcy_flow_mh model.
RegionSet get_region_set(const std::string &set_name) const
void assembly_mh_matrix(MultidimAssembly &assembler)
void postprocess() override
virtual void start_add_assembly()
virtual PetscErrorCode mat_zero_entries()
bool zero_time_term(bool time_global=false) override
static const Input::Type::Record & type_field_descriptor()
void chkerr(unsigned int ierr)
Replacement of new/delete operator in the spirit of xmalloc.
virtual void finish_assembly()=0
const RegionDB & region_db() const
#define ASSERT(expr)
Allow use shorter versions of macro names if these names is not used with external library...
static const std::string field_descriptor_record_description(const string &record_name)
Basic time management functionality for unsteady (and steady) solvers (class Equation).
bool solution_changed_for_scatter
std::shared_ptr< EqData > data_
Basic time management class.
std::unordered_map< unsigned int, unsigned int > edge_new_local_4_mesh_idx_
void prepare_new_time_step() override
postprocess velocity field (add sources)
Assembly explicit Schur complement for the given linear system. Provides method for resolution of the...
VecScatter solution_2_edge_scatter_
PETSC scatter from the solution vector to the parallel edge vector with ghost values.
Global macros to enhance readability and debugging, general constants.
void setup_time_term() override
std::shared_ptr< Balance > balance_
Edge lumped mixed-hybrid solution of unsteady Darcy flow.
const Vec & get_solution()
#define START_TIMER(tag)
Starts a timer with specified tag.
static Input::Type::Abstract & get_input_type()
virtual PetscErrorCode rhs_zero_entries()
std::shared_ptr< EqData > data_
Input::Record input_record_
static const Input::Type::Record & get_input_type()
RichardsLMH(Mesh &mesh, const Input::Record in_rec, TimeGovernor *tm=nullptr)
LocalElementAccessorBase< 3 > accessor(uint local_ele_idx)
void read_initial_condition() override
void assembly_linear_system() override
void set_matrix_changed()
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)
static const Input::Type::Record & get_input_type()
void assembly_source_term() override
Source term is implemented differently in LMH version.
#define FEAL_ASSERT(expr)
Definition of assert for debug and release mode.
Mixed-hybrid model of linear Darcy flow, possibly unsteady.
unsigned int lsize(int proc) const
get local size