OGS 6.2.0-97-g4a610c866
SmallDeformationLocalAssemblerFracture.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #include <vector>
13 
15 
20 
22 #include "SecondaryData.h"
24 
25 namespace ProcessLib
26 {
27 namespace LIE
28 {
29 namespace SmallDeformation
30 {
31 template <typename ShapeFunction, typename IntegrationMethod,
32  int DisplacementDim>
35 {
36 public:
37  using ShapeMatricesType =
43 
49 
51 
56 
58  MeshLib::Element const& e,
59  std::size_t const n_variables,
60  std::size_t const local_matrix_size,
61  std::vector<unsigned> const& dofIndex_to_localIndex,
62  bool const is_axially_symmetric,
63  unsigned const integration_order,
65 
66  void assemble(double const /*t*/, std::vector<double> const& /*local_x*/,
67  std::vector<double>& /*local_M_data*/,
68  std::vector<double>& /*local_K_data*/,
69  std::vector<double>& /*local_b_data*/) override
70  {
71  OGS_FATAL(
72  "SmallDeformationLocalAssembler: assembly without jacobian is not "
73  "implemented.");
74  }
75 
76  void assembleWithJacobian(double const t,
77  Eigen::VectorXd const& local_u,
78  Eigen::VectorXd& local_b,
79  Eigen::MatrixXd& local_J) override;
80 
81  void preTimestepConcrete(std::vector<double> const& /*local_x*/,
82  double const /*t*/,
83  double const /*delta_t*/) override
84  {
85  unsigned const n_integration_points =
86  _integration_method.getNumberOfPoints();
87 
88  for (unsigned ip = 0; ip < n_integration_points; ip++)
89  {
90  _ip_data[ip].pushBackState();
91  }
92  }
93 
95  const double t, Eigen::VectorXd const& local_x) override;
96 
97  Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
98  const unsigned integration_point) const override
99  {
100  auto const& N = _secondary_data.N[integration_point];
101 
102  // assumes N is stored contiguously in memory
103  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
104  }
105 
106  std::vector<double> const& getIntPtSigmaXX(
107  const double /*t*/,
108  GlobalVector const& /*current_solution*/,
109  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
110  std::vector<double>& cache) const override
111  {
112  return getIntPtSigma(cache, 0);
113  }
114 
115  std::vector<double> const& getIntPtSigmaYY(
116  const double /*t*/,
117  GlobalVector const& /*current_solution*/,
118  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
119  std::vector<double>& cache) const override
120  {
121  return getIntPtSigma(cache, 1);
122  }
123 
124  std::vector<double> const& getIntPtSigmaZZ(
125  const double /*t*/,
126  GlobalVector const& /*current_solution*/,
127  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
128  std::vector<double>& cache) const override
129  {
130  return getIntPtSigma(cache, 2);
131  }
132 
133  std::vector<double> const& getIntPtSigmaXY(
134  const double /*t*/,
135  GlobalVector const& /*current_solution*/,
136  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
137  std::vector<double>& cache) const override
138  {
139  return getIntPtSigma(cache, 3);
140  }
141 
142  std::vector<double> const& getIntPtSigmaXZ(
143  const double /*t*/,
144  GlobalVector const& /*current_solution*/,
145  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
146  std::vector<double>& cache) const override
147  {
148  assert(DisplacementDim == 3);
149  return getIntPtSigma(cache, 4);
150  }
151 
152  std::vector<double> const& getIntPtSigmaYZ(
153  const double /*t*/,
154  GlobalVector const& /*current_solution*/,
155  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
156  std::vector<double>& cache) const override
157  {
158  assert(DisplacementDim == 3);
159  return getIntPtSigma(cache, 5);
160  }
161 
162  std::vector<double> const& getIntPtEpsilonXX(
163  const double /*t*/,
164  GlobalVector const& /*current_solution*/,
165  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
166  std::vector<double>& cache) const override
167  {
168  return getIntPtEpsilon(cache, 0);
169  }
170 
171  std::vector<double> const& getIntPtEpsilonYY(
172  const double /*t*/,
173  GlobalVector const& /*current_solution*/,
174  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
175  std::vector<double>& cache) const override
176  {
177  return getIntPtEpsilon(cache, 1);
178  }
179 
180  std::vector<double> const& getIntPtEpsilonZZ(
181  const double /*t*/,
182  GlobalVector const& /*current_solution*/,
183  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
184  std::vector<double>& cache) const override
185  {
186  return getIntPtEpsilon(cache, 2);
187  }
188 
189  std::vector<double> const& getIntPtEpsilonXY(
190  const double /*t*/,
191  GlobalVector const& /*current_solution*/,
192  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
193  std::vector<double>& cache) const override
194  {
195  return getIntPtEpsilon(cache, 3);
196  }
197 
198  std::vector<double> const& getIntPtEpsilonXZ(
199  const double /*t*/,
200  GlobalVector const& /*current_solution*/,
201  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
202  std::vector<double>& cache) const override
203  {
204  assert(DisplacementDim == 3);
205  return getIntPtEpsilon(cache, 4);
206  }
207 
208  std::vector<double> const& getIntPtEpsilonYZ(
209  const double /*t*/,
210  GlobalVector const& /*current_solution*/,
211  NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
212  std::vector<double>& cache) const override
213  {
214  assert(DisplacementDim == 3);
215  return getIntPtEpsilon(cache, 5);
216  }
217 
218 private:
219  std::vector<double> const& getIntPtSigma(
220  std::vector<double>& cache, std::size_t const /*component*/) const
221  {
222  cache.resize(_ip_data.size());
223 
224  return cache;
225  }
226 
227  std::vector<double> const& getIntPtEpsilon(
228  std::vector<double>& cache, std::size_t const /*component*/) const
229  {
230  cache.resize(_ip_data.size());
231 
232  return cache;
233  }
234 
236  std::vector<FractureProperty*> _fracture_props;
237  std::vector<JunctionProperty*> _junction_props;
238  std::unordered_map<int, int> _fracID_to_local;
240 
241  std::vector<IntegrationPointDataFracture<HMatricesType, DisplacementDim>,
242  Eigen::aligned_allocator<IntegrationPointDataFracture<
243  HMatricesType, DisplacementDim>>>
245 
246  IntegrationMethod _integration_method;
247  std::vector<ShapeMatrices, Eigen::aligned_allocator<
252 };
253 
254 } // namespace SmallDeformation
255 } // namespace LIE
256 } // namespace ProcessLib
257 
std::vector< double > const & getIntPtEpsilonZZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtSigma(std::vector< double > &cache, std::size_t const) const
MatrixType< _number_of_dof, _number_of_dof > StiffnessMatrixType
Definition: HMatrixUtils.h:38
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
Definition: SecondaryData.h:29
std::vector< double > const & getIntPtEpsilonYZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
std::vector< ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > _shape_matrices
void assembleWithJacobian(double const t, Eigen::VectorXd const &local_u, Eigen::VectorXd &local_b, Eigen::MatrixXd &local_J) override
void assemble(double const, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
std::vector< double > const & getIntPtEpsilonXZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
void preTimestepConcrete(std::vector< double > const &, double const, double const) override
std::vector< double > const & getIntPtEpsilon(std::vector< double > &cache, std::size_t const) const
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
VectorType< DisplacementDim > ForceVectorType
Definition: HMatrixUtils.h:44
std::vector< double > const & getIntPtSigmaXX(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtEpsilonYY(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
Coordinates mapping matrices at particular location.
Definition: ShapeMatrices.h:45
std::vector< double > const & getIntPtSigmaXY(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
VectorType< ShapeFunction::NPOINTS > NodalVectorType
std::vector< double > const & getIntPtSigmaXZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
MatrixType< DisplacementDim, _number_of_dof > HMatrixType
Definition: HMatrixUtils.h:41
std::vector< double > const & getIntPtSigmaZZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtEpsilonXX(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtEpsilonXY(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
SmallDeformationLocalAssemblerFracture(SmallDeformationLocalAssemblerFracture const &)=delete
VectorType< _number_of_dof > NodalForceVectorType
Definition: HMatrixUtils.h:39
void computeSecondaryVariableConcreteWithVector(const double t, Eigen::VectorXd const &local_x) override
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::vector< double > const & getIntPtSigmaYZ(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
std::vector< IntegrationPointDataFracture< HMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointDataFracture< HMatricesType, DisplacementDim > > > _ip_data
std::vector< double > const & getIntPtSigmaYY(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const override