OGS
TH2MFEM.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <memory>
14#include <vector>
15
34#include "TH2MProcessData.h"
35
36namespace ProcessLib
37{
38namespace TH2M
39{
42template <typename ShapeMatrixType>
44{
45 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>> N_u;
46};
47
48template <typename ShapeFunctionDisplacement, typename ShapeFunctionPressure,
49 int DisplacementDim>
50class TH2MLocalAssembler : public LocalAssemblerInterface<DisplacementDim>
51{
52public:
55
58
59 template <int N>
60 using VectorType =
61 typename ShapeMatricesTypePressure::template VectorType<N>;
62
63 template <int M, int N>
64 using MatrixType =
65 typename ShapeMatricesTypePressure::template MatrixType<M, N>;
66
69
72
73 static int const KelvinVectorSize =
75 using SymmetricTensor = Eigen::Matrix<double, KelvinVectorSize, 1>;
76
77 static constexpr auto& N_u_op = MathLib::eigenBlockMatrixView<
78 DisplacementDim,
80
82
85
87 MeshLib::Element const& e,
88 std::size_t const /*local_matrix_size*/,
89 NumLib::GenericIntegrationMethod const& integration_method,
90 bool const is_axially_symmetric,
92
93private:
96 std::string const& name,
97 double const* values,
98 int const integration_order) override;
99
100 void setInitialConditionsConcrete(std::vector<double> const& local_x_data,
101 double const t,
102 bool const /*use_monolithic_scheme*/,
103 int const /*process_id*/) override;
104
105 void assemble(double const /*t*/, double const /*dt*/,
106 std::vector<double> const& /*local_x*/,
107 std::vector<double> const& /*local_x_prev*/,
108 std::vector<double>& /*local_M_data*/,
109 std::vector<double>& /*local_K_data*/,
110 std::vector<double>& /*local_rhs_data*/) override;
111
112 void assembleWithJacobian(double const t, double const dt,
113 std::vector<double> const& local_x,
114 std::vector<double> const& local_x_prev,
115 std::vector<double>& /*local_M_data*/,
116 std::vector<double>& /*local_K_data*/,
117 std::vector<double>& local_rhs_data,
118 std::vector<double>& local_Jac_data) override;
119
120 void initializeConcrete() override
121 {
122 unsigned const n_integration_points =
124
125 for (unsigned ip = 0; ip < n_integration_points; ip++)
126 {
127 auto& ip_data = _ip_data[ip];
128
129 ParameterLib::SpatialPosition const x_position{
130 std::nullopt, _element.getID(), ip,
133 ShapeFunctionDisplacement,
135
137 if (_process_data.initial_stress != nullptr)
138 {
139 ip_data.sigma_eff =
141 DisplacementDim>((*_process_data.initial_stress)(
142 std::numeric_limits<
143 double>::quiet_NaN() /* time independent */,
144 x_position));
145 }
146
147 double const t = 0; // TODO (naumov) pass t from top
148 ip_data.solid_material.initializeInternalStateVariables(
149 t, x_position, *ip_data.material_state_variables);
150
151 ip_data.pushBackState();
152 }
153 }
154
155 void postTimestepConcrete(Eigen::VectorXd const& /*local_x*/,
156 Eigen::VectorXd const& /*local_x_prev*/,
157 double const /*t*/, double const /*dt*/,
158 bool const /*use_monolithic_scheme*/,
159 int const /*process_id*/) override
160 {
161 unsigned const n_integration_points =
163
164 for (unsigned ip = 0; ip < n_integration_points; ip++)
165 {
166 _ip_data[ip].pushBackState();
167 }
168 }
169
171 double const t, double const dt, Eigen::VectorXd const& local_x,
172 Eigen::VectorXd const& local_x_prev) override;
173
174 Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
175 const unsigned integration_point) const override
176 {
177 auto const& N_u = _secondary_data.N_u[integration_point];
178
179 // assumes N is stored contiguously in memory
180 return Eigen::Map<const Eigen::RowVectorXd>(N_u.data(), N_u.size());
181 }
182
183 std::vector<double> const& getIntPtDarcyVelocityGas(
184 const double t,
185 std::vector<GlobalVector*> const& x,
186 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
187 std::vector<double>& cache) const override;
188 std::vector<double> const& getIntPtDarcyVelocityLiquid(
189 const double t,
190 std::vector<GlobalVector*> const& x,
191 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
192 std::vector<double>& cache) const override;
193 std::vector<double> const& getIntPtDiffusionVelocityVapourGas(
194 const double t,
195 std::vector<GlobalVector*> const& x,
196 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
197 std::vector<double>& cache) const override;
198 std::vector<double> const& getIntPtDiffusionVelocityGasGas(
199 const double t,
200 std::vector<GlobalVector*> const& x,
201 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
202 std::vector<double>& cache) const override;
203 std::vector<double> const& getIntPtDiffusionVelocitySoluteLiquid(
204 const double t,
205 std::vector<GlobalVector*> const& x,
206 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
207 std::vector<double>& cache) const override;
208 std::vector<double> const& getIntPtDiffusionVelocityLiquidLiquid(
209 const double t,
210 std::vector<GlobalVector*> const& x,
211 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_table,
212 std::vector<double>& cache) const override;
213
214 std::vector<ConstitutiveVariables<DisplacementDim>>
215 updateConstitutiveVariables(Eigen::VectorXd const& local_x,
216 Eigen::VectorXd const& local_x_prev,
217 double const t, double const dt);
218
219 std::size_t setSigma(double const* values)
220 {
221 return ProcessLib::setIntegrationPointKelvinVectorData<DisplacementDim>(
222 values, _ip_data, &IpData::sigma_eff);
223 }
224
225 // TODO (naumov) This method is same as getIntPtSigma but for arguments and
226 // the ordering of the cache_mat.
227 // There should be only one.
228 std::vector<double> getSigma() const override
229 {
230 {
231 constexpr int kelvin_vector_size =
233 DisplacementDim);
234
235 return transposeInPlace<kelvin_vector_size>(
236 [this](std::vector<double>& values)
237 { return getIntPtSigma(0, {}, {}, values); });
238 }
239 }
240
241 std::vector<double> const& getIntPtSigma(
242 const double /*t*/,
243 std::vector<GlobalVector*> const& /*x*/,
244 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
245 std::vector<double>& cache) const override
246 {
247 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
248 _ip_data, &IpData::sigma_eff, cache);
249 }
250
251 std::vector<double> getSwellingStress() const override
252 {
253 {
254 constexpr int kelvin_vector_size =
256 DisplacementDim);
257
258 return transposeInPlace<kelvin_vector_size>(
259 [this](std::vector<double>& values)
260 { return getIntPtSwellingStress(0, {}, {}, values); });
261 }
262 }
263
264 std::vector<double> const& getIntPtSwellingStress(
265 const double /*t*/,
266 std::vector<GlobalVector*> const& /*x*/,
267 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
268 std::vector<double>& cache) const override
269 {
270 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
271 _ip_data, &IpData::sigma_sw, cache);
272 }
273
274 std::vector<double> getEpsilon() const override
275 {
276 constexpr int kelvin_vector_size =
278
279 return transposeInPlace<kelvin_vector_size>(
280 [this](std::vector<double>& values)
281 { return getIntPtEpsilon(0, {}, {}, values); });
282 }
283
284 virtual std::vector<double> const& getIntPtEpsilon(
285 const double /*t*/,
286 std::vector<GlobalVector*> const& /*x*/,
287 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
288 std::vector<double>& cache) const override
289 {
290 return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>(
291 _ip_data, &IpData::eps, cache);
292 }
293
294 // TODO: Here is some refactoring potential. All secondary variables could
295 // be stored in some container to avoid defining one method for each
296 // variable.
297
298 virtual std::vector<double> const& getIntPtLiquidDensity(
299 const double /*t*/,
300 std::vector<GlobalVector*> const& /*x*/,
301 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
302 std::vector<double>& cache) const override
303 {
305 &IpData::rhoLR, cache);
306 }
307
308 virtual std::vector<double> const& getIntPtGasDensity(
309 const double /*t*/,
310 std::vector<GlobalVector*> const& /*x*/,
311 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
312 std::vector<double>& cache) const override
313 {
315 &IpData::rhoGR, cache);
316 }
317
318 virtual std::vector<double> const& getIntPtSolidDensity(
319 const double /*t*/,
320 std::vector<GlobalVector*> const& /*x*/,
321 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
322 std::vector<double>& cache) const override
323 {
325 &IpData::rhoSR, cache);
326 }
327
328 virtual std::vector<double> const& getIntPtVapourPressure(
329 const double /*t*/,
330 std::vector<GlobalVector*> const& /*x*/,
331 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
332 std::vector<double>& cache) const override
333 {
335 &IpData::pWGR, cache);
336 }
337
338 virtual std::vector<double> const& getIntPtPorosity(
339 const double /*t*/,
340 std::vector<GlobalVector*> const& /*x*/,
341 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
342 std::vector<double>& cache) const override
343 {
345 cache);
346 }
347
348 std::vector<double> getSaturation() const override
349 {
350 std::vector<double> result;
351 getIntPtSaturation(0, {}, {}, result);
352 return result;
353 }
354
355 virtual std::vector<double> const& getIntPtSaturation(
356 const double /*t*/,
357 std::vector<GlobalVector*> const& /*x*/,
358 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
359 std::vector<double>& cache) const override
360 {
362 cache);
363 }
364
365 virtual std::vector<double> const& getIntPtMoleFractionGas(
366 const double /*t*/,
367 std::vector<GlobalVector*> const& /*x*/,
368 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
369 std::vector<double>& cache) const override
370 {
372 &IpData::xnCG, cache);
373 }
374 virtual std::vector<double> const& getIntPtMassFractionGas(
375 const double /*t*/,
376 std::vector<GlobalVector*> const& /*x*/,
377 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
378 std::vector<double>& cache) const override
379 {
381 &IpData::xmCG, cache);
382 }
383 virtual std::vector<double> const& getIntPtMassFractionLiquid(
384 const double /*t*/,
385 std::vector<GlobalVector*> const& /*x*/,
386 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
387 std::vector<double>& cache) const override
388 {
390 &IpData::xmWL, cache);
391 }
392
393 virtual std::vector<double> const& getIntPtRelativePermeabilityGas(
394 const double /*t*/,
395 std::vector<GlobalVector*> const& /*x*/,
396 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
397 std::vector<double>& cache) const override
398 {
400 _ip_data, &IpData::k_rel_G, cache);
401 }
402 virtual std::vector<double> const& getIntPtRelativePermeabilityLiquid(
403 const double /*t*/,
404 std::vector<GlobalVector*> const& /*x*/,
405 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
406 std::vector<double>& cache) const override
407 {
409 _ip_data, &IpData::k_rel_L, cache);
410 }
411
412 virtual std::vector<double> const& getIntPtIntrinsicPermeability(
413 const double /*t*/,
414 std::vector<GlobalVector*> const& /*x*/,
415 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
416 std::vector<double>& cache) const override
417 {
418 return ProcessLib::getIntegrationPointDimMatrixData<DisplacementDim>(
419 _ip_data, &IpData::k_S, cache);
420 }
421
422 virtual std::vector<double> const& getIntPtEnthalpyGas(
423 const double /*t*/,
424 std::vector<GlobalVector*> const& /*x*/,
425 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
426 std::vector<double>& cache) const override
427 {
429 cache);
430 }
431 virtual std::vector<double> const& getIntPtEnthalpyLiquid(
432 const double /*t*/,
433 std::vector<GlobalVector*> const& /*x*/,
434 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
435 std::vector<double>& cache) const override
436 {
438 cache);
439 }
440 virtual std::vector<double> const& getIntPtEnthalpySolid(
441 const double /*t*/,
442 std::vector<GlobalVector*> const& /*x*/,
443 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
444 std::vector<double>& cache) const override
445 {
447 cache);
448 }
449
450 unsigned getNumberOfIntegrationPoints() const override
451 {
453 }
454
455 int getMaterialID() const override
456 {
457 return _process_data.material_ids == nullptr
458 ? 0
459 : (*_process_data.material_ids)[_element.getID()];
460 }
461
463 std::function<std::span<double>(
465 MaterialStateVariables&)> const& get_values_span,
466 int const& n_components) const override
467 {
470 n_components);
471 }
472
474 DisplacementDim>::MaterialStateVariables const&
475 getMaterialStateVariablesAt(unsigned integration_point) const override
476 {
477 return *_ip_data[integration_point].material_state_variables;
478 }
479
480private:
482
485 using IpData =
487 ShapeMatricesTypePressure, DisplacementDim,
488 ShapeFunctionDisplacement::NPOINTS>;
489 std::vector<IpData, Eigen::aligned_allocator<IpData>> _ip_data;
490
495 typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType>
497
498 // The shape function of pressure has the same form with the shape function
499 // of temperature
500 static const int gas_pressure_index = 0;
501 static const int gas_pressure_size = ShapeFunctionPressure::NPOINTS;
502 static const int capillary_pressure_index = ShapeFunctionPressure::NPOINTS;
503 static const int capillary_pressure_size = ShapeFunctionPressure::NPOINTS;
504 static const int temperature_index = 2 * ShapeFunctionPressure::NPOINTS;
505 static const int temperature_size = ShapeFunctionPressure::NPOINTS;
506 static const int displacement_index = ShapeFunctionPressure::NPOINTS * 3;
507 static const int displacement_size =
508 ShapeFunctionDisplacement::NPOINTS * DisplacementDim;
509
510 static const int C_index = 0;
511 static const int C_size = ShapeFunctionPressure::NPOINTS;
512 static const int W_index = ShapeFunctionPressure::NPOINTS;
513 static const int W_size = ShapeFunctionPressure::NPOINTS;
514};
515
516} // namespace TH2M
517} // namespace ProcessLib
518
519#include "TH2MFEM-impl.h"
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:89
virtual std::vector< double > const & getIntPtEnthalpyGas(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:422
void assemble(double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
void initializeConcrete() override
Definition: TH2MFEM.h:120
void postTimestepConcrete(Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, bool const, int const) override
Definition: TH2MFEM.h:155
static const int gas_pressure_size
Definition: TH2MFEM.h:501
static constexpr auto & N_u_op
Definition: TH2MFEM.h:77
std::vector< double > getMaterialStateVariableInternalState(std::function< std::span< double >(typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables &)> const &get_values_span, int const &n_components) const override
Definition: TH2MFEM.h:462
virtual std::vector< double > const & getIntPtMoleFractionGas(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:365
static const int capillary_pressure_index
Definition: TH2MFEM.h:502
std::vector< ConstitutiveVariables< DisplacementDim > > updateConstitutiveVariables(Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, double const t, double const dt)
Definition: TH2MFEM-impl.h:91
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:241
typename ShapeMatricesTypePressure::template MatrixType< M, N > MatrixType
Definition: TH2MFEM.h:65
Eigen::Matrix< double, KelvinVectorSize, 1 > SymmetricTensor
Definition: TH2MFEM.h:75
std::vector< double > const & getIntPtDiffusionVelocitySoluteLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
std::size_t setSigma(double const *values)
Definition: TH2MFEM.h:219
virtual std::vector< double > const & getIntPtSaturation(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:355
virtual std::vector< double > const & getIntPtRelativePermeabilityLiquid(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:402
virtual std::vector< double > const & getIntPtMassFractionGas(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:374
virtual std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:284
virtual std::vector< double > const & getIntPtEnthalpyLiquid(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:431
ShapeMatrixPolicyType< ShapeFunctionPressure, DisplacementDim > ShapeMatricesTypePressure
Definition: TH2MFEM.h:57
virtual std::vector< double > const & getIntPtRelativePermeabilityGas(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:393
typename ShapeMatricesTypePressure::template VectorType< N > VectorType
Definition: TH2MFEM.h:61
virtual std::vector< double > const & getIntPtSolidDensity(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:318
SecondaryData< typename ShapeMatricesTypeDisplacement::ShapeMatrices::ShapeType > _secondary_data
Definition: TH2MFEM.h:496
virtual std::vector< double > const & getIntPtVapourPressure(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:328
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
TH2MProcessData< DisplacementDim > & _process_data
Definition: TH2MFEM.h:481
MeshLib::Element const & _element
Definition: TH2MFEM.h:492
std::vector< double > const & getIntPtDiffusionVelocityGasGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
Definition: TH2MFEM.h:54
typename ShapeMatricesTypePressure::GlobalDimMatrixType GlobalDimMatrixType
Definition: TH2MFEM.h:68
virtual std::vector< double > const & getIntPtGasDensity(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:308
static const int capillary_pressure_size
Definition: TH2MFEM.h:503
TH2MLocalAssembler(TH2MLocalAssembler &&)=delete
std::vector< double > const & getIntPtDarcyVelocityLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
std::vector< double > const & getIntPtDiffusionVelocityLiquidLiquid(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
static const int temperature_size
Definition: TH2MFEM.h:505
std::vector< double > getSaturation() const override
Definition: TH2MFEM.h:348
unsigned getNumberOfIntegrationPoints() const override
Definition: TH2MFEM.h:450
std::vector< double > getSwellingStress() const override
Definition: TH2MFEM.h:251
virtual std::vector< double > const & getIntPtLiquidDensity(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:298
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
Definition: TH2MFEM.h:174
virtual std::vector< double > const & getIntPtMassFractionLiquid(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:383
virtual std::vector< double > const & getIntPtIntrinsicPermeability(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:412
static const int gas_pressure_index
Definition: TH2MFEM.h:500
std::vector< double > getSigma() const override
Definition: TH2MFEM.h:228
int getMaterialID() const override
Definition: TH2MFEM.h:455
static const int displacement_size
Definition: TH2MFEM.h:507
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
Definition: TH2MFEM.h:489
static int const KelvinVectorSize
Definition: TH2MFEM.h:73
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(unsigned integration_point) const override
Definition: TH2MFEM.h:475
std::vector< double > const & getIntPtDiffusionVelocityVapourGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
BMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > BMatricesType
Definition: TH2MFEM.h:484
virtual std::vector< double > const & getIntPtPorosity(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:338
typename ShapeMatricesTypePressure::GlobalDimVectorType GlobalDimVectorType
Definition: TH2MFEM.h:71
virtual std::vector< double > const & getIntPtEnthalpySolid(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:440
std::size_t setIPDataInitialConditions(std::string const &name, double const *values, int const integration_order) override
Definition: TH2MFEM-impl.h:853
NumLib::GenericIntegrationMethod const & _integration_method
Definition: TH2MFEM.h:491
TH2MLocalAssembler(TH2MLocalAssembler const &)=delete
std::vector< double > getEpsilon() const override
Definition: TH2MFEM.h:274
void setInitialConditionsConcrete(std::vector< double > const &local_x_data, double const t, bool const, int const) override
Definition: TH2MFEM-impl.h:923
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev) override
std::vector< double > const & getIntPtDarcyVelocityGas(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
static const int temperature_index
Definition: TH2MFEM.h:504
std::vector< double > const & getIntPtSwellingStress(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Definition: TH2MFEM.h:264
static const int displacement_index
Definition: TH2MFEM.h:506
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Definition: KelvinVector.h:24
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
Definition: KelvinVector.h:171
constexpr Eigen::CwiseNullaryOp< EigenBlockMatrixViewFunctor< D, M >, typename EigenBlockMatrixViewFunctor< D, M >::Matrix > eigenBlockMatrixView(const Eigen::MatrixBase< M > &matrix)
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
std::vector< double > const & getIntegrationPointScalarData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
std::vector< double > getIntegrationPointDataMaterialStateVariables(IntegrationPointDataVector const &ip_data_vector, MemberType member, std::function< std::span< double >(MaterialStateVariables &)> get_values_span, int const n_components)
MatrixType< GlobalDim, GlobalDim > GlobalDimMatrixType
VectorType< GlobalDim > GlobalDimVectorType
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u
BMatricesType::KelvinVectorType eps
std::unique_ptr< typename MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables > material_state_variables
BMatricesType::KelvinVectorType sigma_eff
BMatricesType::KelvinVectorType sigma_sw
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N_u
Definition: TH2MFEM.h:45