Flow123d  master-f44eb46
richards_lmh.hh
Go to the documentation of this file.
1 /*
2  * richards_lmh.hh
3  *
4  * Created on: Sep 16, 2015
5  * Author: jb
6  */
7 
8 #ifndef SRC_FLOW_RICHARDS_LMH_HH_
9 #define SRC_FLOW_RICHARDS_LMH_HH_
10 
11 
12 #include <memory> // for shared_ptr
13 #include "fields/field.hh" // for Field
14 #include "fields/field_fe.hh" // for FieldFE
15 #include "fields/field_values.hh" // for FieldValue<>::Scalar, FieldValue
16 #include "la/vector_mpi.hh" // for VectorMPI
17 #include "flow/darcy_flow_lmh.hh" // for DarcyLMH, DarcyLMH::EqData
18 #include "flow/darcy_flow_mh_output.hh" // for DarcyFlowMHOutput
19 #include "input/type_base.hh" // for Array
20 #include "input/type_generic.hh" // for Instance
21 
22 class Mesh;
23 class SoilModelBase;
24 namespace Input {
25  class Record;
26  namespace Type {
27  class Record;
28  }
29 }
30 template<unsigned int dim> class InitCondPostprocessAssembly;
31 template<unsigned int dim> class MHMatrixAssemblyRichards;
32 template< template<IntDim...> class DimAssembly> class GenericAssembly;
33 
34 /**
35  * @brief Edge lumped mixed-hybrid solution of unsteady Darcy flow.
36  *
37  * The time term and sources are evenly distributed form an element to its edges.
38  * This applies directly to the second Schur complement. After this system for pressure traces is solved we reconstruct pressures and side flows as follows:
39  *
40  * -# Element pressure is average of edge pressure. This is in fact same as the MH for steady case so we let SchurComplement class do its job.
41  *
42  * -# We let SchurComplement to reconstruct fluxes and then account time term and sources which are evenly distributed from an element to its sides.
43  * It can be proved, that this keeps continuity of the fluxes over the edges.
44  *
45  * This lumping technique preserves discrete maximum principle for any time step provided one use acute mesh. But in practice even worse meshes are tractable.
46  *
47  * Ideas how to unify steady and unsteady flow:
48  * zero_time_step:
49  *
50  * -# Set initial time.
51  * -# Read initial condition. Reconstruct pressures.
52  * -# Assembly system (possibly in matrix free way).
53  * -# Reconstruct velocities (schur complement resolve).
54  * -# Solve iteratively as regular time step if an input flag "steady_initial_time" is set.
55  * -# (Detect that there is no time term. I such case use arbitrary long time step up to next change of data.
56  * Some kind of time step estimator would be nice.
57  *
58  * update solution:
59  *
60  * -# Move to the next time.
61  * -# Update fields
62  * -# Nonlinear solve.
63  * -# In case of slow convergence, use shorter time-step, within estimated limits. Otherwise there is a different problem.
64  */
65 class RichardsLMH : public DarcyLMH
66 {
67 public:
68  /// Class with all fields used in the equation DarcyFlow.
69  /// This is common to all implementations since this provides interface
70  /// to this equation for possible coupling.
71  class EqFields : public DarcyLMH::EqFields {
72  public:
73  EqFields();
74  // input fields
75  Field<3, FieldValue<3>::Scalar > water_content_saturated; // corresponds to the porosity (theta_s = Vw/V = porosity)
79 
80  //output fields
82  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar>> water_content_ptr;
83 
85 // FieldFE<3, FieldValue<3>::Scalar > conductivity_richards;
86  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar>> conductivity_ptr;
87  };
88 
89  class EqData : public DarcyLMH::EqData {
90  public:
91  /// Constructor
92  EqData();
93 
94  // Auxiliary assembly fields.
97 
98  std::shared_ptr<SoilModelBase> soil_model_;
99  };
100 
101  RichardsLMH(Mesh &mesh, const Input::Record in_rec, TimeGovernor *tm = nullptr);
102 
103  static const Input::Type::Record & get_input_type();
104 
105  void accept_time_step() override;
106 
107  virtual ~RichardsLMH() override;
108 
109 protected:
110  /// Registrar of class to factory
111  static const int registrar;
112 
113  bool zero_time_term(bool time_global=false) override;
114 
115  void initialize_specific() override;
116 
117 // void initial_condition_postprocess() override;
118  void assembly_linear_system() override;
119 
120  /// Create and initialize assembly objects
121  void initialize_asm() override;
122 
123  /// Call assemble of read_init_cond_assembly_ and init_cond_postprocess_assembly_
124  void read_init_cond_asm() override;
125 
126 private:
127 
128  std::shared_ptr<EqFields> eq_fields_;
129  std::shared_ptr<EqData> eq_data_;
130 
131  /// general assembly object, hold assembly objects of appropriate dimension
133 
134  static std::string equation_name()
135  { return "Flow_Richards_LMH";}
136 
137 };
138 
139 
140 
141 
142 #endif /* SRC_FLOW_RICHARDS_LMH_HH_ */
Mixed-hybrid model of linear Darcy flow, possibly unsteady.
Mesh & mesh()
Definition: equation.hh:181
Class template representing a field with values dependent on: point, element, and region.
Definition: field.hh:92
Generic class of assemblation.
Accessor to the data with type Type::Record.
Definition: accessors.hh:291
Record type proxy class.
Definition: type_record.hh:182
Definition: mesh.h:362
EqData()
Constructor.
Definition: richards_lmh.cc:79
std::shared_ptr< SoilModelBase > soil_model_
Definition: richards_lmh.hh:98
VectorMPI water_content_previous_time
Definition: richards_lmh.hh:95
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > conductivity_ptr
Definition: richards_lmh.hh:86
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > water_content_ptr
Definition: richards_lmh.hh:82
Field< 3, FieldValue< 3 >::Scalar > genuchten_p_head_scale
Definition: richards_lmh.hh:77
Field< 3, FieldValue< 3 >::Scalar > water_content_saturated
Definition: richards_lmh.hh:75
Field< 3, FieldValue< 3 >::Scalar > water_content
Definition: richards_lmh.hh:81
Field< 3, FieldValue< 3 >::Scalar > water_content_residual
Definition: richards_lmh.hh:76
Field< 3, FieldValue< 3 >::Scalar > conductivity_richards
Definition: richards_lmh.hh:84
Field< 3, FieldValue< 3 >::Scalar > genuchten_n_exponent
Definition: richards_lmh.hh:78
Edge lumped mixed-hybrid solution of unsteady Darcy flow.
Definition: richards_lmh.hh:66
bool zero_time_term(bool time_global=false) override
static const int registrar
Registrar of class to factory.
void read_init_cond_asm() override
Call assemble of read_init_cond_assembly_ and init_cond_postprocess_assembly_.
virtual ~RichardsLMH() override
void accept_time_step() override
postprocess velocity field (add sources)
static const Input::Type::Record & get_input_type()
Definition: richards_lmh.cc:83
GenericAssembly< InitCondPostprocessAssembly > * init_cond_postprocess_assembly_
general assembly object, hold assembly objects of appropriate dimension
void initialize_specific() override
void assembly_linear_system() override
RichardsLMH(Mesh &mesh, const Input::Record in_rec, TimeGovernor *tm=nullptr)
std::shared_ptr< EqFields > eq_fields_
void initialize_asm() override
Create and initialize assembly objects.
std::shared_ptr< EqData > eq_data_
static std::string equation_name()
Basic time management functionality for unsteady (and steady) solvers (class Equation).
Lumped mixed-hybrid model of linear Darcy flow, possibly unsteady.
Output class for darcy_flow_mh model.
unsigned int IntDim
A dimension index type.
Definition: mixed.hh:19
Abstract linear system class.
Definition: balance.hh:40