OGS 6.2.1-97-g73d1aeda3
ThermoMechanicsFEM.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #include <memory>
13 #include <vector>
14 
16 #include "MathLib/KelvinVector.h"
21 #include "ParameterLib/Parameter.h"
25 
28 
29 namespace ProcessLib
30 {
31 namespace ThermoMechanics
32 {
33 template <typename BMatricesType, typename ShapeMatricesType,
34  int DisplacementDim>
36 {
40  : solid_material(solid_material),
42  solid_material.createMaterialStateVariables())
43  {
44  }
45 
46  // total stress
48 
49  // total strain
51 
52  // mechanical strain
54 
55  // Non-equilibrium initial stress.
57  BMatricesType::KelvinVectorType::Zero(
59  DisplacementDim>::value);
60 
62  std::unique_ptr<typename MaterialLib::Solids::MechanicsBase<
63  DisplacementDim>::MaterialStateVariables>
65  double solid_density;
67 
69  typename ShapeMatricesType::NodalRowVectorType N;
70  typename ShapeMatricesType::GlobalDimNodalMatrixType dNdx;
71 
73  {
74  eps_prev = eps;
75  eps_m_prev = eps_m;
76  sigma_prev = sigma;
77  solid_density_prev = solid_density;
78  material_state_variables->pushBackState();
79  }
80 
82 };
83 
86 template <typename ShapeMatrixType>
88 {
89  std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>> N;
90 };
91 
92 template <typename ShapeFunction, typename IntegrationMethod,
93  int DisplacementDim>
96 {
97 public:
98  using ShapeMatricesType =
100 
101  // Types for displacement.
102  // (Higher order elements = ShapeFunction).
105 
107  using RhsVector = typename ShapeMatricesType::template VectorType<
108  ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim>;
109  using JacobianMatrix = typename ShapeMatricesType::template MatrixType<
110  ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim,
111  ShapeFunction::NPOINTS + ShapeFunction::NPOINTS * DisplacementDim>;
112 
114  delete;
116 
118  MeshLib::Element const& e,
119  std::size_t const /*local_matrix_size*/,
120  bool const is_axially_symmetric,
121  unsigned const integration_order,
123 
125  std::size_t setIPDataInitialConditions(
126  std::string const& name,
127  double const* values,
128  int const integration_order) override;
129 
130  void assemble(double const /*t*/, std::vector<double> const& /*local_x*/,
131  std::vector<double>& /*local_M_data*/,
132  std::vector<double>& /*local_K_data*/,
133  std::vector<double>& /*local_rhs_data*/) override
134  {
135  OGS_FATAL(
136  "ThermoMechanicsLocalAssembler: assembly without jacobian is not "
137  "implemented.");
138  }
139 
140  void assembleWithJacobianForStaggeredScheme(
141  double const t, std::vector<double> const& local_xdot,
142  const double dxdot_dx, const double dx_dx,
143  std::vector<double>& local_M_data, std::vector<double>& local_K_data,
144  std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
145  LocalCoupledSolutions const& local_coupled_solutions) override;
146 
147  void assembleWithJacobian(double const t,
148  std::vector<double> const& local_x,
149  std::vector<double> const& local_xdot,
150  const double /*dxdot_dx*/, const double /*dx_dx*/,
151  std::vector<double>& /*local_M_data*/,
152  std::vector<double>& /*local_K_data*/,
153  std::vector<double>& local_rhs_data,
154  std::vector<double>& local_Jac_data) override;
155 
156  void initializeConcrete() override
157  {
158  unsigned const n_integration_points =
159  _integration_method.getNumberOfPoints();
160 
161  for (unsigned ip = 0; ip < n_integration_points; ip++)
162  {
163  _ip_data[ip].pushBackState();
164  }
165  }
166 
167  void postTimestepConcrete(std::vector<double> const& /*local_x*/) override
168  {
169  unsigned const n_integration_points =
170  _integration_method.getNumberOfPoints();
171 
172  for (unsigned ip = 0; ip < n_integration_points; ip++)
173  {
174  _ip_data[ip].pushBackState();
175  }
176  }
177 
178  Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
179  const unsigned integration_point) const override
180  {
181  auto const& N = _secondary_data.N[integration_point];
182 
183  // assumes N is stored contiguously in memory
184  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
185  }
186 
187 private:
215  void assembleWithJacobianForDeformationEquations(
216  double const t, std::vector<double> const& local_xdot,
217  const double dxdot_dx, const double dx_dx,
218  std::vector<double>& local_M_data, std::vector<double>& local_K_data,
219  std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
220  LocalCoupledSolutions const& local_coupled_solutions);
221 
247  void assembleWithJacobianForHeatConductionEquations(
248  double const t, std::vector<double> const& local_xdot,
249  const double dxdot_dx, const double dx_dx,
250  std::vector<double>& local_M_data, std::vector<double>& local_K_data,
251  std::vector<double>& local_b_data, std::vector<double>& local_Jac_data,
252  LocalCoupledSolutions const& local_coupled_solutions);
253 
254  std::size_t setSigma(double const* values);
255 
256  // TODO (naumov) This method is same as getIntPtSigma but for arguments and
257  // the ordering of the cache_mat.
258  // There should be only one.
259  std::vector<double> getSigma() const override;
260 
261  std::vector<double> const& getIntPtSigma(
262  const double /*t*/,
263  GlobalVector const& /*current_solution*/,
264  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
265  std::vector<double>& cache) const override;
266 
267  std::size_t setEpsilon(double const* values);
268 
269  std::vector<double> getEpsilon() const override;
270 
271  std::vector<double> const& getIntPtEpsilon(
272  const double /*t*/,
273  GlobalVector const& /*current_solution*/,
274  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
275  std::vector<double>& cache) const override;
276 
277  std::size_t setEpsilonMechanical(double const* values);
278 
279  std::vector<double> getEpsilonMechanical() const override;
280 
282 
283  std::vector<
285  Eigen::aligned_allocator<IntegrationPointData<
286  BMatricesType, ShapeMatricesType, DisplacementDim>>>
288 
289  IntegrationMethod _integration_method;
293 
294  static const int temperature_index = 0;
295  static const int temperature_size = ShapeFunction::NPOINTS;
296  static const int displacement_index = ShapeFunction::NPOINTS;
297  static const int displacement_size =
298  ShapeFunction::NPOINTS * DisplacementDim;
299 };
300 
301 } // namespace ThermoMechanics
302 } // namespace ProcessLib
303 
304 #include "ThermoMechanicsFEM-impl.h"
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:23
IntegrationPointData(MaterialLib::Solids::MechanicsBase< DisplacementDim > const &solid_material)
ShapeMatricesType::GlobalDimNodalMatrixType dNdx
std::unique_ptr< typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables > material_state_variables
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.
Definition: BMatrixPolicy.h:43
ShapeMatricesType::NodalRowVectorType N
ThermoMechanicsProcessData< DisplacementDim > & _process_data
MaterialLib::Solids::MechanicsBase< DisplacementDim > const & solid_material
void postTimestepConcrete(std::vector< double > const &) override
typename ShapeMatricesType::template MatrixType< ShapeFunction::NPOINTS+ShapeFunction::NPOINTS *DisplacementDim, ShapeFunction::NPOINTS+ShapeFunction::NPOINTS *DisplacementDim > JacobianMatrix
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
Coordinates mapping matrices at particular location.
Definition: ShapeMatrices.h:45
void assemble(double const, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
typename BMatricesType::NodalForceVectorType NodalForceVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::vector< IntegrationPointData< BMatricesType, ShapeMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointData< BMatricesType, ShapeMatricesType, DisplacementDim > > > _ip_data
typename ShapeMatricesType::template VectorType< ShapeFunction::NPOINTS+ShapeFunction::NPOINTS *DisplacementDim > RhsVector
Eigen::Matrix< double, KelvinVectorDimensions< DisplacementDim >::value, 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:49