OGS
ThermoMechanicsFEM-impl.h
Go to the documentation of this file.
1
13#pragma once
14
21
22namespace ProcessLib
23{
24namespace ThermoMechanics
25{
26template <typename ShapeFunction, int DisplacementDim>
29 MeshLib::Element const& e,
30 std::size_t const /*local_matrix_size*/,
31 NumLib::GenericIntegrationMethod const& integration_method,
32 bool const is_axially_symmetric,
34 : _process_data(process_data),
35 _integration_method(integration_method),
36 _element(e),
37 _is_axially_symmetric(is_axially_symmetric)
38{
39 unsigned const n_integration_points =
41
42 _ip_data.reserve(n_integration_points);
43 _secondary_data.N.resize(n_integration_points);
44
45 auto const shape_matrices =
47 DisplacementDim>(e, is_axially_symmetric,
49
51 _process_data.solid_materials, _process_data.material_ids, e.getID());
52
53 for (unsigned ip = 0; ip < n_integration_points; ip++)
54 {
55 _ip_data.emplace_back(solid_material);
56 auto& ip_data = _ip_data[ip];
57 ip_data.integration_weight =
59 shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
60
61 static const int kelvin_vector_size =
63 // Initialize current time step values
64 ip_data.sigma.setZero(kelvin_vector_size);
65 ip_data.eps.setZero(kelvin_vector_size);
66
67 // Previous time step values are not initialized and are set later.
68 ip_data.sigma_prev.resize(kelvin_vector_size);
69 ip_data.eps_prev.resize(kelvin_vector_size);
70
71 ip_data.eps_m.setZero(kelvin_vector_size);
72 ip_data.eps_m_prev.setZero(kelvin_vector_size);
74 x_position.setElementID(_element.getID());
75 ip_data.N = shape_matrices[ip].N;
76 ip_data.dNdx = shape_matrices[ip].dNdx;
77
78 _secondary_data.N[ip] = shape_matrices[ip].N;
79 }
80}
81
82template <typename ShapeFunction, int DisplacementDim>
84 setInitialConditionsConcrete(Eigen::VectorXd const local_x,
85 double const /*t*/,
86 int const /*process_id*/)
87{
88 unsigned const n_integration_points =
89 this->_integration_method.getNumberOfPoints();
90 auto const local_u =
91 local_x.template segment<displacement_size>(displacement_index);
92
93 for (unsigned ip = 0; ip < n_integration_points; ip++)
94 {
95 auto const& N = _ip_data[ip].N;
96 auto const& dNdx = _ip_data[ip].dNdx;
97
98 ParameterLib::SpatialPosition const x_position{
99 std::nullopt, _element.getID(),
102 _element, N))};
103
104 auto const x_coord = x_position.getCoordinates().value()[0];
105
106 auto const& B =
107 LinearBMatrix::computeBMatrix<DisplacementDim,
108 ShapeFunction::NPOINTS,
110 dNdx, N, x_coord, _is_axially_symmetric);
111
112 _ip_data[ip].eps.noalias() = B * local_u;
113 }
114}
115
116template <typename ShapeFunction, int DisplacementDim>
118 setIPDataInitialConditions(std::string_view const name,
119 double const* values,
120 int const integration_order)
121{
122 if (integration_order !=
123 static_cast<int>(_integration_method.getIntegrationOrder()))
124 {
125 OGS_FATAL(
126 "Setting integration point initial conditions; The integration "
127 "order of the local assembler for element {:d} is different from "
128 "the integration order in the initial condition.",
129 _element.getID());
130 }
131
132 if (name == "sigma")
133 {
134 return setSigma(values);
135 }
136 if (name == "epsilon")
137 {
138 return setEpsilon(values);
139 }
140 if (name == "epsilon_m")
141 {
142 return setEpsilonMechanical(values);
143 }
144
145 return 0;
146}
147
148template <typename ShapeFunction, int DisplacementDim>
150 assembleWithJacobian(double const t, double const dt,
151 std::vector<double> const& local_x,
152 std::vector<double> const& local_x_prev,
153 std::vector<double>& local_rhs_data,
154 std::vector<double>& local_Jac_data)
155{
156 auto const local_matrix_size = local_x.size();
157 assert(local_matrix_size == temperature_size + displacement_size);
158
159 auto T = Eigen::Map<typename ShapeMatricesType::template VectorType<
160 temperature_size> const>(local_x.data() + temperature_index,
161 temperature_size);
162
163 auto u = Eigen::Map<typename ShapeMatricesType::template VectorType<
164 displacement_size> const>(local_x.data() + displacement_index,
165 displacement_size);
166 bool const is_u_non_zero = u.norm() > 0.0;
167
168 auto T_prev = Eigen::Map<typename ShapeMatricesType::template VectorType<
169 temperature_size> const>(local_x_prev.data() + temperature_index,
170 temperature_size);
171
173 local_Jac_data, local_matrix_size, local_matrix_size);
174
175 auto local_rhs = MathLib::createZeroedVector<RhsVector>(local_rhs_data,
176 local_matrix_size);
177
178 typename ShapeMatricesType::template MatrixType<displacement_size,
179 temperature_size>
180 KuT;
181 KuT.setZero(displacement_size, temperature_size);
182
184 KTT.setZero(temperature_size, temperature_size);
185
187 DTT.setZero(temperature_size, temperature_size);
188
189 unsigned const n_integration_points =
190 _integration_method.getNumberOfPoints();
191
192 MPL::VariableArray variables;
193 MPL::VariableArray variables_prev;
195 x_position.setElementID(_element.getID());
196
197 auto const& medium = _process_data.media_map.getMedium(_element.getID());
198 auto const& solid_phase = medium->phase("Solid");
199
200 for (unsigned ip = 0; ip < n_integration_points; ip++)
201 {
202 auto const& w = _ip_data[ip].integration_weight;
203 auto const& N = _ip_data[ip].N;
204 auto const& dNdx = _ip_data[ip].dNdx;
205
206 ParameterLib::SpatialPosition const x_position{
207 std::nullopt, _element.getID(),
210 _element, N))};
211
212 auto const x_coord = x_position.getCoordinates().value()[0];
213 auto const& B =
214 LinearBMatrix::computeBMatrix<DisplacementDim,
215 ShapeFunction::NPOINTS,
217 dNdx, N, x_coord, _is_axially_symmetric);
218
219 auto& sigma = _ip_data[ip].sigma;
220 auto const& sigma_prev = _ip_data[ip].sigma_prev;
221
222 auto& eps = _ip_data[ip].eps;
223 auto const& eps_prev = _ip_data[ip].eps_prev;
224
225 auto& eps_m = _ip_data[ip].eps_m;
226 auto const& eps_m_prev = _ip_data[ip].eps_m_prev;
227
228 auto& state = _ip_data[ip].material_state_variables;
229
230 const double T_ip = N.dot(T); // T at integration point
231 double const T_prev_ip = N.dot(T_prev);
232
233 // Consider also anisotropic thermal expansion.
234 auto const solid_linear_thermal_expansivity_vector =
236 solid_phase
237 .property(
239 .value(variables, x_position, t, dt));
240
242 dthermal_strain =
243 solid_linear_thermal_expansivity_vector * (T_ip - T_prev_ip);
244
245 //
246 // displacement equation, displacement part
247 //
248 // For the restart computation, the displacement may not be
249 // reloaded but the initial strains are always available. For such case,
250 // the following computation is skipped.
251 if (is_u_non_zero)
252 {
253 eps.noalias() = B * u;
254 }
255
256 eps_m.noalias() = eps_m_prev + eps - eps_prev - dthermal_strain;
257
258 variables_prev.stress
260 sigma_prev);
261 variables_prev.mechanical_strain
263 eps_m_prev);
264 variables_prev.temperature = T_ip;
265 variables.mechanical_strain
267 eps_m);
268 variables.temperature = T_ip;
269
270 auto&& solution = _ip_data[ip].solid_material.integrateStress(
271 variables_prev, variables, t, x_position, dt, *state);
272
273 if (!solution)
274 {
275 OGS_FATAL("Computation of local constitutive relation failed.");
276 }
277
279 std::tie(sigma, state, C) = std::move(*solution);
280
281 local_Jac
282 .template block<displacement_size, displacement_size>(
283 displacement_index, displacement_index)
284 .noalias() += B.transpose() * C * B * w;
285
286 typename ShapeMatricesType::template MatrixType<DisplacementDim,
287 displacement_size>
288 N_u = ShapeMatricesType::template MatrixType<
289 DisplacementDim, displacement_size>::Zero(DisplacementDim,
290 displacement_size);
291
292 for (int i = 0; i < DisplacementDim; ++i)
293 {
294 N_u.template block<1, displacement_size / DisplacementDim>(
295 i, i * displacement_size / DisplacementDim)
296 .noalias() = N;
297 }
298
299 auto const rho_s =
300 solid_phase.property(MPL::PropertyType::density)
301 .template value<double>(variables, x_position, t, dt);
302
303 auto const& b = _process_data.specific_body_force;
304 local_rhs.template block<displacement_size, 1>(displacement_index, 0)
305 .noalias() -=
306 (B.transpose() * sigma - N_u.transpose() * rho_s * b) * w;
307
308 //
309 // displacement equation, temperature part
310 // The computation of KuT can be ignored.
311 auto const alpha_T_tensor =
313 solid_linear_thermal_expansivity_vector);
314 KuT.noalias() += B.transpose() * (C * alpha_T_tensor) * N * w;
315
316 if (_ip_data[ip].solid_material.getConstitutiveModel() ==
318 {
319 auto const s = Invariants::deviatoric_projection * sigma;
320 double const norm_s = Invariants::FrobeniusNorm(s);
321 const double creep_coefficient =
322 _ip_data[ip].solid_material.getTemperatureRelatedCoefficient(
323 t, dt, x_position, T_ip, norm_s);
324 KuT.noalias() += creep_coefficient * B.transpose() * s * N * w;
325 }
326
327 //
328 // temperature equation, temperature part;
329 //
330 auto const lambda =
331 solid_phase
332 .property(
334 .value(variables, x_position, t, dt);
335
336 GlobalDimMatrixType const thermal_conductivity =
338
339 KTT.noalias() += dNdx.transpose() * thermal_conductivity * dNdx * w;
340
341 auto const c =
342 solid_phase
343 .property(
345 .template value<double>(variables, x_position, t, dt);
346 DTT.noalias() += N.transpose() * rho_s * c * N * w;
347 }
348
349 // temperature equation, temperature part
350 local_Jac
351 .template block<temperature_size, temperature_size>(temperature_index,
352 temperature_index)
353 .noalias() += KTT + DTT / dt;
354
355 // displacement equation, temperature part
356 local_Jac
357 .template block<displacement_size, temperature_size>(displacement_index,
358 temperature_index)
359 .noalias() -= KuT;
360
361 local_rhs.template block<temperature_size, 1>(temperature_index, 0)
362 .noalias() -= KTT * T + DTT * (T - T_prev) / dt;
363}
364
365template <typename ShapeFunction, int DisplacementDim>
367 assembleWithJacobianForStaggeredScheme(const double t, double const dt,
368 Eigen::VectorXd const& local_x,
369 Eigen::VectorXd const& local_x_prev,
370 int const process_id,
371 std::vector<double>& local_b_data,
372 std::vector<double>& local_Jac_data)
373{
374 // For the equations with pressure
375 if (process_id == _process_data.heat_conduction_process_id)
376 {
377 assembleWithJacobianForHeatConductionEquations(
378 t, dt, local_x, local_x_prev, local_b_data, local_Jac_data);
379 return;
380 }
381
382 // For the equations with deformation
383 assembleWithJacobianForDeformationEquations(t, dt, local_x, local_x_prev,
384 local_b_data, local_Jac_data);
385}
386
387template <typename ShapeFunction, int DisplacementDim>
390 const double t, double const dt, Eigen::VectorXd const& local_x,
391 Eigen::VectorXd const& local_x_prev, std::vector<double>& local_b_data,
392 std::vector<double>& local_Jac_data)
393{
394 auto const local_T =
395 local_x.template segment<temperature_size>(temperature_index);
396
397 auto const local_T_prev =
398 local_x_prev.template segment<temperature_size>(temperature_index);
399
400 auto const local_u =
401 local_x.template segment<displacement_size>(displacement_index);
402 bool const is_u_non_zero = local_u.norm() > 0.0;
403
404 auto local_Jac = MathLib::createZeroedMatrix<
405 typename ShapeMatricesType::template MatrixType<displacement_size,
406 displacement_size>>(
407 local_Jac_data, displacement_size, displacement_size);
408
409 auto local_rhs = MathLib::createZeroedVector<
410 typename ShapeMatricesType::template VectorType<displacement_size>>(
411 local_b_data, displacement_size);
412
413 unsigned const n_integration_points =
414 _integration_method.getNumberOfPoints();
415
416 MPL::VariableArray variables;
417 MPL::VariableArray variables_prev;
418 auto const& medium = _process_data.media_map.getMedium(_element.getID());
419 auto const& solid_phase = medium->phase("Solid");
420
421 for (unsigned ip = 0; ip < n_integration_points; ip++)
422 {
423 auto const& w = _ip_data[ip].integration_weight;
424 auto const& N = _ip_data[ip].N;
425 auto const& dNdx = _ip_data[ip].dNdx;
426
427 ParameterLib::SpatialPosition const x_position{
428 std::nullopt, _element.getID(),
431 _element, N))};
432
433 auto const x_coord = x_position.getCoordinates().value()[0];
434
435 auto const& B =
436 LinearBMatrix::computeBMatrix<DisplacementDim,
437 ShapeFunction::NPOINTS,
439 dNdx, N, x_coord, _is_axially_symmetric);
440
441 auto& sigma = _ip_data[ip].sigma;
442 auto const& sigma_prev = _ip_data[ip].sigma_prev;
443
444 auto& eps = _ip_data[ip].eps;
445 auto const& eps_prev = _ip_data[ip].eps_prev;
446
447 auto& eps_m = _ip_data[ip].eps_m;
448 auto const& eps_m_prev = _ip_data[ip].eps_m_prev;
449
450 auto& state = _ip_data[ip].material_state_variables;
451
452 const double T_ip = N.dot(local_T); // T at integration point
453 variables.temperature = T_ip;
454 double const dT_ip = T_ip - N.dot(local_T_prev);
455
456 //
457 // displacement equation, displacement part
458 //
459 // For the restart computation, the displacement may not be
460 // reloaded but the initial strains are always available. For such case,
461 // the following computation is skipped.
462 if (is_u_non_zero)
463 {
464 eps.noalias() = B * local_u;
465 }
466
467 // Consider also anisotropic thermal expansion.
468 auto const solid_linear_thermal_expansivity_vector =
470 solid_phase
471 .property(
473 .value(variables, x_position, t, dt));
474
476 dthermal_strain = solid_linear_thermal_expansivity_vector * dT_ip;
477
478 eps_m.noalias() = eps_m_prev + eps - eps_prev - dthermal_strain;
479
480 variables_prev.stress
482 sigma_prev);
483 variables_prev.mechanical_strain
485 eps_m_prev);
486 variables_prev.temperature = T_ip;
487 variables.mechanical_strain
489 eps_m);
490
491 auto&& solution = _ip_data[ip].solid_material.integrateStress(
492 variables_prev, variables, t, x_position, dt, *state);
493
494 if (!solution)
495 {
496 OGS_FATAL("Computation of local constitutive relation failed.");
497 }
498
500 std::tie(sigma, state, C) = std::move(*solution);
501
502 local_Jac.noalias() += B.transpose() * C * B * w;
503
504 typename ShapeMatricesType::template MatrixType<DisplacementDim,
505 displacement_size>
506 N_u = ShapeMatricesType::template MatrixType<
507 DisplacementDim, displacement_size>::Zero(DisplacementDim,
508 displacement_size);
509
510 for (int i = 0; i < DisplacementDim; ++i)
511 {
512 N_u.template block<1, displacement_size / DisplacementDim>(
513 i, i * displacement_size / DisplacementDim)
514 .noalias() = N;
515 }
516
517 auto const rho_s =
518 solid_phase.property(MPL::PropertyType::density)
519 .template value<double>(variables, x_position, t, dt);
520
521 auto const& b = _process_data.specific_body_force;
522 local_rhs.noalias() -=
523 (B.transpose() * sigma - N_u.transpose() * rho_s * b) * w;
524 }
525}
526
527template <typename ShapeFunction, int DisplacementDim>
530 const double t, double const dt, Eigen::VectorXd const& local_x,
531 Eigen::VectorXd const& local_x_prev, std::vector<double>& local_b_data,
532 std::vector<double>& local_Jac_data)
533{
534 auto const local_T =
535 local_x.template segment<temperature_size>(temperature_index);
536
537 auto const local_T_prev =
538 local_x_prev.template segment<temperature_size>(temperature_index);
539
540 auto local_Jac = MathLib::createZeroedMatrix<
541 typename ShapeMatricesType::template MatrixType<temperature_size,
542 temperature_size>>(
543 local_Jac_data, temperature_size, temperature_size);
544
545 auto local_rhs = MathLib::createZeroedVector<
546 typename ShapeMatricesType::template VectorType<temperature_size>>(
547 local_b_data, temperature_size);
548
550 mass.setZero(temperature_size, temperature_size);
551
552 typename ShapeMatricesType::NodalMatrixType laplace;
553 laplace.setZero(temperature_size, temperature_size);
554
555 unsigned const n_integration_points =
556 _integration_method.getNumberOfPoints();
557
558 auto const& medium = _process_data.media_map.getMedium(_element.getID());
559 auto const& solid_phase = medium->phase("Solid");
560 MPL::VariableArray variables;
561
562 for (unsigned ip = 0; ip < n_integration_points; ip++)
563 {
564 auto const& w = _ip_data[ip].integration_weight;
565 auto const& N = _ip_data[ip].N;
566 auto const& dNdx = _ip_data[ip].dNdx;
567
568 ParameterLib::SpatialPosition const x_position{
569 std::nullopt, _element.getID(),
572 _element, N))};
573
574 const double T_ip = N.dot(local_T); // T at integration point
575 variables.temperature = T_ip;
576
577 auto const rho_s =
578 solid_phase.property(MPL::PropertyType::density)
579 .template value<double>(variables, x_position, t, dt);
580 auto const c_p =
581 solid_phase.property(MPL::PropertyType::specific_heat_capacity)
582 .template value<double>(variables, x_position, t, dt);
583
584 mass.noalias() += N.transpose() * rho_s * c_p * N * w;
585
586 auto const lambda =
587 solid_phase
588 .property(
590 .value(variables, x_position, t, dt);
591
592 GlobalDimMatrixType const thermal_conductivity =
594
595 laplace.noalias() += dNdx.transpose() * thermal_conductivity * dNdx * w;
596 }
597 local_Jac.noalias() += laplace + mass / dt;
598
599 local_rhs.noalias() -=
600 laplace * local_T + mass * (local_T - local_T_prev) / dt;
601}
602
603template <typename ShapeFunction, int DisplacementDim>
604std::size_t
611
612template <typename ShapeFunction, int DisplacementDim>
613std::vector<double>
615{
616 constexpr int kelvin_vector_size =
618
620 [this](std::vector<double>& values)
621 { return getIntPtSigma(0, {}, {}, values); });
622}
623
624template <typename ShapeFunction, int DisplacementDim>
625std::vector<double> const&
627 const double /*t*/,
628 std::vector<GlobalVector*> const& /*x*/,
629 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
630 std::vector<double>& cache) const
631{
633 _ip_data, &IpData::sigma, cache);
634}
635
636template <typename ShapeFunction, int DisplacementDim>
637std::size_t
644
645template <typename ShapeFunction, int DisplacementDim>
646std::vector<double> ThermoMechanicsLocalAssembler<
647 ShapeFunction, DisplacementDim>::getEpsilon() const
648{
649 constexpr int kelvin_vector_size =
651
653 [this](std::vector<double>& values)
654 { return getIntPtEpsilon(0, {}, {}, values); });
655}
656
657template <typename ShapeFunction, int DisplacementDim>
658std::vector<double> const&
660 const double /*t*/,
661 std::vector<GlobalVector*> const& /*x*/,
662 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
663 std::vector<double>& cache) const
664{
666 _ip_data, &IpData::eps, cache);
667}
668
669template <typename ShapeFunction, int DisplacementDim>
670std::vector<double> const&
673 const double /*t*/,
674 std::vector<GlobalVector*> const& /*x*/,
675 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
676 std::vector<double>& cache) const
677{
679 _ip_data, &IpData::eps_m, cache);
680}
681
682template <typename ShapeFunction, int DisplacementDim>
684 ShapeFunction, DisplacementDim>::setEpsilonMechanical(double const* values)
685{
687 values, _ip_data, &IpData::eps_m);
688}
689template <typename ShapeFunction, int DisplacementDim>
690std::vector<double> ThermoMechanicsLocalAssembler<
691 ShapeFunction, DisplacementDim>::getEpsilonMechanical() const
692{
693 constexpr int kelvin_vector_size =
695
697 [this](std::vector<double>& values)
698 { return getIntPtEpsilonMechanical(0, {}, {}, values); });
699}
700
701template <typename ShapeFunction, int DisplacementDim>
703 ShapeFunction, DisplacementDim>::getNumberOfIntegrationPoints() const
704{
705 return _integration_method.getNumberOfPoints();
706}
707
708template <typename ShapeFunction, int DisplacementDim>
709int ThermoMechanicsLocalAssembler<ShapeFunction,
710 DisplacementDim>::getMaterialID() const
711{
712 return _process_data.material_ids == nullptr
713 ? 0
714 : (*_process_data.material_ids)[_element.getID()];
715}
716
717template <typename ShapeFunction, int DisplacementDim>
719 DisplacementDim>::MaterialStateVariables const&
721 getMaterialStateVariablesAt(unsigned integration_point) const
722{
723 return *_ip_data[integration_point].material_state_variables;
724}
725} // namespace ThermoMechanics
726} // namespace ProcessLib
#define OGS_FATAL(...)
Definition Error.h:26
constexpr double getWeight() const
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:91
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
void setElementID(std::size_t element_id)
std::optional< MathLib::Point3d > const getCoordinates() const
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
std::vector< double > const & getIntPtEpsilonMechanical(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
std::size_t setIPDataInitialConditions(std::string_view const name, double const *values, int const integration_order) override
Returns number of read integration points.
std::vector< double > const & getIntPtSigma(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
ThermoMechanicsLocalAssembler(ThermoMechanicsLocalAssembler const &)=delete
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(unsigned integration_point) const override
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
ThermoMechanicsProcessData< DisplacementDim > & _process_data
NumLib::GenericIntegrationMethod const & _integration_method
void assembleWithJacobianForHeatConductionEquations(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
void assembleWithJacobianForDeformationEquations(const double t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
void assembleWithJacobianForStaggeredScheme(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_x_prev, int const process_id, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
typename ShapeMatricesType::GlobalDimMatrixType GlobalDimMatrixType
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
void setInitialConditionsConcrete(Eigen::VectorXd const local_x, double const t, int const process_id) override
std::vector< double > const & getIntPtEpsilon(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
auto & selectSolidConstitutiveRelation(SolidMaterialsMap const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
MathLib::KelvinVector::KelvinVectorType< GlobalDim > formKelvinVector(MaterialPropertyLib::PropertyDataType const &values)
A function to form a Kelvin vector from strain or stress alike property like thermal expansivity for ...
Eigen::Matrix< double, GlobalDim, GlobalDim > formEigenTensor(MaterialPropertyLib::PropertyDataType const &values)
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
Eigen::Map< Matrix > createZeroedMatrix(std::vector< double > &data, Eigen::MatrixXd::Index rows, Eigen::MatrixXd::Index cols)
std::vector< typename ShapeMatricesType::ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > initShapeMatrices(MeshLib::Element const &e, bool const is_axially_symmetric, IntegrationMethod const &integration_method)
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
BMatrixType computeBMatrix(DNDX_Type const &dNdx, N_Type const &N, const double radius, const bool is_axially_symmetric)
Fills a B-matrix based on given shape function dN/dx values.
std::vector< double > transposeInPlace(StoreValuesFunction const &store_values_function)
std::vector< double > const & getIntegrationPointKelvinVectorData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
std::size_t setIntegrationPointKelvinVectorData(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType IpData::*const member)
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N