Flow123d  master-f44eb46
elasticity.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 elasticity.hh
15  * @brief FEM for linear elasticity.
16  * @author Jan Stebel
17  */
18 
19 #ifndef ELASTICITY_HH_
20 #define ELASTICITY_HH_
21 
22 #include "fields/bc_field.hh"
23 #include "fields/field.hh"
24 #include "fields/field_fe.hh"
25 #include "fields/multi_field.hh"
26 #include "la/linsys.hh"
27 #include "la/vector_mpi.hh"
29 #include "coupling/equation.hh"
30 #include "fem/fe_values_views.hh"
31 #include "tools/mixed.hh"
33 
34 class Distribution;
35 class OutputTime;
36 class DOFHandlerMultiDim;
37 template<unsigned int dim> class FiniteElement;
38 class Elasticity;
39 template<unsigned int dim> class StiffnessAssemblyElasticity;
40 template<unsigned int dim> class RhsAssemblyElasticity;
41 template<unsigned int dim> class ConstraintAssemblyElasticity;
42 template<unsigned int dim> class OutpuFieldsAssemblyElasticity;
43 template< template<IntDim...> class DimAssembly> class GenericAssembly;
44 
45 
46 
47 
48 class Elasticity : public EquationBase
49 {
50 public:
51 
52  class EqFields : public FieldSet {
53  public:
54 
55  enum Bc_types {
60  };
61 
62  EqFields();
63 
65 
73  Field<3, FieldValue<3>::Scalar> fracture_sigma; ///< Transition parameter for diffusive transfer on fractures.
75 
76  /// Pointer to DarcyFlow field cross_section
79  Field<3, FieldValue<3>::Scalar > potential_load; ///< Potential of an additional (external) load.
80  Field<3, FieldValue<3>::Scalar > ref_potential_load; ///< Potential of reference external load on boundary. TODO: Switch to BCField when possible.
83 
90 
91  /// @name Instances of FieldModel used in assembly methods
92  // @{
93 
97 
98  // @}
99 
100  std::shared_ptr<FieldFE<3, FieldValue<3>::VectorFixed> > output_field_ptr;
101  std::shared_ptr<FieldFE<3, FieldValue<3>::TensorFixed> > output_stress_ptr;
102  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar> > output_von_mises_stress_ptr;
103  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar> > output_mean_stress_ptr;
104  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar> > output_cross_section_ptr;
105  std::shared_ptr<FieldFE<3, FieldValue<3>::Scalar> > output_div_ptr;
106 
108 
109  };
110 
111  class EqData {
112  public:
113 
115  : ls(nullptr), constraint_matrix(nullptr), constraint_vec(nullptr) {}
116 
118  if (ls!=nullptr) delete ls;
119  if (constraint_matrix!=nullptr) MatDestroy(&constraint_matrix);
120  if (constraint_vec!=nullptr) VecDestroy(&constraint_vec);
121  }
122 
123  /// Create DOF handler objects
124  void create_dh(Mesh * mesh, unsigned int fe_order);
125 
126  /// Objects for distribution of dofs.
127  std::shared_ptr<DOFHandlerMultiDim> dh_;
128  std::shared_ptr<DOFHandlerMultiDim> dh_scalar_;
129  std::shared_ptr<DOFHandlerMultiDim> dh_tensor_;
130 
131  /// @name Solution of algebraic system
132  // @{
133 
134  /// Linear algebraic system.
136 
139 
140  // map local element -> constraint index
142 
143  // @}
144 
145  /// Shared Balance object
146  std::shared_ptr<Balance> balance_;
147 
148  };
149 
150 
151  /**
152  * @brief Constructor.
153  * @param init_mesh computational mesh
154  * @param in_rec input record
155  * @param tm time governor (if nullptr then it is created from input record)
156  */
157  Elasticity(Mesh &init_mesh, const Input::Record in_rec, TimeGovernor *tm = nullptr);
158  /**
159 
160  * @brief Declare input record type for the equation TransportDG.
161  */
162  static const Input::Type::Record & get_input_type();
163 
164  /**
165  * @brief Initialize solution in the zero time.
166  */
167  void zero_time_step() override;
168 
169  /**
170  * @brief Computes the solution in one time instant.
171  */
172  void update_solution() override;
173 
174  /// Pass to next time and update equation data.
175  void next_time();
176 
177  /// Solve without updating time step and without output.
178  void solve_linear_system();
179 
180  /**
181  * @brief Postprocesses the solution and writes to output file.
182  */
183  void output_data();
184 
185  /**
186  * @brief Destructor.
187  */
188  ~Elasticity();
189 
190  void initialize() override;
191 
192  // Recompute fields for output (stress, divergence etc.)
193  void update_output_fields();
194 
195  void set_potential_load(const Field<3, FieldValue<3>::Scalar> &potential,
196  const Field<3, FieldValue<3>::Scalar> &ref_potential)
197  {
198  eq_fields_->potential_load = potential;
199  eq_fields_->ref_potential_load = ref_potential;
200  }
201 
203 
204  const Vec &get_solution()
205  { return eq_data_->ls->get_solution(); }
206 
207  inline EqFields &eq_fields() { return *eq_fields_; }
208 
209  inline EqData &eq_data() { return *eq_data_; }
210 
211 
212 
214 
215 
216 
217 
218 private:
219  /// Registrar of class to factory
220  static const int registrar;
221 
222  void preallocate();
223 
224 
226 
227 
228  /// @name Physical parameters
229  // @{
230 
231  /// Fields for model parameters.
232  std::shared_ptr<EqFields> eq_fields_;
233 
234  /// Data for model parameters.
235  std::shared_ptr<EqData> eq_data_;
236 
237  /// Indicator of contact conditions on fractures.
239 
240 
241  // @}
242 
243 
244 
245  /// @name Output to file
246  // @{
247 
248  std::shared_ptr<OutputTime> output_stream_;
249 
250  /// Record with input specification.
252 
253 
254  // @}
255 
256 
257  static constexpr const char * name_ = "Mechanics_LinearElasticity";
258 
259 
260  /// general assembly objects, hold assembly objects of appropriate dimension
265 
266 };
267 
268 
269 /*
270  * TODO Remove these two methods after implementation new assembly algorithm in HM_Iterative class.
271  */
272 double lame_mu(double young, double poisson);
273 double lame_lambda(double young, double poisson);
274 
275 
276 
277 
278 
279 
280 #endif /* ELASTICITY_HH_ */
Provides the numbering of the finite element degrees of freedom on the computational mesh.
Definition: dofhandler.hh:151
LinSys * ls
Linear algebraic system.
Definition: elasticity.hh:135
std::shared_ptr< DOFHandlerMultiDim > dh_scalar_
Definition: elasticity.hh:128
void create_dh(Mesh *mesh, unsigned int fe_order)
Create DOF handler objects.
Definition: elasticity.cc:282
std::shared_ptr< DOFHandlerMultiDim > dh_tensor_
Definition: elasticity.hh:129
std::shared_ptr< DOFHandlerMultiDim > dh_
Objects for distribution of dofs.
Definition: elasticity.hh:127
std::map< LongIdx, LongIdx > constraint_idx
Definition: elasticity.hh:141
std::shared_ptr< Balance > balance_
Shared Balance object.
Definition: elasticity.hh:146
Field< 3, FieldValue< 3 >::Scalar > potential_load
Potential of an additional (external) load.
Definition: elasticity.hh:79
Field< 3, FieldValue< 3 >::TensorFixed > output_stress
Definition: elasticity.hh:85
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::TensorFixed > > output_stress_ptr
Definition: elasticity.hh:101
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > output_von_mises_stress_ptr
Definition: elasticity.hh:102
Field< 3, FieldValue< 3 >::Scalar > cross_section_min
Definition: elasticity.hh:78
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > output_cross_section_ptr
Definition: elasticity.hh:104
Field< 3, FieldValue< 3 >::Scalar > subdomain
Definition: elasticity.hh:82
Field< 3, FieldValue< 3 >::Scalar > region_id
Definition: elasticity.hh:81
EquationOutput output_fields
Definition: elasticity.hh:107
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::VectorFixed > > output_field_ptr
Definition: elasticity.hh:100
BCField< 3, FieldValue< 3 >::Enum > bc_type
Definition: elasticity.hh:66
Field< 3, FieldValue< 3 >::Scalar > lame_lambda
Definition: elasticity.hh:95
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > output_div_ptr
Definition: elasticity.hh:105
Field< 3, FieldValue< 3 >::Scalar > output_mean_stress
Definition: elasticity.hh:87
BCField< 3, FieldValue< 3 >::VectorFixed > bc_displacement
Definition: elasticity.hh:67
Field< 3, FieldValue< 3 >::Scalar > poisson_ratio
Definition: elasticity.hh:72
Field< 3, FieldValue< 3 >::VectorFixed > output_field
Definition: elasticity.hh:84
Field< 3, FieldValue< 3 >::TensorFixed > initial_stress
Definition: elasticity.hh:74
static const Input::Type::Selection & get_bc_type_selection()
Definition: elasticity.cc:117
BCField< 3, FieldValue< 3 >::VectorFixed > bc_traction
Definition: elasticity.hh:68
Field< 3, FieldValue< 3 >::Scalar > output_von_mises_stress
Definition: elasticity.hh:86
Field< 3, FieldValue< 3 >::Scalar > ref_potential_load
Potential of reference external load on boundary. TODO: Switch to BCField when possible.
Definition: elasticity.hh:80
Field< 3, FieldValue< 3 >::VectorFixed > load
Definition: elasticity.hh:70
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > output_mean_stress_ptr
Definition: elasticity.hh:103
Field< 3, FieldValue< 3 >::Scalar > output_divergence
Definition: elasticity.hh:89
Field< 3, FieldValue< 3 >::Scalar > fracture_sigma
Transition parameter for diffusive transfer on fractures.
Definition: elasticity.hh:73
Field< 3, FieldValue< 3 >::Scalar > young_modulus
Definition: elasticity.hh:71
Field< 3, FieldValue< 3 >::Scalar > dirichlet_penalty
Definition: elasticity.hh:96
Field< 3, FieldValue< 3 >::Scalar > cross_section
Pointer to DarcyFlow field cross_section.
Definition: elasticity.hh:77
Field< 3, FieldValue< 3 >::Scalar > output_cross_section
Definition: elasticity.hh:88
BCField< 3, FieldValue< 3 >::TensorFixed > bc_stress
Definition: elasticity.hh:69
Field< 3, FieldValue< 3 >::Scalar > lame_mu
Definition: elasticity.hh:94
std::shared_ptr< EqData > eq_data_
Data for model parameters.
Definition: elasticity.hh:235
static const int registrar
Registrar of class to factory.
Definition: elasticity.hh:220
void update_solution() override
Computes the solution in one time instant.
Definition: elasticity.cc:548
GenericAssembly< RhsAssemblyElasticity > * rhs_assembly_
Definition: elasticity.hh:262
void update_output_fields()
Definition: elasticity.cc:471
void assemble_constraint_matrix()
Definition: elasticity.cc:641
GenericAssembly< OutpuFieldsAssemblyElasticity > * output_fields_assembly_
Definition: elasticity.hh:264
static constexpr const char * name_
Definition: elasticity.hh:257
EqFields & eq_fields()
Definition: elasticity.hh:207
GenericAssembly< StiffnessAssemblyElasticity > * stiffness_assembly_
general assembly objects, hold assembly objects of appropriate dimension
Definition: elasticity.hh:261
void solve_linear_system()
Solve without updating time step and without output.
Definition: elasticity.cc:571
std::shared_ptr< OutputTime > output_stream_
Definition: elasticity.hh:248
void initialize() override
Definition: elasticity.cc:349
void set_potential_load(const Field< 3, FieldValue< 3 >::Scalar > &potential, const Field< 3, FieldValue< 3 >::Scalar > &ref_potential)
Definition: elasticity.hh:195
Input::Record input_rec
Record with input specification.
Definition: elasticity.hh:251
void calculate_cumulative_balance()
Definition: elasticity.cc:632
void output_data()
Postprocesses the solution and writes to output file.
Definition: elasticity.cc:611
~Elasticity()
Destructor.
Definition: elasticity.cc:456
static const Input::Type::Record & get_input_type()
Declare input record type for the equation TransportDG.
Definition: elasticity.cc:48
bool has_contact_
Indicator of contact conditions on fractures.
Definition: elasticity.hh:238
GenericAssembly< ConstraintAssemblyElasticity > * constraint_assembly_
Definition: elasticity.hh:263
Elasticity(Mesh &init_mesh, const Input::Record in_rec, TimeGovernor *tm=nullptr)
Constructor.
Definition: elasticity.cc:307
void zero_time_step() override
Initialize solution in the zero time.
Definition: elasticity.cc:501
const Vec & get_solution()
Definition: elasticity.hh:204
Elasticity FactoryBaseType
Definition: elasticity.hh:213
std::shared_ptr< EqFields > eq_fields_
Fields for model parameters.
Definition: elasticity.hh:232
EqData & eq_data()
Definition: elasticity.hh:209
void preallocate()
Definition: elasticity.cc:534
void next_time()
Pass to next time and update equation data.
Definition: elasticity.cc:562
Mesh & mesh()
Definition: equation.hh:181
Container for various descendants of FieldCommonBase.
Definition: field_set.hh:159
Class template representing a field with values dependent on: point, element, and region.
Definition: field.hh:92
Abstract class for the description of a general finite element on a reference simplex in dim dimensio...
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
Template for classes storing finite set of named values.
Definition: mesh.h:362
The class for outputting data during time.
Definition: output_time.hh:51
Basic time management functionality for unsteady (and steady) solvers (class Equation).
double lame_mu(double young, double poisson)
Definition: elasticity.cc:80
double lame_lambda(double young, double poisson)
Definition: elasticity.cc:86
Abstract base class for equation clasess.
Wrappers for linear systems based on MPIAIJ and MATIS format.
unsigned int IntDim
A dimension index type.
Definition: mixed.hh:19
Definitions of particular quadrature rules on simplices.