OGS
ThermoRichardsMechanicsFEM-impl.h
Go to the documentation of this file.
1
12#pragma once
13
14#include <spdlog/fmt/bundled/format.h>
15
16#include <cassert>
17#include <type_traits>
18
31
32namespace ProcessLib
33{
34namespace ThermoRichardsMechanics
35{
36template <typename ShapeFunctionDisplacement, typename ShapeFunction,
37 int DisplacementDim, typename ConstitutiveTraits>
38ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement, ShapeFunction,
39 DisplacementDim, ConstitutiveTraits>::
40 ThermoRichardsMechanicsLocalAssembler(
41 MeshLib::Element const& e,
42 std::size_t const /*local_matrix_size*/,
43 NumLib::GenericIntegrationMethod const& integration_method,
44 bool const is_axially_symmetric,
46 process_data)
47 : LocalAssemblerInterface<DisplacementDim, ConstitutiveTraits>(
48 e, integration_method, is_axially_symmetric, process_data)
49{
50 unsigned const n_integration_points =
52
53 ip_data_.resize(n_integration_points);
54
55 auto const shape_matrices_u =
56 NumLib::initShapeMatrices<ShapeFunctionDisplacement,
58 DisplacementDim>(e, is_axially_symmetric,
59 integration_method);
60
61 auto const shape_matrices =
63 DisplacementDim>(e, is_axially_symmetric,
64 integration_method);
65
66 for (unsigned ip = 0; ip < n_integration_points; ip++)
67 {
68 auto& ip_data = ip_data_[ip];
69 auto const& sm_u = shape_matrices_u[ip];
70 ip_data_[ip].integration_weight =
71 integration_method.getWeightedPoint(ip).getWeight() *
72 sm_u.integralMeasure * sm_u.detJ;
73
74 ip_data.N_u = sm_u.N;
75 ip_data.dNdx_u = sm_u.dNdx;
76
77 // ip_data.N_p and ip_data.dNdx_p are used for both p and T variables
78 ip_data.N_p = shape_matrices[ip].N;
79 ip_data.dNdx_p = shape_matrices[ip].dNdx;
80 }
81}
82
83template <typename ShapeFunctionDisplacement, typename ShapeFunction,
84 int DisplacementDim, typename ConstitutiveTraits>
86 ShapeFunctionDisplacement, ShapeFunction, DisplacementDim,
87 ConstitutiveTraits>::setInitialConditionsConcrete(Eigen::VectorXd const
88 local_x,
89 double const t,
90 int const /*process_id*/)
91{
92 assert(local_x.size() ==
93 temperature_size + pressure_size + displacement_size);
94
95 auto const p_L = local_x.template segment<pressure_size>(pressure_index);
96 auto const T =
97 local_x.template segment<temperature_size>(temperature_index);
98
99 constexpr double dt = std::numeric_limits<double>::quiet_NaN();
100 auto const& medium =
101 *this->process_data_.media_map.getMedium(this->element_.getID());
102
103 MediaData const media_data{medium};
104
105 typename ConstitutiveTraits::ConstitutiveSetting const constitutive_setting;
106 auto models = ConstitutiveTraits::createConstitutiveModels(
107 this->process_data_, this->solid_material_);
108
109 unsigned const n_integration_points =
110 this->integration_method_.getNumberOfPoints();
111 for (unsigned ip = 0; ip < n_integration_points; ip++)
112 {
113 // N is used for both T and p variables.
114 auto const& N = ip_data_[ip].N_p;
115
116 ParameterLib::SpatialPosition const x_position{
117 std::nullopt, this->element_.getID(), ip,
119 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
121 this->element_, ip_data_[ip].N_u))};
122
123 double T_ip;
125
126 double p_cap_ip;
127 NumLib::shapeFunctionInterpolate(-p_L, N, p_cap_ip);
128
129 MPL::VariableArray variables;
130 variables.capillary_pressure = p_cap_ip;
131 variables.liquid_phase_pressure = -p_cap_ip;
132 // setting pG to 1 atm
133 // TODO : rewrite equations s.t. p_L = pG-p_cap
134 variables.gas_phase_pressure = 1.0e5;
135 variables.temperature = T_ip;
136
137 double const S_L =
138 medium.property(MPL::PropertyType::saturation)
139 .template value<double>(variables, x_position, t, dt);
140 std::get<PrevState<SaturationData>>(this->prev_states_[ip])->S_L = S_L;
141
142 constitutive_setting.init(models, t, dt, x_position, media_data,
143 {T_ip, 0, {}}, this->current_states_[ip],
144 this->prev_states_[ip]);
145
146 if (this->process_data_.initial_stress.value)
147 {
148 variables.liquid_saturation = S_L;
149 convertInitialStressType(ip, t, x_position, medium, variables,
150 -p_cap_ip);
151 }
152 }
153}
154
155template <typename ShapeFunctionDisplacement, typename ShapeFunction,
156 int DisplacementDim, typename ConstitutiveTraits>
157void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
158 ShapeFunction, DisplacementDim,
159 ConstitutiveTraits>::
160 convertInitialStressType(unsigned const ip,
161 double const t,
162 ParameterLib::SpatialPosition const x_position,
163 MaterialPropertyLib::Medium const& medium,
164 MPL::VariableArray const& variables,
165 double const p_at_ip)
166{
167 bool constexpr is_strain_temperature_constitutive =
169 DisplacementDim>,
170 ConstitutiveTraits>::value;
171 if (is_strain_temperature_constitutive &&
172 this->process_data_.initial_stress.type ==
174 {
175 return;
176 }
177
178 if (!is_strain_temperature_constitutive &&
179 this->process_data_.initial_stress.type == InitialStress::Type::Total)
180 {
181 return;
182 }
183
184 double const alpha_b =
185 medium.property(MPL::PropertyType::biot_coefficient)
186 .template value<double>(variables, x_position, t, 0.0 /*dt*/);
187
188 double const bishop =
189 medium.property(MPL::PropertyType::bishops_effective_stress)
190 .template value<double>(variables, x_position, t, 0.0 /*dt*/);
191
192 ConstitutiveTraits::ConstitutiveSetting::convertInitialStressType(
193 this->current_states_[ip], this->prev_states_[ip],
194 bishop * alpha_b * p_at_ip * Invariants::identity2);
195}
196
197template <typename ShapeFunctionDisplacement, typename ShapeFunction,
198 int DisplacementDim, typename ConstitutiveTraits>
199void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
200 ShapeFunction, DisplacementDim,
201 ConstitutiveTraits>::
202 assembleWithJacobian(double const t, double const dt,
203 std::vector<double> const& local_x,
204 std::vector<double> const& local_x_prev,
205 std::vector<double>& /*local_M_data*/,
206 std::vector<double>& /*local_K_data*/,
207 std::vector<double>& local_rhs_data,
208 std::vector<double>& local_Jac_data)
209{
210 auto& medium =
211 *this->process_data_.media_map.getMedium(this->element_.getID());
212
213 LocalMatrices loc_mat;
214 loc_mat.setZero();
215 LocalMatrices loc_mat_current_ip;
216 loc_mat_current_ip.setZero(); // only to set the right matrix sizes
217
218 typename ConstitutiveTraits::ConstitutiveSetting constitutive_setting;
219
220 for (unsigned ip = 0; ip < this->integration_method_.getNumberOfPoints();
221 ++ip)
222 {
223 ParameterLib::SpatialPosition const x_position{
224 std::nullopt, this->element_.getID(), ip,
226 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
228 this->element_, ip_data_[ip].N_u))};
229
230 assembleWithJacobianSingleIP(
231 t, dt, x_position, //
232 local_x, local_x_prev, //
233 ip_data_[ip], constitutive_setting,
234 medium, //
235 loc_mat_current_ip, //
236 this->current_states_[ip], this->prev_states_[ip],
237 this->material_states_[ip], this->output_data_[ip]);
238 loc_mat += loc_mat_current_ip;
239 }
240
241 massLumping(loc_mat);
242
243 addToLocalMatrixData(dt, local_x, local_x_prev, loc_mat, local_rhs_data,
244 local_Jac_data);
245}
246
247template <typename ShapeFunctionDisplacement, typename ShapeFunction,
248 int DisplacementDim, typename ConstitutiveTraits>
249void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
250 ShapeFunction, DisplacementDim,
251 ConstitutiveTraits>::
252 massLumping(typename ThermoRichardsMechanicsLocalAssembler<
253 ShapeFunctionDisplacement, ShapeFunction, DisplacementDim,
254 ConstitutiveTraits>::LocalMatrices& loc_mat) const
255{
256 if (this->process_data_.apply_mass_lumping)
257 {
258 loc_mat.storage_p_a_p =
259 loc_mat.storage_p_a_p.colwise().sum().eval().asDiagonal();
260 loc_mat.storage_p_a_S =
261 loc_mat.storage_p_a_S.colwise().sum().eval().asDiagonal();
262 loc_mat.storage_p_a_S_Jpp =
263 loc_mat.storage_p_a_S_Jpp.colwise().sum().eval().asDiagonal();
264 }
265}
266
267template <typename ShapeFunctionDisplacement, typename ShapeFunction,
268 int DisplacementDim, typename ConstitutiveTraits>
269void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
270 ShapeFunction, DisplacementDim,
271 ConstitutiveTraits>::
272 addToLocalMatrixData(
273 double const dt,
274 std::vector<double> const& local_x,
275 std::vector<double> const& local_x_prev,
277 ShapeFunctionDisplacement, ShapeFunction, DisplacementDim,
278 ConstitutiveTraits>::LocalMatrices const& loc_mat,
279 std::vector<double>& local_rhs_data,
280 std::vector<double>& local_Jac_data) const
281{
282 constexpr auto local_matrix_dim =
283 displacement_size + pressure_size + temperature_size;
284 assert(local_x.size() == local_matrix_dim);
285
286 auto local_Jac = MathLib::createZeroedMatrix<
287 typename ShapeMatricesTypeDisplacement::template MatrixType<
288 local_matrix_dim, local_matrix_dim>>(
289 local_Jac_data, local_matrix_dim, local_matrix_dim);
290
291 auto local_rhs =
292 MathLib::createZeroedVector<typename ShapeMatricesTypeDisplacement::
293 template VectorType<local_matrix_dim>>(
294 local_rhs_data, local_matrix_dim);
295
296 local_Jac.noalias() = loc_mat.Jac;
297 local_rhs.noalias() = -loc_mat.res;
298
299 //
300 // -- Jacobian
301 //
302 block_TT(local_Jac).noalias() += loc_mat.M_TT / dt + loc_mat.K_TT;
303 block_Tp(local_Jac).noalias() +=
304 loc_mat.K_Tp + loc_mat.dK_TT_dp + loc_mat.M_Tp / dt;
305
306 block_pT(local_Jac).noalias() += loc_mat.M_pT / dt + loc_mat.K_pT;
307 block_pp(local_Jac).noalias() +=
308 loc_mat.K_pp + loc_mat.storage_p_a_p / dt + loc_mat.storage_p_a_S_Jpp;
309 block_pu(local_Jac).noalias() = loc_mat.M_pu / dt;
310
311 //
312 // -- Residual
313 //
314 auto const [T, p_L, u] = localDOF(local_x);
315 auto const [T_prev, p_L_prev, u_prev] = localDOF(local_x_prev);
316
317 block_T(local_rhs).noalias() -= loc_mat.M_TT * (T - T_prev) / dt +
318 loc_mat.K_TT * T + loc_mat.K_Tp * p_L +
319 loc_mat.M_Tp * (p_L - p_L_prev) / dt;
320 block_p(local_rhs).noalias() -=
321 loc_mat.K_pp * p_L + loc_mat.K_pT * T +
322 (loc_mat.storage_p_a_p + loc_mat.storage_p_a_S) * (p_L - p_L_prev) /
323 dt +
324 loc_mat.M_pu * (u - u_prev) / dt + loc_mat.M_pT * (T - T_prev) / dt;
325}
326
327template <typename ShapeFunctionDisplacement, typename ShapeFunction,
328 int DisplacementDim, typename ConstitutiveTraits>
329void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
330 ShapeFunction, DisplacementDim,
331 ConstitutiveTraits>::
332 assembleWithJacobianSingleIP(
333 double const t, double const dt,
334 ParameterLib::SpatialPosition const& x_position,
335 std::vector<double> const& local_x,
336 std::vector<double> const& local_x_prev,
338 ShapeFunctionDisplacement, ShapeFunction, DisplacementDim,
339 ConstitutiveTraits>::IpData const& ip_data,
340 typename ConstitutiveTraits::ConstitutiveSetting& CS,
343 ShapeFunctionDisplacement, ShapeFunction, DisplacementDim,
344 ConstitutiveTraits>::LocalMatrices& out,
345 typename ConstitutiveTraits::StatefulData& current_state,
346 typename ConstitutiveTraits::StatefulDataPrev const& prev_state,
348 typename ConstitutiveTraits::OutputData& output_data) const
349{
350 auto const& N_u = ip_data.N_u;
351 auto const& dNdx_u = ip_data.dNdx_u;
352
353 // N and dNdx are used for both p and T variables
354 auto const& N = ip_data.N_p;
355 auto const& dNdx = ip_data.dNdx_p;
356
357 auto const B =
358 LinearBMatrix::computeBMatrix<DisplacementDim,
359 ShapeFunctionDisplacement::NPOINTS,
361 dNdx_u, N_u, (*x_position.getCoordinates())[0],
362 this->is_axially_symmetric_);
363
364 typename ConstitutiveTraits::ConstitutiveData CD;
365
366 auto const [T, p_L, u] = localDOF(local_x);
367 auto const [T_prev, p_L_prev, u_prev] = localDOF(local_x_prev);
368
369 {
370 auto models = ConstitutiveTraits::createConstitutiveModels(
371 this->process_data_, this->solid_material_);
372 typename ConstitutiveTraits::ConstitutiveTempData tmp;
373
374 double const T_ip = N * T;
375 double const T_prev_ip = N * T_prev;
376 GlobalDimVectorType const grad_T_ip = dNdx * T;
377
378 double const p_cap_ip = -N * p_L;
379 double const p_cap_prev_ip = -N * p_L_prev;
380 GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
381
382 KelvinVectorType eps = B * u;
383
384 CS.eval(models, t, dt, x_position, //
385 medium, //
386 {T_ip, T_prev_ip, grad_T_ip}, //
387 {p_cap_ip, p_cap_prev_ip, grad_p_cap_ip}, //
388 eps, current_state, prev_state, mat_state, tmp, output_data,
389 CD);
390 }
391
392 using NodalMatrix = typename ShapeMatricesType::NodalMatrixType;
393 NodalMatrix const NTN = N.transpose() * N;
394 NodalMatrix const dNTdN = dNdx.transpose() * dNdx;
395
396 // TODO is identity2.transpose() * B something like divergence?
397 auto const& identity2 = MathLib::KelvinVector::Invariants<
399 DisplacementDim)>::identity2;
400 typename ShapeMatricesTypeDisplacement::template MatrixType<
401 displacement_size, pressure_size>
402 BTI2N = B.transpose() * identity2 * N;
403
404 /*
405 * Conventions:
406 *
407 * * use positive signs exclusively, any negative signs must be included in
408 * the coefficients coming from the constitutive setting
409 * * the used coefficients from the constitutive setting are named after the
410 * terms they appear in, e.g. K_TT_X_dNTdN means it is a scalar (X) that
411 * will be multiplied by dNdx.transpose() * dNdx. Placefolders for the
412 * coefficients are:
413 * * X -> scalar
414 * * V -> vector
415 * * K -> Kelvin vector
416 * * the Laplace terms have a special name, e.g., K_TT_Laplace
417 * * there shall be only one contribution to each of the LocalMatrices,
418 * assigned with = assignment; this point might be relaxed in the future
419 * * this method will overwrite the data in the passed LocalMatrices& out
420 * argument, not add to it
421 */
422
423 // residual, order T, p, u
424 block_p(out.res).noalias() =
425 dNdx.transpose() * std::get<EqPData<DisplacementDim>>(CD).rhs_p_dNT_V;
426 block_u(out.res).noalias() =
427 B.transpose() *
428 ProcessLib::Graph::get<TotalStressData<DisplacementDim>>(
429 CD, current_state)
430 .sigma_total -
431 static_cast<int>(this->process_data_.apply_body_force_for_deformation) *
432 N_u_op(N_u).transpose() *
433 std::get<GravityData<DisplacementDim>>(CD).volumetric_body_force;
434
435 // Storage matrices
436 out.storage_p_a_p.noalias() =
437 std::get<EqPData<DisplacementDim>>(CD).storage_p_a_p_X_NTN * NTN;
438 out.storage_p_a_S.noalias() =
439 std::get<TRMStorageData>(CD).storage_p_a_S_X_NTN * NTN;
440 out.storage_p_a_S_Jpp.noalias() =
441 std::get<TRMStorageData>(CD).storage_p_a_S_Jpp_X_NTN * NTN;
442
443 // M matrices, order T, p, u
444 out.M_TT.noalias() =
445 std::get<EqTData<DisplacementDim>>(CD).M_TT_X_NTN * NTN;
446 out.M_Tp.noalias() =
447 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).M_Tp_X_NTN * NTN;
448
449 out.M_pT.noalias() =
450 std::get<EqPData<DisplacementDim>>(CD).M_pT_X_NTN * NTN;
451 out.M_pu.noalias() =
452 std::get<EqPData<DisplacementDim>>(CD).M_pu_X_BTI2N * BTI2N.transpose();
453
454 // K matrices, order T, p, u
455 out.K_TT.noalias() =
456 dNdx.transpose() *
457 std::get<TRMHeatStorageAndFluxData<DisplacementDim>>(CD)
458 .K_TT_Laplace *
459 dNdx +
460 N.transpose() *
461 (std::get<EqTData<DisplacementDim>>(CD).K_TT_NT_V_dN.transpose() *
462 dNdx) +
464 dNTdN;
465
466 out.dK_TT_dp.noalias() =
467 N.transpose() *
468 (std::get<TRMHeatStorageAndFluxData<DisplacementDim>>(CD)
469 .K_Tp_NT_V_dN.transpose() *
470 dNdx) +
472 NTN;
473 out.K_Tp.noalias() =
474 dNdx.transpose() *
475 std::get<ThermoOsmosisData<DisplacementDim>>(CD).K_Tp_Laplace *
476 dNdx +
478 dNTdN;
479
480 out.K_pp.noalias() =
481 dNdx.transpose() * std::get<EqPData<DisplacementDim>>(CD).K_pp_Laplace *
482 dNdx +
484 dNTdN;
485 out.K_pT.noalias() =
486 dNdx.transpose() *
487 std::get<ThermoOsmosisData<DisplacementDim>>(CD).K_pT_Laplace * dNdx;
488
489 // direct Jacobian contributions, order T, p, u
490 block_pT(out.Jac).noalias() =
491 std::get<TRMVaporDiffusionData<DisplacementDim>>(CD).J_pT_X_dNTdN *
492 dNTdN;
493 block_pp(out.Jac).noalias() =
494 std::get<TRMStorageData>(CD).J_pp_X_NTN * NTN +
495 std::get<EqPData<DisplacementDim>>(CD).J_pp_X_BTI2NT_u_dot_N *
496 BTI2N.transpose() * (u - u_prev) / dt *
497 N // TODO something with volumetric strain rate?
498 + dNdx.transpose() *
499 std::get<EqPData<DisplacementDim>>(CD).J_pp_dNT_V_N * N;
500
501 block_uT(out.Jac).noalias() =
502 B.transpose() *
503 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD).J_uT_BT_K_N *
504 N;
505 block_up(out.Jac).noalias() =
506 B.transpose() *
507 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD)
508 .J_up_BT_K_N *
509 N +
510 N_u_op(N_u).transpose() *
511 std::get<GravityData<DisplacementDim>>(CD).J_up_HT_V_N * N;
512 block_uu(out.Jac).noalias() =
513 B.transpose() *
514 std::get<SolidMechanicsDataStateless<DisplacementDim>>(CD)
515 .stiffness_tensor *
516 B;
517
518 out *= ip_data.integration_weight;
519}
520
521template <typename ShapeFunctionDisplacement, typename ShapeFunction,
522 int DisplacementDim, typename ConstitutiveTraits>
523void ThermoRichardsMechanicsLocalAssembler<ShapeFunctionDisplacement,
524 ShapeFunction, DisplacementDim,
525 ConstitutiveTraits>::
526 computeSecondaryVariableConcrete(double const t, double const dt,
527 Eigen::VectorXd const& local_x,
528 Eigen::VectorXd const& local_x_prev)
529{
530 auto const T = block_T(local_x);
531 auto const p_L = block_p(local_x);
532 auto const u = block_u(local_x);
533
534 auto const T_prev = block_T(local_x_prev);
535 auto const p_L_prev = block_p(local_x_prev);
536
537 auto const e_id = this->element_.getID();
538 auto const& process_data = this->process_data_;
539 auto& medium = *process_data.media_map.getMedium(e_id);
540
541 unsigned const n_integration_points =
542 this->integration_method_.getNumberOfPoints();
543
544 typename ConstitutiveTraits::ConstitutiveSetting constitutive_setting;
545
546 auto models = ConstitutiveTraits::createConstitutiveModels(
547 process_data, this->solid_material_);
548 typename ConstitutiveTraits::ConstitutiveTempData tmp;
549 typename ConstitutiveTraits::ConstitutiveData CD;
550
551 for (unsigned ip = 0; ip < n_integration_points; ip++)
552 {
553 auto& current_state = this->current_states_[ip];
554 auto& output_data = this->output_data_[ip];
555
556 auto const& ip_data = ip_data_[ip];
557
558 // N is used for both p and T variables
559 auto const& N = ip_data.N_p;
560 auto const& N_u = ip_data.N_u;
561 auto const& dNdx_u = ip_data.dNdx_u;
562 auto const& dNdx = ip_data.dNdx_p;
563
564 ParameterLib::SpatialPosition const x_position{
565 std::nullopt, this->element_.getID(), ip,
567 NumLib::interpolateCoordinates<ShapeFunctionDisplacement,
569 this->element_, N_u))};
570 auto const x_coord =
571 NumLib::interpolateXCoordinate<ShapeFunctionDisplacement,
573 this->element_, N_u);
574 auto const B =
575 LinearBMatrix::computeBMatrix<DisplacementDim,
576 ShapeFunctionDisplacement::NPOINTS,
578 dNdx_u, N_u, x_coord, this->is_axially_symmetric_);
579
580 double const T_ip = N * T;
581 double const T_prev_ip = N * T_prev;
582 GlobalDimVectorType const grad_T_ip = dNdx * T;
583
584 double const p_cap_ip = -N * p_L;
585 double const p_cap_prev_ip = -N * p_L_prev;
586 GlobalDimVectorType const grad_p_cap_ip = -dNdx * p_L;
587
588 KelvinVectorType eps = B * u;
589
590 constitutive_setting.eval(models, //
591 t, dt, x_position, //
592 medium, //
593 {T_ip, T_prev_ip, grad_T_ip}, //
594 {p_cap_ip, p_cap_prev_ip, grad_p_cap_ip}, //
595 eps, current_state, this->prev_states_[ip],
596 this->material_states_[ip], tmp, output_data,
597 CD);
598 }
599
601 ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
602 DisplacementDim>(this->element_, this->is_axially_symmetric_, p_L,
603 *process_data.pressure_interpolated);
605 ShapeFunction, typename ShapeFunctionDisplacement::MeshElement,
606 DisplacementDim>(this->element_, this->is_axially_symmetric_, T,
607 *process_data.temperature_interpolated);
608}
609} // namespace ThermoRichardsMechanics
610} // namespace ProcessLib
Property const & property(PropertyType const &p) const
Definition Medium.cpp:54
double getWeight() const
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
std::optional< MathLib::Point3d > const & getCoordinates() const
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
ShapeMatrixPolicyType< ShapeFunctionDisplacement, DisplacementDim > ShapeMatricesTypeDisplacement
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
constexpr int kelvin_vector_dimensions(int const displacement_dim)
Kelvin vector dimensions for given displacement dimension.
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)
void shapeFunctionInterpolate(const NodalValues &, const ShapeMatrix &)
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
void interpolateToHigherOrderNodes(MeshLib::Element const &element, bool const is_axially_symmetric, Eigen::MatrixBase< EigenMatrixType > const &node_values, MeshLib::PropertyVector< double > &interpolated_values_global_vector)
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.
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType