Flow123d  build_with_4.0.3-86a16ad
first_order_reaction_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 first_order_reaction_base.hh
15  * @brief
16  */
17 
18 #ifndef FIRST_ORDER_REACTION_BASE_H_
19 #define FIRST_ORDER_REACTION_BASE_H_
20 
21 #include <memory> // for shared_ptr
22 #include <string> // for string
23 #include <vector>
25 
26 #include "armadillo"
27 
28 class Mesh;
29 class LinearODESolver;
30 namespace Input { class Record; }
31 
32 /** @brief Base class for linear reactions and decay chain.
33  *
34  * The class implements common interface for linear reactions and decay chains.
35  * One step of the linear reaction or the decay is represented as a product of a reaction matrix and
36  * a vector of concentrations of transported substances on a single element.
37  *
38  * It uses armadillo to compute the reaction matrix which then multiplies to concetration vector.
39  * This class also resolves the choice of the numerical method which is used to compute the reaction matrix.
40  */
42 {
43 public:
44  /// Constructor.
45  FirstOrderReactionBase(Mesh &init_mesh, Input::Record in_rec);
46 
47  /// Destructor.
49 
50  /// Prepares the object to usage.
51  /**
52  * Allocating memory, reading input, initialization of fields.
53  */
54  void initialize() override;
55 
56  /// Moves the model to zero time.
57  /** The assembly of the system matrix is called here.
58  */
59  void zero_time_step() override;
60 
61  /// Updates the solution.
62  /**
63  * Goes through local distribution of elements and calls @p compute_reaction.
64  */
65  void update_solution(void) override;
66 
67 protected:
68  /// Assembles the matrix of the ODEs.
69  /**
70  * We solve the system of \f$N\f$ equations
71  * \f[
72  * \frac{\textrm{d} c_i}{\textrm{d}t}=-\sum\limits_{j=1}^N \frac{M_i}{M_j} \lambda_{i} b_j c_j, \qquad \textrm i=1,\ldots,N
73  * \f]
74  * where \f$M_i, M_j\f$ are the molar masses of the parent substances and products, \f$\lambda_i\f$ are the
75  * reaction rate constants (in case of decays converted from half_lives) and \f$b_j\f$ are the branching ratios.
76  * The constant coefficients \f$\mathbf{A}_{ij}=\frac{M_i}{M_j} \lambda_{i} b_j\f$ are the elements of the system matrix.
77  *
78  * It is pure virtual and must be implemented in descendants.
79  */
80  virtual void assemble_ode_matrix(void) = 0;
81 
82  /// Computes the reaction on a specified element.
83  virtual void compute_reaction(const DHCellAccessor& dh_cell) override;
84 
85  /// Initializes private members of sorption from the input record.
86  virtual void initialize_from_input() = 0;
87 
88  /** Help function to create mapping of substance indices.
89  * Finds a position of a string in specified array.
90  */
91  unsigned int find_subst_name(const std::string &name);
92 
93  /**
94  * Sequence of integers describing an order of isotopes.
95  * substance_ids_[reactant][local_product_idx] = global_substance_idx
96  */
98 
99  /**
100  * Two dimensional array contains mass percentage of every single decay bifurcation on every single row.
101  */
103 
104  /// Number of all transported substances. It is the dimension of the reaction matrix.
105  unsigned int n_substances_;
106 
107  arma::mat reaction_matrix_; ///< Reaction matrix.
108  arma::vec prev_conc_; ///< Column vector storing previous concetrations on an element.
109 
110  arma::mat molar_matrix_; ///< Diagonal matrix with molar masses of substances.
111  arma::mat molar_mat_inverse_; ///< Inverse of @p molar_matrix_.
112 
113  std::shared_ptr<LinearODESolver> linear_ode_solver_;
114 };
115 
116 #endif // FIRST_ORDER_REACTION_BASE_H_
Cell accessor allow iterate over DOF handler cells.
Base class for linear reactions and decay chain.
arma::mat molar_mat_inverse_
Inverse of molar_matrix_.
arma::vec prev_conc_
Column vector storing previous concetrations on an element.
std::shared_ptr< LinearODESolver > linear_ode_solver_
unsigned int n_substances_
Number of all transported substances. It is the dimension of the reaction matrix.
void zero_time_step() override
Moves the model to zero time.
arma::mat molar_matrix_
Diagonal matrix with molar masses of substances.
virtual void compute_reaction(const DHCellAccessor &dh_cell) override
Computes the reaction on a specified element.
std::vector< std::vector< unsigned int > > substance_ids_
arma::mat reaction_matrix_
Reaction matrix.
std::vector< std::vector< double > > bifurcation_
unsigned int find_subst_name(const std::string &name)
virtual void assemble_ode_matrix(void)=0
Assembles the matrix of the ODEs.
~FirstOrderReactionBase(void)
Destructor.
void update_solution(void) override
Updates the solution.
void initialize() override
Prepares the object to usage.
FirstOrderReactionBase(Mesh &init_mesh, Input::Record in_rec)
Constructor.
virtual void initialize_from_input()=0
Initializes private members of sorption from the input record.
Accessor to the data with type Type::Record.
Definition: accessors.hh:291
Class for linear ODE solver.
Definition: mesh.h:362
ArmaMat< double, N, M > mat
Definition: armor.hh:888
ArmaVec< double, N > vec
Definition: armor.hh:885
Abstract linear system class.
Definition: balance.hh:40
Class ReactionTerm is an abstract class representing reaction term in transport.