Loading [MathJax]/extensions/MathMenu.js
OGS
SmallDeformationNonlocalFEM.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <algorithm>
14#include <limits>
15#include <memory>
16#include <vector>
17
18#include "Damage.h"
36
37namespace ProcessLib
38{
39namespace SmallDeformationNonlocal
40{
41namespace MPL = MaterialPropertyLib;
42
45template <typename ShapeMatrixType>
47{
48 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>> N;
49};
50
51template <typename ShapeFunction, int DisplacementDim>
54{
55public:
63
69 using IpData =
71
76
78 MeshLib::Element const& e,
79 std::size_t const /*local_matrix_size*/,
80 NumLib::GenericIntegrationMethod const& integration_method,
81 bool const is_axially_symmetric,
83 : _process_data(process_data),
84 _integration_method(integration_method),
85 _element(e),
86 _is_axially_symmetric(is_axially_symmetric)
87 {
88 unsigned const n_integration_points =
90
91 _ip_data.reserve(n_integration_points);
92 _secondary_data.N.resize(n_integration_points);
93
94 auto const shape_matrices =
96 DisplacementDim>(e, is_axially_symmetric,
98
99 auto& solid_material =
101 _process_data.solid_materials,
102 _process_data.material_ids,
103 e.getID());
104 auto* ehlers_solid_material = dynamic_cast<
106 &solid_material);
107 if (ehlers_solid_material == nullptr)
108 {
109 OGS_FATAL(
110 "The SmallDeformationNonlocalLocal process supports only "
111 "Ehlers material at the moment. For other materials the "
112 "interface must be extended first.");
113 }
114
115 for (unsigned ip = 0; ip < n_integration_points; ip++)
116 {
117 _ip_data.emplace_back(*ehlers_solid_material);
118 auto& ip_data = _ip_data[ip];
119 auto const& sm = shape_matrices[ip];
120 _ip_data[ip].integration_weight =
122 sm.integralMeasure * sm.detJ;
123
124 ip_data.N = sm.N;
125 ip_data.dNdx = sm.dNdx;
126
127 // Initialize current time step values
128 ip_data.sigma.setZero(
130 DisplacementDim));
132 DisplacementDim));
133
134 // Previous time step values are not initialized and are set later.
135 ip_data.sigma_prev.resize(
137 DisplacementDim));
138 ip_data.eps_prev.resize(
140 DisplacementDim));
141
142 _secondary_data.N[ip] = shape_matrices[ip].N;
143
144 ip_data.coordinates = getSingleIntegrationPointCoordinates(ip);
145 }
146 }
147
148 std::size_t setIPDataInitialConditions(std::string_view const name,
149 double const* values,
150 int const integration_order) override
151 {
152 if (integration_order !=
153 static_cast<int>(_integration_method.getIntegrationOrder()))
154 {
155 OGS_FATAL(
156 "Setting integration point initial conditions; The integration "
157 "order of the local assembler for element {:d} is different "
158 "from the integration order in the initial condition.",
159 _element.getID());
160 }
161
162 if (name == "sigma")
163 {
164 return setSigma(values);
165 }
166
167 if (name == "kappa_d")
168 {
171 }
172
173 return 0;
174 }
175
177 std::string const& name, std::vector<double> const& value) override
178 {
179 if (name == "kappa_d_ip")
180 {
181 if (value.size() != 1)
182 {
183 OGS_FATAL(
184 "CellData for kappa_d initial conditions has wrong number "
185 "of components. 1 expected, got {:d}.",
186 value.size());
187 }
188 setKappaD(value[0]);
189 }
190 }
191
192 double alpha_0(double const distance2) const
193 {
194 double const internal_length2 = _process_data.internal_length_squared;
195 return (distance2 > internal_length2)
196 ? 0
197 : (1 - distance2 / (internal_length2)) *
198 (1 - distance2 / (internal_length2));
199 }
200
202 std::size_t const /*mesh_item_id*/,
203 std::vector<
205 DisplacementDim>>> const& local_assemblers) override
206 {
207 auto const search_element_ids = MeshLib::findElementsWithinRadius(
208 _element, _process_data.internal_length_squared);
209
210 unsigned const n_integration_points =
212
213 std::vector<double> distances; // Cache for ip-ip distances.
214 //
215 // For every integration point in this element collect the neighbouring
216 // integration points falling in given radius (internal length) and
217 // compute the alpha_kl weight.
218 //
219 for (unsigned k = 0; k < n_integration_points; k++)
220 {
221 //
222 // Collect the integration points.
223 //
224
225 auto const& xyz = _ip_data[k].coordinates;
226
227 // For all neighbors of element
228 for (auto const search_element_id : search_element_ids)
229 {
230 auto const& la = local_assemblers[search_element_id];
231 la->getIntegrationPointCoordinates(xyz, distances);
232 for (int ip = 0; ip < static_cast<int>(distances.size()); ++ip)
233 {
234 if (distances[ip] >= _process_data.internal_length_squared)
235 {
236 continue;
237 }
238 // save into current ip_k
239 _ip_data[k].non_local_assemblers.push_back(
240 {la->getIPDataPtr(ip),
241 std::numeric_limits<double>::quiet_NaN(),
242 distances[ip]});
243 }
244 }
245 if (_ip_data[k].non_local_assemblers.size() == 0)
246 {
247 OGS_FATAL("no neighbours found!");
248 }
249
250 double a_k_sum_m = 0;
251 for (auto const& tuple : _ip_data[k].non_local_assemblers)
252 {
253 double const distance2_m = tuple.distance2;
254
255 auto const& w_m = tuple.ip_l_pointer->integration_weight;
256
257 a_k_sum_m += w_m * alpha_0(distance2_m);
258 }
259
260 //
261 // Calculate alpha_kl =
262 // alpha_0(|x_k - x_l|) / int_{m \in ip} alpha_0(|x_k - x_m|)
263 //
264 for (auto& tuple : _ip_data[k].non_local_assemblers)
265 {
266 double const distance2_l = tuple.distance2;
267 double const a_kl = alpha_0(distance2_l) / a_k_sum_m;
268
269 // Store the a_kl already multiplied with the integration
270 // weight of that l integration point.
271 auto const w_l = tuple.ip_l_pointer->integration_weight;
272 tuple.alpha_kl_times_w_l = a_kl * w_l;
273 }
274 }
275 }
276
278 int integration_point) const
279 {
280 auto const& N = _secondary_data.N[integration_point];
281
282 Eigen::Vector3d xyz = Eigen::Vector3d::Zero(); // Resulting coordinates
283 auto* nodes = _element.getNodes();
284 for (int i = 0; i < N.size(); ++i)
285 {
286 auto const& node_coordinates{nodes[i]->asEigenVector3d()};
287 xyz += node_coordinates * N[i];
288 }
289 return xyz;
290 }
291
296 Eigen::Vector3d const& coords,
297 std::vector<double>& distances) const override
298 {
299 unsigned const n_integration_points =
301
302 distances.resize(n_integration_points);
303
304 for (unsigned ip = 0; ip < n_integration_points; ip++)
305 {
306 auto const& xyz = _ip_data[ip].coordinates;
307 distances[ip] = (xyz - coords).squaredNorm();
308 }
309 }
310
311 void assemble(double const /*t*/, double const /*dt*/,
312 std::vector<double> const& /*local_x*/,
313 std::vector<double> const& /*local_x_prev*/,
314 std::vector<double>& /*local_M_data*/,
315 std::vector<double>& /*local_K_data*/,
316 std::vector<double>& /*local_b_data*/) override
317 {
318 OGS_FATAL(
319 "SmallDeformationNonlocalLocalAssembler: assembly without jacobian "
320 "is not "
321 "implemented.");
322 }
323
324 void preAssemble(double const t, double const dt,
325 std::vector<double> const& local_x) override
326 {
327 auto const n_integration_points =
329
330 MPL::VariableArray variables;
331 MPL::VariableArray variables_prev;
333 x_position.setElementID(_element.getID());
334
335 for (unsigned ip = 0; ip < n_integration_points; ip++)
336 {
337 auto const& N = _ip_data[ip].N;
338 auto const& dNdx = _ip_data[ip].dNdx;
339
340 auto const x_coord =
341 NumLib::interpolateXCoordinate<ShapeFunction,
343 auto const B = LinearBMatrix::computeBMatrix<
344 DisplacementDim, ShapeFunction::NPOINTS,
345 typename BMatricesType::BMatrixType>(dNdx, N, x_coord,
347 auto const& eps_prev = _ip_data[ip].eps_prev;
348 auto const& sigma_prev = _ip_data[ip].sigma_prev;
349
350 auto& eps = _ip_data[ip].eps;
351 auto& sigma = _ip_data[ip].sigma;
352 auto& C = _ip_data[ip].C;
353 auto& state = _ip_data[ip].material_state_variables;
354 double const& damage_prev = _ip_data[ip].damage_prev;
355
356 eps.noalias() =
357 B *
358 Eigen::Map<typename BMatricesType::NodalForceVectorType const>(
359 local_x.data(), ShapeFunction::NPOINTS * DisplacementDim);
360
361 // sigma is for plastic part only.
362 std::unique_ptr<
364 new_C;
365 std::unique_ptr<typename MaterialLib::Solids::MechanicsBase<
366 DisplacementDim>::MaterialStateVariables>
367 new_state;
368
369 // Compute sigma_eff from damage total stress sigma
370 using KelvinVectorType = typename BMatricesType::KelvinVectorType;
371 KelvinVectorType const sigma_eff_prev =
372 sigma_prev /
373 (1. - damage_prev); // damage_prev is in [0,1) range. See
374 // calculateDamage() function.
375
376 variables_prev.stress.emplace<
378 sigma_eff_prev);
379 variables_prev.mechanical_strain.emplace<
381 eps_prev);
382 variables_prev.temperature = _process_data.reference_temperature;
383 variables.mechanical_strain.emplace<
385 variables.temperature = _process_data.reference_temperature;
386
387 auto&& solution = _ip_data[ip].solid_material.integrateStress(
388 variables_prev, variables, t, x_position, dt, *state);
389
390 if (!solution)
391 {
392 OGS_FATAL("Computation of local constitutive relation failed.");
393 }
394
395 std::tie(sigma, state, C) = std::move(*solution);
396
398 {
399 auto const& ehlers_material =
401 DisplacementDim> const&>(_ip_data[ip].solid_material);
402 auto const damage_properties =
403 ehlers_material.evaluatedDamageProperties(t, x_position);
404 auto const material_properties =
405 ehlers_material.evaluatedMaterialProperties(t, x_position);
406
407 // Ehlers material state variables
408 auto& state_vars =
410 DisplacementDim>&>(
411 *_ip_data[ip].material_state_variables);
412
413 double const eps_p_eff_diff =
414 state_vars.eps_p.eff - state_vars.eps_p_prev.eff;
415
417 eps_p_eff_diff, sigma, _ip_data[ip].kappa_d_prev,
418 damage_properties.h_d, material_properties);
419
420 if (!_ip_data[ip].active_self)
421 {
422 _ip_data[ip].active_self |= _ip_data[ip].kappa_d > 0;
423 if (_ip_data[ip].active_self)
424 {
425 for (auto const& tuple :
426 _ip_data[ip].non_local_assemblers)
427 {
428 // Activate the integration point.
429 tuple.ip_l_pointer->activated = true;
430 }
431 }
432 }
433 }
434 }
435 }
436
437 void assembleWithJacobian(double const t, double const /*dt*/,
438 std::vector<double> const& local_x,
439 std::vector<double> const& /*local_x_prev*/,
440 std::vector<double>& local_b_data,
441 std::vector<double>& local_Jac_data) override
442 {
443 auto const local_matrix_size = local_x.size();
444
446 local_Jac_data, local_matrix_size, local_matrix_size);
447
449 local_b_data, local_matrix_size);
450
451 unsigned const n_integration_points =
453
455 x_position.setElementID(_element.getID());
456
457 // Non-local integration.
458 for (unsigned ip = 0; ip < n_integration_points; ip++)
459 {
460 auto const& w = _ip_data[ip].integration_weight;
461
462 auto const& N = _ip_data[ip].N;
463 auto const& dNdx = _ip_data[ip].dNdx;
464
465 auto const x_coord =
466 NumLib::interpolateXCoordinate<ShapeFunction,
468 auto const B = LinearBMatrix::computeBMatrix<
469 DisplacementDim, ShapeFunction::NPOINTS,
470 typename BMatricesType::BMatrixType>(dNdx, N, x_coord,
472
473 auto& sigma = _ip_data[ip].sigma;
474 auto& C = _ip_data[ip].C;
475 double& damage = _ip_data[ip].damage;
476
477 {
478 double nonlocal_kappa_d = 0;
479
480 if (_ip_data[ip].active_self || _ip_data[ip].activated)
481 {
482 for (auto const& tuple : _ip_data[ip].non_local_assemblers)
483 {
484 // Get local variable for the integration point l.
485 double const kappa_d_l = tuple.ip_l_pointer->kappa_d;
486 double const a_kl_times_w_l = tuple.alpha_kl_times_w_l;
487 nonlocal_kappa_d += a_kl_times_w_l * kappa_d_l;
488 }
489 }
490
491 auto const& ehlers_material =
493 DisplacementDim> const&>(_ip_data[ip].solid_material);
494
495 //
496 // Overnonlocal formulation
497 //
498 // See (Di Luzio & Bazant 2005, IJSS) for details.
499 // The implementation would go here and would be for a given
500 // gamma_nonlocal:
501 //
502 // Update nonlocal damage with local damage (scaled with 1 -
503 // \gamma_{nonlocal}) for the current integration point and the
504 // nonlocal integral part.
505 // nonlocal_kappa_d = (1. - gamma_nonlocal) * kappa_d +
506 // gamma_nonlocal * nonlocal_kappa_d;
507
508 nonlocal_kappa_d = std::max(0., nonlocal_kappa_d);
509
510 // Update damage based on nonlocal kappa_d
511 {
512 auto const damage_properties =
513 ehlers_material.evaluatedDamageProperties(t,
514 x_position);
515 damage = calculateDamage(nonlocal_kappa_d,
516 damage_properties.alpha_d,
517 damage_properties.beta_d);
518 damage = std::max(0., damage);
519 }
520 sigma = sigma * (1. - damage);
521 }
522
523 local_b.noalias() -= B.transpose() * sigma * w;
524 local_Jac.noalias() += B.transpose() * C * (1. - damage) * B * w;
525 }
526 }
527
528 void initializeConcrete() override
529 {
530 unsigned const n_integration_points =
532
533 for (unsigned ip = 0; ip < n_integration_points; ip++)
534 {
535 _ip_data[ip].pushBackState();
536 }
537 }
538
539 void postTimestepConcrete(Eigen::VectorXd const& /*local_x*/,
540 Eigen::VectorXd const& /*local_x_prev*/,
541 double const /*t*/, double const /*dt*/,
542 int const /*process_id*/) override
543 {
544 unsigned const n_integration_points =
546
547 for (unsigned ip = 0; ip < n_integration_points; ip++)
548 {
549 _ip_data[ip].pushBackState();
550 }
551 }
552
553 void computeCrackIntegral(std::size_t mesh_item_id,
554 NumLib::LocalToGlobalIndexMap const& dof_table,
555 GlobalVector const& x,
556 double& crack_volume) override
557 {
558 auto const indices = NumLib::getIndices(mesh_item_id, dof_table);
559 auto local_x = x.get(indices);
560
561 auto u = Eigen::Map<typename BMatricesType::NodalForceVectorType const>(
562 local_x.data(), ShapeFunction::NPOINTS * DisplacementDim);
563
564 int const n_integration_points =
566
567 for (int ip = 0; ip < n_integration_points; ip++)
568 {
569 auto const& dNdx = _ip_data[ip].dNdx;
570 auto const& d = _ip_data[ip].damage;
571 auto const& w = _ip_data[ip].integration_weight;
572
573 double const div_u =
574 (dNdx.template topRows<DisplacementDim>() *
575 u.reshaped(ShapeFunction::NPOINTS, DisplacementDim))
576 .trace();
577 crack_volume += div_u * d * w;
578 }
579 }
580
581 Eigen::Map<const Eigen::RowVectorXd> getShapeMatrix(
582 const unsigned integration_point) const override
583 {
584 auto const& N = _secondary_data.N[integration_point];
585
586 // assumes N is stored contiguously in memory
587 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
588 }
589
590 std::vector<double> const& getNodalValues(
591 std::vector<double>& nodal_values) const override
592 {
593 nodal_values.clear();
595 nodal_values, ShapeFunction::NPOINTS * DisplacementDim);
596
597 unsigned const n_integration_points =
599
600 for (unsigned ip = 0; ip < n_integration_points; ip++)
601 {
602 auto const& w = _ip_data[ip].integration_weight;
603
604 auto const& N = _ip_data[ip].N;
605 auto const& dNdx = _ip_data[ip].dNdx;
606
607 auto const x_coord =
608 NumLib::interpolateXCoordinate<ShapeFunction,
610 auto const B = LinearBMatrix::computeBMatrix<
611 DisplacementDim, ShapeFunction::NPOINTS,
612 typename BMatricesType::BMatrixType>(dNdx, N, x_coord,
614 auto& sigma = _ip_data[ip].sigma;
615
616 local_b.noalias() += B.transpose() * sigma * w;
617 }
618
619 return nodal_values;
620 }
621
622 std::vector<double> const& getIntPtFreeEnergyDensity(
623 const double /*t*/,
624 std::vector<GlobalVector*> const& /*x*/,
625 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
626 std::vector<double>& cache) const override
627 {
628 cache.clear();
629 cache.reserve(_ip_data.size());
630
631 transform(cbegin(_ip_data), cend(_ip_data), back_inserter(cache),
632 [](auto const& ip_data)
633 { return ip_data.free_energy_density; });
634
635 return cache;
636 }
637
638 std::vector<double> const& getIntPtEpsPV(
639 const double /*t*/,
640 std::vector<GlobalVector*> const& /*x*/,
641 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
642 std::vector<double>& cache) const override
643 {
644 cache.clear();
645 cache.reserve(_ip_data.size());
646
647 transform(cbegin(_ip_data), cend(_ip_data), back_inserter(cache),
648 [](auto const& ip_data) { return *ip_data.eps_p_V; });
649
650 return cache;
651 }
652
653 std::vector<double> const& getIntPtEpsPDXX(
654 const double /*t*/,
655 std::vector<GlobalVector*> const& /*x*/,
656 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
657 std::vector<double>& cache) const override
658 {
659 cache.clear();
660 cache.reserve(_ip_data.size());
661
662 transform(cbegin(_ip_data), cend(_ip_data), back_inserter(cache),
663 [](auto const& ip_data) { return *ip_data.eps_p_D_xx; });
664
665 return cache;
666 }
667
668 std::vector<double> const& getIntPtSigma(
669 const double /*t*/,
670 std::vector<GlobalVector*> const& /*x*/,
671 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
672 std::vector<double>& cache) const override
673 {
675 _ip_data, &IpData::sigma, cache);
676 }
677
678 std::vector<double> const& getIntPtEpsilon(
679 const double /*t*/,
680 std::vector<GlobalVector*> const& /*x*/,
681 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
682 std::vector<double>& cache) const override
683 {
685 _ip_data, &IpData::eps, cache);
686 }
687
688 std::size_t setSigma(double const* values)
689 {
691 values, _ip_data, &IpData::sigma);
692 }
693
694 // TODO (naumov) This method is same as getIntPtSigma but for arguments and
695 // the ordering of the cache_mat.
696 // There should be only one.
697 std::vector<double> getSigma() const override
698 {
701 }
702
703 void setKappaD(double value)
704 {
705 for (auto& ip_data : _ip_data)
706 {
707 ip_data.kappa_d = value;
708 }
709 }
710 std::vector<double> getKappaD() const override
711 {
712 unsigned const n_integration_points =
714
715 std::vector<double> result_values;
716 result_values.resize(n_integration_points);
717
718 for (unsigned ip = 0; ip < n_integration_points; ++ip)
719 {
720 result_values[ip] = _ip_data[ip].kappa_d;
721 }
722
723 return result_values;
724 }
725
726 std::vector<double> const& getIntPtDamage(
727 const double /*t*/,
728 std::vector<GlobalVector*> const& /*x*/,
729 std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/,
730 std::vector<double>& cache) const override
731 {
733 _ip_data, &IpData::damage, cache);
734 }
735
736 unsigned getNumberOfIntegrationPoints() const override
737 {
739 }
740
742 DisplacementDim>::MaterialStateVariables const&
743 getMaterialStateVariablesAt(int const integration_point) const override
744 {
745 return *_ip_data[integration_point].material_state_variables;
746 }
747
748private:
749 std::vector<double> const& getIntPtSigma(std::vector<double>& cache,
750 std::size_t const component) const
751 {
752 cache.clear();
753 cache.reserve(_ip_data.size());
754
755 for (auto const& ip_data : _ip_data)
756 {
757 if (component < 3)
758 { // xx, yy, zz components
759 cache.push_back(ip_data.sigma[component]);
760 }
761 else
762 { // mixed xy, yz, xz components
763 cache.push_back(ip_data.sigma[component] / std::sqrt(2));
764 }
765 }
766
767 return cache;
768 }
769
770 std::vector<double> const& getIntPtEpsilon(
771 std::vector<double>& cache, std::size_t const component) const
772 {
773 cache.clear();
774 cache.reserve(_ip_data.size());
775
776 for (auto const& ip_data : _ip_data)
777 {
778 if (component < 3) // xx, yy, zz components
779 cache.push_back(ip_data.eps[component]);
780 else // mixed xy, yz, xz components
781 cache.push_back(ip_data.eps[component] / std::sqrt(2));
782 }
783
784 return cache;
785 }
786
788 {
789 return &_ip_data[ip];
790 }
791
792private:
794
795 std::vector<IpData, Eigen::aligned_allocator<IpData>> _ip_data;
796
801
802 static const int displacement_size =
803 ShapeFunction::NPOINTS * DisplacementDim;
804};
805
806} // namespace SmallDeformationNonlocal
807} // namespace ProcessLib
#define OGS_FATAL(...)
Definition Error.h:26
DamageProperties evaluatedDamageProperties(double const t, ParameterLib::SpatialPosition const &x) const
Definition Ehlers.h:343
Global vector based on Eigen vector.
Definition EigenVector.h:25
double get(IndexType rowId) const
get entry
Definition EigenVector.h:58
constexpr double getWeight() const
virtual Node *const * getNodes() const =0
Get array of element nodes.
std::size_t getID() const
Returns the ID of the element.
Definition Element.h:89
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
void setElementID(std::size_t element_id)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
MatrixType< _number_of_dof, _number_of_dof > StiffnessMatrixType
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.
VectorType< _kelvin_vector_size > KelvinVectorType
void postTimestepConcrete(Eigen::VectorXd const &, Eigen::VectorXd const &, double const, double const, int const) override
std::size_t setIPDataInitialConditions(std::string_view const name, double const *values, int const integration_order) override
void computeCrackIntegral(std::size_t mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double &crack_volume) override
std::vector< double > const & getIntPtEpsilon(std::vector< double > &cache, std::size_t const component) const
void assemble(double const, double const, std::vector< double > const &, std::vector< double > const &, std::vector< double > &, std::vector< double > &, std::vector< double > &) override
IntegrationPointDataNonlocalInterface * getIPDataPtr(int const ip) override
SmallDeformationNonlocalLocalAssembler(MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, SmallDeformationNonlocalProcessData< DisplacementDim > &process_data)
void setIPDataInitialConditionsFromCellData(std::string const &name, std::vector< double > const &value) override
void nonlocal(std::size_t const, std::vector< std::unique_ptr< SmallDeformationNonlocalLocalAssemblerInterface< DisplacementDim > > > const &local_assemblers) override
std::vector< double > const & getIntPtSigma(std::vector< double > &cache, std::size_t const component) const
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
std::vector< double > const & getIntPtEpsPDXX(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
MaterialLib::Solids::MechanicsBase< DisplacementDim >::MaterialStateVariables const & getMaterialStateVariablesAt(int const integration_point) const override
std::vector< double > const & getIntPtEpsPV(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
void preAssemble(double const t, double const dt, std::vector< double > const &local_x) override
std::vector< double > const & getIntPtFreeEnergyDensity(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
SmallDeformationNonlocalLocalAssembler(SmallDeformationNonlocalLocalAssembler const &)=delete
std::vector< double > const & getIntPtSigma(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtEpsilon(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
void getIntegrationPointCoordinates(Eigen::Vector3d const &coords, std::vector< double > &distances) const override
std::vector< double > const & getNodalValues(std::vector< double > &nodal_values) const override
void assembleWithJacobian(double const t, double const, std::vector< double > const &local_x, std::vector< double > const &, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
SmallDeformationNonlocalLocalAssembler(SmallDeformationNonlocalLocalAssembler &&)=delete
std::vector< double > const & getIntPtDamage(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
auto & selectSolidConstitutiveRelation(SolidMaterialsMap const &constitutive_relations, MeshLib::PropertyVector< int > const *const material_ids, std::size_t const element_id)
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< std::size_t > findElementsWithinRadius(Element const &start_element, double const radius_squared)
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
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)
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.
double calculateDamage(double const kappa_d, double const alpha_d, double const beta_d)
Definition Damage.h:25
double calculateDamageKappaD(double const eps_p_eff_diff, KelvinVectorType const &sigma, double const kappa_d_prev, double const h_d, MaterialLib::Solids::Ehlers::MaterialProperties const &mp)
Definition Damage.h:39
std::vector< double > const & getIntegrationPointScalarData(IntegrationPointDataVector const &ip_data_vector, MemberType IpData::*const member, std::vector< double > &cache)
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)
std::size_t setIntegrationPointScalarData(double const *values, IntegrationPointDataVector &ip_data_vector, MemberType IpData::*const member)
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
VectorType< GlobalDim > GlobalDimVectorType
VectorType< ShapeFunction::NPOINTS > NodalVectorType
PlasticStrain< KelvinVector > eps_p
plastic part of the state.
Definition Ehlers.h:244
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N