Loading [MathJax]/extensions/tex2jax.js
OGS
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, int DisplacementDim>
class ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >

Definition at line 37 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

#include <SmallDeformationLocalAssemblerMatrixNearFracture.h>

Inheritance diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >:
[legend]

Public Types

using ShapeMatricesType
using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
using BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>
using BBarMatrixType = typename BMatricesType::BBarMatrixType
using BMatrixType = typename BMatricesType::BMatrixType
using StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
using NodalDisplacementVectorType

Public Member Functions

 SmallDeformationLocalAssemblerMatrixNearFracture (SmallDeformationLocalAssemblerMatrixNearFracture const &)=delete
 SmallDeformationLocalAssemblerMatrixNearFracture (SmallDeformationLocalAssemblerMatrixNearFracture &&)=delete
 SmallDeformationLocalAssemblerMatrixNearFracture (MeshLib::Element const &e, std::size_t const n_variables, std::size_t const local_matrix_size, std::vector< unsigned > const &dofIndex_to_localIndex, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, SmallDeformationProcessData< DisplacementDim > &process_data)
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 assembleWithJacobian (double const t, double const dt, Eigen::VectorXd const &local_u, Eigen::VectorXd &local_b, Eigen::MatrixXd &local_J) override
void preTimestepConcrete (std::vector< double > const &, double const, double const) override
void computeSecondaryVariableConcreteWithVector (double const t, Eigen::VectorXd const &local_x) override
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 & getIntPtSigma (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 & getIntPtEpsilon (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 & getIntPtFractureStress (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
std::vector< double > const & getIntPtFractureAperture (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
Public Member Functions inherited from ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface
 SmallDeformationLocalAssemblerInterface ()
 SmallDeformationLocalAssemblerInterface (std::size_t n_local_size, std::vector< unsigned > dofIndex_to_localIndex)
void assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x_, std::vector< double > const &, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data) override
void computeSecondaryVariableConcrete (double const t, double const, Eigen::VectorXd const &local_x, Eigen::VectorXd const &) override
Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual ~LocalAssemblerInterface ()=default
virtual void setInitialConditions (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, int const process_id)
virtual void initialize (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
virtual void preAssemble (double const, double const, std::vector< double > const &)
virtual void assembleForStaggeredScheme (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_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
virtual 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)
virtual void computeSecondaryVariable (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_prev, int const process_id)
virtual void preTimestep (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void postTimestep (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
void postNonLinearSolver (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id)
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< double > const &) const
virtual Eigen::Vector3d getFlux (MathLib::Point3d const &, double const, std::vector< std::vector< double > > const &) const
 Fits to staggered scheme.
virtual std::optional< VectorSegmentgetVectorDeformationSegment () const
Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default

Private Types

using IntegrationPointDataType

Private Member Functions

std::optional< BBarMatrixTypegetDilatationalBBarMatrix () const

Private Attributes

SmallDeformationProcessData< DisplacementDim > & _process_data
std::vector< FractureProperty * > _fracture_props
std::vector< JunctionProperty * > _junction_props
std::unordered_map< int, int > _fracID_to_local
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
NumLib::GenericIntegrationMethod const & _integration_method
MeshLib::Element const & _element
bool const _is_axially_symmetric
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data

Member Typedef Documentation

◆ BBarMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::BBarMatrixType = typename BMatricesType::BBarMatrixType

◆ BMatricesType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>

◆ BMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::BMatrixType = typename BMatricesType::BMatrixType

◆ IntegrationPointDataType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::IntegrationPointDataType
private

◆ NodalDisplacementVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::NodalDisplacementVectorType
Initial value:
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.

Definition at line 52 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

◆ NodalForceVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename BMatricesType::NodalForceVectorType

◆ NodalMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::NodalMatrixType = typename ShapeMatricesType::NodalMatrixType

◆ NodalVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType

◆ ShapeMatrices

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

◆ ShapeMatricesType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::ShapeMatricesType
Initial value:
EigenFixedShapeMatrixPolicy< ShapeFunction, GlobalDim > ShapeMatrixPolicyType

Definition at line 41 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

◆ StiffnessMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType

Constructor & Destructor Documentation

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [1/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrixNearFracture ( SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim > const & )
delete

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [2/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrixNearFracture ( SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim > && )
delete

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [3/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerMatrixNearFracture ( MeshLib::Element const & e,
std::size_t const n_variables,
std::size_t const local_matrix_size,
std::vector< unsigned > const & dofIndex_to_localIndex,
NumLib::GenericIntegrationMethod const & integration_method,
bool const is_axially_symmetric,
SmallDeformationProcessData< DisplacementDim > & process_data )

Definition at line 47 of file SmallDeformationLocalAssemblerMatrixNearFracture-impl.h.

61 _element(e),
63{
70
71 unsigned const n_integration_points =
72 _integration_method.getNumberOfPoints();
73
76
78 _process_data.solid_materials, _process_data.material_ids, e.getID());
79
80 for (unsigned ip = 0; ip < n_integration_points; ip++)
81 {
82 _ip_data.emplace_back(solid_material);
83 auto& ip_data = _ip_data[ip];
84 auto const& sm = shape_matrices[ip];
85 ip_data.N_u = sm.N;
86 ip_data.dNdx_u = sm.dNdx;
87 ip_data.integration_weight =
88 _integration_method.getWeightedPoint(ip).getWeight() *
89 sm.integralMeasure * sm.detJ;
90
91 // Initialize current time step values
92 static const int kelvin_vector_size =
94 ip_data._sigma.setZero(kelvin_vector_size);
95 ip_data._eps.setZero(kelvin_vector_size);
96
97 // Previous time step values are not initialized and are set later.
98 ip_data._sigma_prev.resize(kelvin_vector_size);
99 ip_data._eps_prev.resize(kelvin_vector_size);
100
102
103 _secondary_data.N[ip] = sm.N;
104 }
105
106 for (auto fid : process_data.vec_ele_connected_fractureIDs[e.getID()])
107 {
108 _fracID_to_local.insert({fid, _fracture_props.size()});
109 _fracture_props.push_back(&_process_data.fracture_properties[fid]);
110 }
111
112 _junction_props = process_data.vec_ele_connected_junctionIDs[e.getID()] |
114 [&](auto const jid)
115 { return &_process_data.junction_properties[jid]; }) |
117}
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
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.
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface::SmallDeformationLocalAssemblerInterface(), _element, _fracID_to_local, _fracture_props, _integration_method, _ip_data, _is_axially_symmetric, _junction_props, _process_data, _secondary_data, MeshLib::Element::getID(), NumLib::initShapeMatrices(), MathLib::KelvinVector::kelvin_vector_dimensions(), MaterialLib::Solids::selectSolidConstitutiveRelation(), ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::vec_ele_connected_fractureIDs, and ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::vec_ele_connected_junctionIDs.

Member Function Documentation

◆ assemble()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::assemble ( double const ,
double const ,
std::vector< double > const & ,
std::vector< double > const & ,
std::vector< double > & ,
std::vector< double > & ,
std::vector< double > &  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 69 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

75 {
77 "SmallDeformationLocalAssembler: assembly without jacobian is not "
78 "implemented.");
79 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobian()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::assembleWithJacobian ( double const t,
double const dt,
Eigen::VectorXd const & local_u,
Eigen::VectorXd & local_b,
Eigen::MatrixXd & local_J )
overridevirtual

Reimplemented from ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 122 of file SmallDeformationLocalAssemblerMatrixNearFracture-impl.h.

126{
127 assert(_element.getDimension() == DisplacementDim);
128
130 auto const n_fractures = _fracture_props.size();
131 auto const n_junctions = _junction_props.size();
132 auto const n_enrich_var = n_fractures + n_junctions;
133
134 using BlockVectorType =
136 using BlockMatrixType =
138
139 //--------------------------------------------------------------------------------------
140 // prepare sub vectors, matrices for regular displacement (u) and
141 // displacement jumps (g)
142 //
143 // example with two fractures with one intersection:
144 // |b(u)|
145 // b = |b(g1)|
146 // |b(g2)|
147 // |b(j1)|
148 //
149 // |J(u,u) J(u,g1) J(u,g2) J(u,j1) |
150 // J = |J(g1,u) J(g1,g1) J(g1,g2) J(g1,j1)|
151 // |J(g2,u) J(g2,g1) J(g2,g2) J(g2,j1)|
152 // |J(j1,u) J(j1,g1) J(j1,g2) J(j1,j1)|
153 //--------------------------------------------------------------------------------------
154 auto local_b_u = local_b.segment<N_DOF_PER_VAR>(0);
156 for (unsigned i = 0; i < n_enrich_var; i++)
157 {
158 vec_local_b_g.push_back(
159 local_b.segment<N_DOF_PER_VAR>(N_DOF_PER_VAR * (i + 1)));
160 }
161
162 auto local_J_uu = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(0, 0);
166 for (unsigned i = 0; i < n_enrich_var; i++)
167 {
169 0, N_DOF_PER_VAR * (i + 1));
170 vec_local_J_ug.push_back(sub_ug);
171
173 N_DOF_PER_VAR * (i + 1), 0);
174 vec_local_J_gu.push_back(sub_gu);
175
176 for (unsigned j = 0; j < n_enrich_var; j++)
177 {
179 N_DOF_PER_VAR * (i + 1), N_DOF_PER_VAR * (j + 1));
180 vec_local_J_gg[i].push_back(sub_gg);
181 }
182 }
183
184 auto const nodal_u = local_u.segment<N_DOF_PER_VAR>(0);
186 for (unsigned i = 0; i < n_enrich_var; i++)
187 {
189 N_DOF_PER_VAR * (i + 1));
190 vec_nodal_g.push_back(sub);
191 }
192
193 //------------------------------------------------
194 // integration
195 //------------------------------------------------
196 unsigned const n_integration_points =
197 _integration_method.getNumberOfPoints();
198
202 x_position.setElementID(_element.getID());
203
205
206 for (unsigned ip = 0; ip < n_integration_points; ip++)
207 {
208 auto& ip_data = _ip_data[ip];
209 auto const& w = _ip_data[ip].integration_weight;
210
211 auto const& N = ip_data.N_u;
212 auto const& dNdx = ip_data.dNdx_u;
213
214 // levelset functions
220
221 // u = u^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x) *
222 // [u]_i)
224 for (unsigned i = 0; i < n_enrich_var; i++)
225 {
227 }
228
229 auto const x_coord =
231 _element, N);
236
237 // strain, stress
238 auto const& eps_prev = ip_data._eps_prev;
239 auto const& sigma_prev = ip_data._sigma_prev;
240
241 auto& eps = ip_data._eps;
242 auto& sigma = ip_data._sigma;
243 auto& state = ip_data._material_state_variables;
244
245 eps.noalias() = B * nodal_total_u;
246
247 variables.mechanical_strain
249 eps);
250
251 variables_prev.stress
253 sigma_prev);
254 variables_prev.mechanical_strain
256 eps_prev);
257 variables_prev.temperature = _process_data.reference_temperature;
258
259 auto&& solution = _ip_data[ip]._solid_material.integrateStress(
261
262 if (!solution)
263 {
264 OGS_FATAL("Computation of local constitutive relation failed.");
265 }
266
269
270 // r_u = B^T * Sigma = B^T * C * B * (u+phi*[u])
271 // r_[u] = (phi*B)^T * Sigma = (phi*B)^T * C * B * (u+phi*[u])
272 local_b_u.noalias() -= B.transpose() * sigma * w;
273 for (unsigned i = 0; i < n_enrich_var; i++)
274 {
275 vec_local_b_g[i].noalias() -=
276 levelsets[i] * B.transpose() * sigma * w;
277 }
278
279 // J_uu += B^T * C * B
280 local_J_uu.noalias() += B.transpose() * C * B * w;
281
282 for (unsigned i = 0; i < n_enrich_var; i++)
283 {
284 // J_u[u] += B^T * C * (levelset * B)
285 vec_local_J_ug[i].noalias() +=
286 B.transpose() * C * (levelsets[i] * B) * w;
287
288 // J_[u]u += (levelset * B)^T * C * B
289 vec_local_J_gu[i].noalias() +=
290 (levelsets[i] * B.transpose()) * C * B * w;
291
292 for (unsigned j = 0; j < n_enrich_var; j++)
293 {
294 // J_[u][u] += (levelset * B)^T * C * (levelset * B)
295 vec_local_J_gg[i][j].noalias() +=
296 (levelsets[i] * B.transpose()) * C * (levelsets[j] * B) * w;
297 }
298 }
299 }
300}
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
Eigen::Vector3d computePhysicalCoordinates(MeshLib::Element const &e, Eigen::MatrixBase< Derived > const &shape)
Definition Utils.h:24
std::vector< double > uGlobalEnrichments(std::vector< FractureProperty * > const &frac_props, std::vector< JunctionProperty * > const &junction_props, std::unordered_map< int, int > const &fracID_to_local, Eigen::Vector3d const &x)

References _element, _fracID_to_local, _fracture_props, _integration_method, _ip_data, _is_axially_symmetric, _junction_props, _process_data, assembleWithJacobian(), ProcessLib::LinearBMatrix::computeBMatrixPossiblyWithBbar(), ProcessLib::LIE::computePhysicalCoordinates(), getDilatationalBBarMatrix(), NumLib::interpolateXCoordinate(), MaterialPropertyLib::VariableArray::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), MaterialPropertyLib::VariableArray::stress, MaterialPropertyLib::VariableArray::temperature, and ProcessLib::LIE::uGlobalEnrichments().

Referenced by assembleWithJacobian().

◆ computeSecondaryVariableConcreteWithVector()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::computeSecondaryVariableConcreteWithVector ( double const t,
Eigen::VectorXd const & local_x )
overridevirtual

Implements ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface.

Definition at line 307 of file SmallDeformationLocalAssemblerMatrixNearFracture-impl.h.

310{
311 // Compute average value per element
314 auto const e_id = _element.getID();
315
316 unsigned const n_integration_points =
317 _integration_method.getNumberOfPoints();
318 for (unsigned ip = 0; ip < n_integration_points; ip++)
319 {
320 sigma_avg += _ip_data[ip]._sigma;
321 }
323
325 &(*_process_data.element_stresses)[e_id * KV::RowsAtCompileTime]) =
327}
Eigen::Matrix< double, 4, 1 > kelvinVectorToSymmetricTensor(Eigen::Matrix< double, 4, 1, Eigen::ColMajor, 4, 1 > const &v)

References _element, _integration_method, _ip_data, _process_data, and MathLib::KelvinVector::kelvinVectorToSymmetricTensor().

◆ getDilatationalBBarMatrix()

template<typename ShapeFunction, int DisplacementDim>
std::optional< BBarMatrixType > ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::getDilatationalBBarMatrix ( ) const
inlineprivate

◆ getIntPtEpsilon()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::getIntPtEpsilon ( const double t,
std::vector< GlobalVector * > const & x,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & dof_table,
std::vector< double > & cache ) const
overridevirtual

◆ getIntPtFractureAperture()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::getIntPtFractureAperture ( const double ,
std::vector< GlobalVector * > const & ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & ,
std::vector< double > & cache ) const
inlineoverridevirtual

◆ getIntPtFractureStress()

template<typename ShapeFunction, int DisplacementDim>
std::vector< double > const & ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::getIntPtFractureStress ( const double ,
std::vector< GlobalVector * > const & ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const & ,
std::vector< double > & cache ) const
inlineoverridevirtual

◆ getIntPtSigma()

◆ getShapeMatrix()

template<typename ShapeFunction, int DisplacementDim>
Eigen::Map< const Eigen::RowVectorXd > ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::getShapeMatrix ( const unsigned integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 102 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

104 {
105 auto const& N = _secondary_data.N[integration_point];
106
107 // assumes N is stored contiguously in memory
108 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
109 }

References _secondary_data.

◆ preTimestepConcrete()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::preTimestepConcrete ( std::vector< double > const & ,
double const ,
double const  )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 86 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

89 {
90 unsigned const n_integration_points =
91 _integration_method.getNumberOfPoints();
92
93 for (unsigned ip = 0; ip < n_integration_points; ip++)
94 {
95 _ip_data[ip].pushBackState();
96 }
97 }

References _integration_method, and _ip_data.

Member Data Documentation

◆ _element

◆ _fracID_to_local

template<typename ShapeFunction, int DisplacementDim>
std::unordered_map<int, int> ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::_fracID_to_local
private

◆ _fracture_props

template<typename ShapeFunction, int DisplacementDim>
std::vector<FractureProperty*> ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::_fracture_props
private

◆ _integration_method

◆ _ip_data

◆ _is_axially_symmetric

template<typename ShapeFunction, int DisplacementDim>
bool const ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::_is_axially_symmetric
private

◆ _junction_props

template<typename ShapeFunction, int DisplacementDim>
std::vector<JunctionProperty*> ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::_junction_props
private

◆ _process_data

◆ _secondary_data

template<typename ShapeFunction, int DisplacementDim>
SecondaryData<typename ShapeMatrices::ShapeType> ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, DisplacementDim >::_secondary_data
private

The documentation for this class was generated from the following files: