Flow123d  JS_before_hm-2087-gb06a8f5d8
sorption_base.hh
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file sorption_base.hh
15  * @brief Class SorptionBase is abstract class representing model of sorption in transport.
16  *
17  * The sorption is described by several types of isotherms - linear, Freundlich or Langmuir.
18  * Limited solubility can be considered.
19  *
20  * Interpolation tables are used to speed up evaluation of isotherms.
21  *
22  */
23 
24 #ifndef SORPTION_BASE_H
25 #define SORPTION_BASE_H
26 
27 
28 #include <memory> // for shared_ptr
29 #include <string> // for string
30 #include <vector>
31 #include "reaction/reaction_term.hh" // for ReactionTerm
32 #include "fields/field.hh" // for Field
33 #include "fields/field_values.hh" // for FieldValue<>::Scalar, FieldVa...
34 #include "fields/field_set.hh"
35 #include "fields/multi_field.hh"
36 #include "la/vector_mpi.hh"
38 #include "input/input_exception.hh" // for DECLARE_INPUT_EXCEPTION, Exce...
39 #include "input/type_base.hh" // for Array
40 #include "input/type_generic.hh" // for Instance
41 #include "petscvec.h" // for Vec, VecScatter, _p_VecScatter
42 //
43 #include "system/exceptions.hh" // for operator<<, ExcStream, EI
44 
45 class Isotherm;
46 class Mesh;
47 namespace Input {
48  class Record;
49  namespace Type {
50  class Record;
51  class Selection;
52  }
53 }
54 template <int spacedim> class ElementAccessor;
55 
56 
57 
58 
60 {
61 public:
62  TYPEDEF_ERR_INFO( EI_ArrayName, std::string);
63  TYPEDEF_ERR_INFO( EI_Subst, unsigned int);
64  DECLARE_INPUT_EXCEPTION( ExcSubstanceCountMatch, << "The size of the input array " << EI_ArrayName::qval
65  << " does not match the number of substances.");
66  DECLARE_EXCEPTION( ExcNotPositiveScaling,
67  << "Scaling parameter in sorption is not positive. Check the input for rock density and molar mass of " << EI_Subst::val << ". substance.\n" );
68 
69  /**
70  * Static variable for new input data types input
71  */
72  static const Input::Type::Record & get_input_type();
73 
74  static Input::Type::Instance make_output_type(const string &equation_name, const string &output_field_name, const string &output_field_desc )
75  {
76  return EqData(output_field_name, output_field_desc).output_fields.make_output_type(equation_name, "");
77  }
78 
79  class EqData : public FieldSet
80  {
81  public:
82  /**
83  * Sorption type specifies a kind of equilibrial description of adsorption.
84  */
86 
87  /// Collect all fields
88  EqData(const string &output_field_name, const string &output_field_desc);
89 
90  MultiField<3, FieldValue<3>::Enum > sorption_type; ///< Discrete need Selection for initialization.
91  Field<3, FieldValue<3>::Scalar > rock_density; ///< Rock matrix density.
92 
93  /// Multiplication coefficients (k, omega) for all types of isotherms.
94  /** Langmuir: c_s = omega * (alpha*c_a)/(1- alpha*c_a), Linear: c_s = k*c_a */
96  /// Langmuir sorption coeficients alpha (in fraction c_s = omega * (alpha*c_a)/(1- alpha*c_a)).
98 
99  MultiField<3, FieldValue<3>::Scalar> init_conc_solid; ///< Initial sorbed concentrations.
100  Field<3, FieldValue<3>::Scalar > porosity; ///< Porosity field copied from transport.
101 
102  MultiField<3, FieldValue<3>::Scalar> conc_solid; ///< Calculated sorbed concentrations, for output only.
103  FieldFEScalarVec conc_solid_fe; ///< Underlaying FieldFE for each substance of conc_solid.
104 
105  /// Input data set - fields in this set are read from the input file.
107 
108  /// Fields indended for output, i.e. all input fields plus those representing solution.
110 
111  };
112 
113  /**
114  * Constructor with parameter for initialization of a new declared class member
115  */
116  SorptionBase(Mesh &init_mesh, Input::Record in_rec);
117  /**
118  * Destructor.
119  */
120  virtual ~SorptionBase(void);
121 
122  /// Prepares the object to usage.
123  /**
124  * Allocating memory, reading input, initialization of fields.
125  */
126  void initialize() override;
127 
128  /**
129  * Does first computation after initialization process.
130  * The time is set and initial condition is set and output.
131  */
132  void zero_time_step() override;
133 
134  /// Updates the solution.
135  /**
136  * Goes through local distribution of elements and calls @p compute_reaction.
137  */
138  void update_solution(void) override;
139 
140  void output_data(void) override;
141 
142 
143 protected:
144  /**
145  * This method disables to use constructor without parameters.
146  */
147  SorptionBase();
148 
149  /** Initializes possible following reactions from input record.
150  * It should be called after setting mesh, time_governor, distribution and concentration_matrix
151  * if there are some setting methods for reactions called (they are not at the moment, so it could be part of init_from_input).
152  */
153  void make_reactions();
154 
155  /// Reads names of substances from input and creates indexing to global vector of substance.
156  /** Also creates the local vector of molar masses. */
158 
159  /// Initializes private members of sorption from the input record.
160  void initialize_from_input();
161 
162  /// Initializes field sets.
163  void initialize_fields();
164 
165  ///Reads and sets initial condition for concentration in solid.
166  void set_initial_condition();
167 
168  /// Compute reaction on a single element.
169  void compute_reaction(const DHCellAccessor& dh_cell) override;
170 
171  /// Reinitializes the isotherm.
172  /**
173  * On data change the isotherm is recomputed, possibly new interpolation table is made.
174  * NOTE: Be sure to update common element data (porosity, rock density etc.)
175  * by @p compute_common_ele_data(), before calling reinitialization!
176  */
177  void isotherm_reinit(unsigned int i_subst, const ElementAccessor<3> &elm);
178 
179  /// Calls @p isotherm_reinit for all isotherms.
180  void isotherm_reinit_all(const ElementAccessor<3> &elm);
181 
182  /**
183  * Creates interpolation table for isotherms.
184  */
185  void make_tables(void);
186 
187  /// Computes maximal aqueous concentration at the current step.
188  void update_max_conc();
189 
190  /// Sets max conc to zeros on all regins.
191  void clear_max_conc();
192 
193  /// Pointer to equation data. The object is constructed in descendants.
195 
196  /**
197  * Temporary nr_of_points can be computed using step_length. Should be |nr_of_region x nr_of_substances| matrix later.
198  */
200  /**
201  * Density of the solvent.
202  * TODO: Could be done region dependent, easily.
203  */
205  /**
206  * Critical concentrations of species dissolved in water.
207  */
209  /**
210  * Concentration table limits of species dissolved in water.
211  */
213  /**
214  * Maximum concentration per region.
215  * It is used for optimization of interpolation table.
216  */
218  /**
219  * Three dimensional array contains intersections between isotherms and mass balance lines.
220  * It describes behaviour of sorbents in mobile pores of various rock matrix enviroments.
221  * Up to |nr_of_region x nr_of_substances x n_points| doubles. Because of equidistant step
222  * lenght in cocidered system of coordinates, just function values are stored.
223  */
225 
226  unsigned int n_substances_; //< number of substances that take part in the sorption mode
227 
228  /// Mapping from local indexing of substances to global.
230 
231  /**
232  * Reaction model that follows the sorption.
233  */
234  std::shared_ptr<ReactionTerm> reaction_liquid;
235  std::shared_ptr<ReactionTerm> reaction_solid;
236 
237 
238  /** Structure for data respectful to element, but indepedent of actual isotherm.
239  * Reads mobile/immobile porosity, rock density and then computes concentration scaling parameters.
240  * Is kind of optimization, so that these data are computed only when necessary.
241  */
243  double scale_aqua;
244  double scale_sorbed;
246  } common_ele_data;
247 
248  /** Computes @p CommonElementData.
249  * Is pure virtual, implemented differently for simple/mobile/immobile sorption class.
250  */
251  virtual void compute_common_ele_data(const ElementAccessor<3> &elem) = 0;
252 };
253 
254 #endif //SORPTION_BASE_H
SorptionBase::table_limit_
std::vector< double > table_limit_
Definition: sorption_base.hh:212
SorptionBase::get_input_type
static const Input::Type::Record & get_input_type()
Definition: sorption_base.cc:54
SorptionBase::data_
EqData * data_
Pointer to equation data. The object is constructed in descendants.
Definition: sorption_base.hh:194
SorptionBase::reaction_liquid
std::shared_ptr< ReactionTerm > reaction_liquid
Definition: sorption_base.hh:234
SorptionBase::initialize_substance_ids
void initialize_substance_ids()
Reads names of substances from input and creates indexing to global vector of substance.
Definition: sorption_base.cc:207
SorptionBase::reaction_solid
std::shared_ptr< ReactionTerm > reaction_solid
Definition: sorption_base.hh:235
vector_mpi.hh
SorptionBase::DECLARE_INPUT_EXCEPTION
DECLARE_INPUT_EXCEPTION(ExcSubstanceCountMatch,<< "The size of the input array "<< EI_ArrayName::qval<< " does not match the number of substances.")
Input
Abstract linear system class.
Definition: balance.hh:40
SorptionBase::clear_max_conc
void clear_max_conc()
Sets max conc to zeros on all regins.
Definition: sorption_base.cc:438
SorptionBase::initialize_from_input
void initialize_from_input()
Initializes private members of sorption from the input record.
Definition: sorption_base.cc:253
SorptionBase::EqData::output_fields
EquationOutput output_fields
Fields indended for output, i.e. all input fields plus those representing solution.
Definition: sorption_base.hh:109
SorptionBase::initialize_fields
void initialize_fields()
Initializes field sets.
Definition: sorption_base.cc:302
SorptionBase::n_interpolation_steps_
unsigned int n_interpolation_steps_
Definition: sorption_base.hh:199
SorptionBase::update_max_conc
void update_max_conc()
Computes maximal aqueous concentration at the current step.
Definition: sorption_base.cc:449
field_set.hh
std::vector< std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > >
ElementAccessor
Definition: dh_cell_accessor.hh:32
SorptionBase::solubility_vec_
std::vector< double > solubility_vec_
Definition: sorption_base.hh:208
SorptionBase::EqData::sorption_type
MultiField< 3, FieldValue< 3 >::Enum > sorption_type
Discrete need Selection for initialization.
Definition: sorption_base.hh:90
SorptionBase::CommonElementData::no_sorbing_surface_cond
double no_sorbing_surface_cond
Definition: sorption_base.hh:245
type_base.hh
exceptions.hh
SorptionBase::isotherm_reinit
void isotherm_reinit(unsigned int i_subst, const ElementAccessor< 3 > &elm)
Reinitializes the isotherm.
Definition: sorption_base.cc:400
SorptionBase::compute_reaction
void compute_reaction(const DHCellAccessor &dh_cell) override
Compute reaction on a single element.
Definition: sorption_base.cc:530
SorptionBase::output_data
void output_data(void) override
Output method.
Definition: sorption_base.cc:583
SorptionBase
Definition: sorption_base.hh:59
SorptionBase::compute_common_ele_data
virtual void compute_common_ele_data(const ElementAccessor< 3 > &elem)=0
EquationOutput::make_output_type
const Input::Type::Instance & make_output_type(const string &equation_name, const string &aditional_description="")
Definition: equation_output.cc:109
SorptionBase::EqData::init_conc_solid
MultiField< 3, FieldValue< 3 >::Scalar > init_conc_solid
Initial sorbed concentrations.
Definition: sorption_base.hh:99
equation_output.hh
SorptionBase::make_reactions
void make_reactions()
Definition: sorption_base.cc:138
SorptionBase::EqData::EqData
EqData(const string &output_field_name, const string &output_field_desc)
Collect all fields.
Definition: sorption_base.cc:78
SorptionBase::EqData::rock_density
Field< 3, FieldValue< 3 >::Scalar > rock_density
Rock matrix density.
Definition: sorption_base.hh:91
SorptionBase::SorptionBase
SorptionBase()
SorptionBase::make_output_type
static Input::Type::Instance make_output_type(const string &equation_name, const string &output_field_name, const string &output_field_desc)
Definition: sorption_base.hh:74
SorptionBase::EqData
Definition: sorption_base.hh:79
Input::Record
Accessor to the data with type Type::Record.
Definition: accessors.hh:291
SorptionBase::set_initial_condition
void set_initial_condition()
Reads and sets initial condition for concentration in solid.
Definition: sorption_base.cc:364
type_generic.hh
EquationOutput
Definition: equation_output.hh:44
SorptionBase::zero_time_step
void zero_time_step() override
Definition: sorption_base.cc:338
SorptionBase::EqData::get_sorption_type_selection
static const Input::Type::Selection & get_sorption_type_selection()
Definition: sorption_base.cc:40
SorptionBase::EqData::porosity
Field< 3, FieldValue< 3 >::Scalar > porosity
Porosity field copied from transport.
Definition: sorption_base.hh:100
field_values.hh
SorptionBase::CommonElementData::scale_sorbed
double scale_sorbed
Definition: sorption_base.hh:244
Input::Type::Instance
Helper class that stores data of generic types.
Definition: type_generic.hh:89
SorptionBase::~SorptionBase
virtual ~SorptionBase(void)
Definition: sorption_base.cc:133
input_exception.hh
FieldSet
Container for various descendants of FieldCommonBase.
Definition: field_set.hh:159
Input::Type::Selection
Template for classes storing finite set of named values.
Definition: type_selection.hh:65
SorptionBase::EqData::conc_solid_fe
FieldFEScalarVec conc_solid_fe
Underlaying FieldFE for each substance of conc_solid.
Definition: sorption_base.hh:103
SorptionBase::isotherms
std::vector< std::vector< Isotherm > > isotherms
Definition: sorption_base.hh:224
Input::Type::Record
Record type proxy class.
Definition: type_record.hh:182
SorptionBase::update_solution
void update_solution(void) override
Updates the solution.
Definition: sorption_base.cc:380
SorptionBase::solvent_density_
double solvent_density_
Definition: sorption_base.hh:204
SorptionBase::substance_global_idx_
std::vector< unsigned int > substance_global_idx_
Mapping from local indexing of substances to global.
Definition: sorption_base.hh:229
SorptionBase::DECLARE_EXCEPTION
DECLARE_EXCEPTION(ExcNotPositiveScaling,<< "Scaling parameter in sorption is not positive. Check the input for rock density and molar mass of "<< EI_Subst::val<< ". substance.\n")
SorptionBase::make_tables
void make_tables(void)
Definition: sorption_base.cc:465
Mesh
Definition: mesh.h:355
SorptionBase::TYPEDEF_ERR_INFO
TYPEDEF_ERR_INFO(EI_ArrayName, std::string)
SorptionBase::CommonElementData
Definition: sorption_base.hh:242
multi_field.hh
DHCellAccessor
Cell accessor allow iterate over DOF handler cells.
Definition: dh_cell_accessor.hh:43
SorptionBase::EqData::input_data_set_
FieldSet input_data_set_
Input data set - fields in this set are read from the input file.
Definition: sorption_base.hh:106
SorptionBase::n_substances_
unsigned int n_substances_
Definition: sorption_base.hh:226
MultiField
Class for representation of a vector of fields of the same physical quantity.
Definition: multi_field.hh:87
Isotherm
Definition: isotherm.hh:158
SorptionBase::isotherm_reinit_all
void isotherm_reinit_all(const ElementAccessor< 3 > &elm)
Calls isotherm_reinit for all isotherms.
Definition: sorption_base.cc:430
SorptionBase::max_conc
std::vector< std::vector< double > > max_conc
Definition: sorption_base.hh:217
SorptionBase::EqData::conc_solid
MultiField< 3, FieldValue< 3 >::Scalar > conc_solid
Calculated sorbed concentrations, for output only.
Definition: sorption_base.hh:102
SorptionBase::CommonElementData::scale_aqua
double scale_aqua
Definition: sorption_base.hh:243
SorptionBase::initialize
void initialize() override
Prepares the object to usage.
Definition: sorption_base.cc:165
Field
Class template representing a field with values dependent on: point, element, and region.
Definition: field.hh:93
SorptionBase::EqData::isotherm_other
MultiField< 3, FieldValue< 3 >::Scalar > isotherm_other
Langmuir sorption coeficients alpha (in fraction c_s = omega * (alpha*c_a)/(1- alpha*c_a)).
Definition: sorption_base.hh:97
reaction_term.hh
Class ReactionTerm is an abstract class representing reaction term in transport.
ReactionTerm
Definition: reaction_term.hh:45
SorptionBase::EqData::distribution_coefficient
MultiField< 3, FieldValue< 3 >::Scalar > distribution_coefficient
Multiplication coefficients (k, omega) for all types of isotherms.
Definition: sorption_base.hh:95
field.hh
SorptionBase::common_ele_data
struct SorptionBase::CommonElementData common_ele_data