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

Detailed Description

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

Definition at line 32 of file SmallDeformationLocalAssemblerFracture.h.

#include <SmallDeformationLocalAssemblerFracture.h>

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

Public Types

using ShapeMatricesType
using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
using HMatricesType = HMatrixPolicyType<ShapeFunction, DisplacementDim>
using HMatrixType = typename HMatricesType::HMatrixType
using StiffnessMatrixType = typename HMatricesType::StiffnessMatrixType
using NodalForceVectorType = typename HMatricesType::NodalForceVectorType
using NodalDisplacementVectorType
using ForceVectorType = typename HMatricesType::ForceVectorType
using GlobalDimVectorType = Eigen::Matrix<double, DisplacementDim, 1>

Public Member Functions

 SmallDeformationLocalAssemblerFracture (SmallDeformationLocalAssemblerFracture const &)=delete
 SmallDeformationLocalAssemblerFracture (SmallDeformationLocalAssemblerFracture &&)=delete
 SmallDeformationLocalAssemblerFracture (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 (const double t, Eigen::VectorXd const &local_u) 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, 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
std::vector< double > const & getIntPtFractureStress (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 & getIntPtFractureAperture (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, 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 Attributes

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

Member Typedef Documentation

◆ ForceVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::ForceVectorType = typename HMatricesType::ForceVectorType

Definition at line 49 of file SmallDeformationLocalAssemblerFracture.h.

◆ GlobalDimVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::GlobalDimVectorType = Eigen::Matrix<double, DisplacementDim, 1>

Definition at line 50 of file SmallDeformationLocalAssemblerFracture.h.

◆ HMatricesType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::HMatricesType = HMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 41 of file SmallDeformationLocalAssemblerFracture.h.

◆ HMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::HMatrixType = typename HMatricesType::HMatrixType

Definition at line 43 of file SmallDeformationLocalAssemblerFracture.h.

◆ IntegrationPointDataType

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

◆ NodalDisplacementVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::NodalDisplacementVectorType
Initial value:
VectorType< _number_of_dof > NodalForceVectorType

Definition at line 46 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalForceVectorType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename HMatricesType::NodalForceVectorType

Definition at line 45 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalMatrixType

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

Definition at line 38 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalVectorType

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

Definition at line 39 of file SmallDeformationLocalAssemblerFracture.h.

◆ ShapeMatrices

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

Definition at line 40 of file SmallDeformationLocalAssemblerFracture.h.

◆ ShapeMatricesType

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

Definition at line 36 of file SmallDeformationLocalAssemblerFracture.h.

◆ StiffnessMatrixType

template<typename ShapeFunction, int DisplacementDim>
using ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::StiffnessMatrixType = typename HMatricesType::StiffnessMatrixType

Definition at line 44 of file SmallDeformationLocalAssemblerFracture.h.

Constructor & Destructor Documentation

◆ SmallDeformationLocalAssemblerFracture() [1/3]

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

◆ SmallDeformationLocalAssemblerFracture() [2/3]

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

◆ SmallDeformationLocalAssemblerFracture() [3/3]

template<typename ShapeFunction, int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::SmallDeformationLocalAssemblerFracture ( 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 31 of file SmallDeformationLocalAssemblerFracture-impl.h.

49 _element(e)
50{
51 assert(_element.getDimension() == DisplacementDim - 1);
52
53 unsigned const n_integration_points =
54 _integration_method.getNumberOfPoints();
55
58
59 auto mat_id = (*_process_data.mesh_prop_materialIDs)[e.getID()];
60 auto frac_id = _process_data.map_materialID_to_fractureID[mat_id];
61 _fracture_property = &_process_data.fracture_properties[frac_id];
62 for (auto fid : process_data.vec_ele_connected_fractureIDs[e.getID()])
63 {
64 _fracID_to_local.insert({fid, _fracture_props.size()});
65 _fracture_props.push_back(&_process_data.fracture_properties[fid]);
66 }
67
68 _junction_props = process_data.vec_ele_connected_junctionIDs[e.getID()] |
70 [&](auto const jid)
71 { return &_process_data.junction_properties[jid]; }) |
73
74 for (unsigned ip = 0; ip < n_integration_points; ip++)
75 {
76 _ip_data.emplace_back(*_process_data.fracture_model);
77 auto const& sm = _shape_matrices[ip];
78 auto& ip_data = _ip_data[ip];
79
81 std::nullopt, _element.getID(),
84 _element, sm.N))};
85
86 ip_data.integration_weight =
87 _integration_method.getWeightedPoint(ip).getWeight() *
88 sm.integralMeasure * sm.detJ;
89 ip_data.h_matrices.setZero(DisplacementDim,
91
94 HMatrixType>(sm.N, ip_data.h_matrices);
95
96 // Initialize current time step values
97 ip_data.w.setZero(DisplacementDim);
98 ip_data.sigma.setZero(DisplacementDim);
99
100 // Previous time step values are not initialized and are set later.
101 ip_data.sigma_prev.resize(DisplacementDim);
102 ip_data.w_prev.resize(DisplacementDim);
103
105
106 ip_data.aperture0 = _fracture_property->aperture0(0, x_position)[0];
107 ip_data.aperture_prev = ip_data.aperture0;
108
109 _secondary_data.N[ip] = sm.N;
110 }
111}
std::vector< IntegrationPointDataType, Eigen::aligned_allocator< IntegrationPointDataType > > _ip_data
std::vector< ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > _shape_matrices
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerInterface::SmallDeformationLocalAssemblerInterface(), _element, _fracID_to_local, _fracture_property, _fracture_props, _integration_method, _ip_data, _junction_props, _process_data, _secondary_data, _shape_matrices, ProcessLib::computeHMatrix(), MeshLib::Element::getID(), NumLib::interpolateCoordinates(), 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::SmallDeformationLocalAssemblerFracture< 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 66 of file SmallDeformationLocalAssemblerFracture.h.

72 {
74 "SmallDeformationLocalAssembler: assembly without jacobian is not "
75 "implemented.");
76 }
#define OGS_FATAL(...)
Definition Error.h:26

References OGS_FATAL.

◆ assembleWithJacobian()

template<typename ShapeFunction, int DisplacementDim>
void ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< 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 114 of file SmallDeformationLocalAssemblerFracture-impl.h.

118{
120 auto const n_fractures = _fracture_props.size();
121 auto const n_junctions = _junction_props.size();
122 auto const n_enrich_var = n_fractures + n_junctions;
123
124 //--------------------------------------------------------------------------------------
125 // prepare sub vectors, matrices for regular displacement (u) and
126 // displacement jumps (g)
127 //
128 // example with two fractures with one intersection:
129 // b = |b(g1)|
130 // |b(g2)|
131 //
132 // J = |J(g1,g1) J(g1,g2)|
133 // |J(g2,g1) J(g2,g2)|
134 //--------------------------------------------------------------------------------------
135
136 using BlockVectorType =
138 using BlockMatrixType =
140
142 for (unsigned i = 0; i < n_enrich_var; i++)
143 {
144 vec_local_b_g.push_back(
146 }
148 for (unsigned i = 0; i < n_enrich_var; i++)
149 {
150 for (unsigned j = 0; j < n_enrich_var; j++)
151 {
154 vec_local_J_gg[i].push_back(sub_gg);
155 }
156 }
157
159 for (unsigned i = 0; i < n_enrich_var; i++)
160 {
162 N_DOF_PER_VAR * i);
163 vec_nodal_g.push_back(sub);
164 }
165
166 //------------------------------------------------
167 // integration
168 //------------------------------------------------
169 // the index of a normal (normal to a fracture plane) component
170 // in a displacement vector
171 int const index_normal = DisplacementDim - 1;
172 auto const& R = _fracture_property->R;
173
174 unsigned const n_integration_points =
175 _integration_method.getNumberOfPoints();
176
178 x_position.setElementID(_element.getID());
179
180 for (unsigned ip = 0; ip < n_integration_points; ip++)
181 {
182 auto& ip_data = _ip_data[ip];
183 auto const& integration_weight = ip_data.integration_weight;
184 auto const& H = ip_data.h_matrices;
185 auto& mat = ip_data.fracture_material;
186 auto& sigma = ip_data.sigma;
187 auto const& sigma_prev = ip_data.sigma_prev;
188 auto& w = ip_data.w;
189 auto const& w_prev = ip_data.w_prev;
190 auto& C = ip_data.C;
191 auto& state = *ip_data.material_state_variables;
192 auto const& N = _secondary_data.N[ip];
193
198
199 // du = du^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x)
200 // * [u]_i)
202 nodal_gap.setZero();
203 for (unsigned i = 0; i < n_enrich_var; i++)
204 {
206 }
207
208 // displacement jumps
209 w.noalias() = R * H * nodal_gap;
210
211 // total aperture
212 ip_data.aperture = ip_data.aperture0 + w[index_normal];
213
214 // local C, local stress
215 mat.computeConstitutiveRelation(
216 t, x_position, ip_data.aperture0,
218 // Replace with
219 // initial
220 // stress values
222
223 // r_[u] += H^T*Stress
224 for (unsigned i = 0; i < n_enrich_var; i++)
225 {
226 vec_local_b_g[i].noalias() -= levelsets[i] * H.transpose() *
227 R.transpose() * sigma *
229 }
230
231 // J_[u][u] += H^T*C*H
232 for (unsigned i = 0; i < n_enrich_var; i++)
233 {
234 for (unsigned j = 0; j < n_enrich_var; j++)
235 {
236 // J_[u][u] += (levelset * B)^T * C * (levelset * B)
237 vec_local_J_gg[i][j].noalias() +=
238 (levelsets[i] * H.transpose() * R.transpose()) * C *
240 }
241 }
242 }
243}
std::vector< double > duGlobalEnrichments(std::size_t this_frac_id, 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)
Eigen::Vector3d computePhysicalCoordinates(MeshLib::Element const &e, Eigen::MatrixBase< Derived > const &shape)
Definition Utils.h:24

References _element, _fracID_to_local, _fracture_property, _fracture_props, _integration_method, _ip_data, _junction_props, _secondary_data, ProcessLib::LIE::computePhysicalCoordinates(), ProcessLib::LIE::duGlobalEnrichments(), and ParameterLib::SpatialPosition::setElementID().

◆ computeSecondaryVariableConcreteWithVector()

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

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

Definition at line 246 of file SmallDeformationLocalAssemblerFracture-impl.h.

249{
251 auto const n_fractures = _fracture_props.size();
252 auto const n_junctions = _junction_props.size();
253 auto const n_enrich_var = n_fractures + n_junctions;
254
255 auto const& R = _fracture_property->R;
256
257 // the index of a normal (normal to a fracture plane) component
258 // in a displacement vector
259 int const index_normal = DisplacementDim - 1;
260
261 unsigned const n_integration_points =
262 _integration_method.getNumberOfPoints();
263
265 auto const e_id = _element.getID();
266 x_position.setElementID(e_id);
267
269 for (unsigned i = 0; i < n_enrich_var; i++)
270 {
272 N_DOF_PER_VAR * i);
273 vec_nodal_g.push_back(sub);
274 }
275
276 for (unsigned ip = 0; ip < n_integration_points; ip++)
277 {
278 auto& ip_data = _ip_data[ip];
279 auto const& H = ip_data.h_matrices;
280 auto& mat = ip_data.fracture_material;
281 auto& sigma = ip_data.sigma;
282 auto const& sigma_prev = ip_data.sigma_prev;
283 auto& w = ip_data.w;
284 auto const& w_prev = ip_data.w_prev;
285 auto& C = ip_data.C;
286 auto& state = *ip_data.material_state_variables;
287 auto& b_m = ip_data.aperture;
288 auto const& N = _secondary_data.N[ip];
289
294
295 // du = du^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x)
296 // * [u]_i)
298 nodal_gap.setZero();
299 for (unsigned i = 0; i < n_enrich_var; i++)
300 {
302 }
303
304 // displacement jumps in local coordinates
305 w.noalias() = R * H * nodal_gap;
306
307 // aperture
308 b_m = ip_data.aperture0 + w[index_normal];
309 if (b_m < 0.0)
310 {
311 OGS_FATAL(
312 "Element {:d}, gp {:d}: Fracture aperture is {:g}, but it must "
313 "be non-negative.",
314 _element.getID(), ip, b_m);
315 }
316
317 // local C, local stress
318 mat.computeConstitutiveRelation(
319 t, x_position, ip_data.aperture0,
321 // Replace with
322 // initial
323 // stress values
325 }
326
327 double ele_b = 0;
330
331 for (unsigned ip = 0; ip < n_integration_points; ip++)
332 {
333 auto const& ip_data = _ip_data[ip];
334 ele_b += ip_data.aperture;
335 ele_w += ip_data.w;
336 ele_sigma += ip_data.sigma;
337 }
341 (*_process_data.mesh_prop_b)[_element.getID()] = ele_b;
342
344 &(*_process_data.element_fracture_stresses)[e_id * DisplacementDim]) =
345 ele_sigma;
346
348 &(*_process_data.element_local_jumps)[e_id * DisplacementDim]) = ele_w;
349}

References _element, _fracID_to_local, _fracture_property, _fracture_props, _integration_method, _ip_data, _junction_props, _process_data, _secondary_data, ProcessLib::LIE::computePhysicalCoordinates(), ProcessLib::LIE::duGlobalEnrichments(), OGS_FATAL, and ParameterLib::SpatialPosition::setElementID().

◆ getIntPtEpsilon()

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

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

Definition at line 118 of file SmallDeformationLocalAssemblerFracture.h.

123 {
124 cache.resize(0);
125 return cache;
126 }

◆ getIntPtFractureAperture()

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

◆ getIntPtFractureStress()

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

◆ getIntPtSigma()

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

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

Definition at line 108 of file SmallDeformationLocalAssemblerFracture.h.

113 {
114 cache.resize(0);
115 return cache;
116 }

◆ getShapeMatrix()

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

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 99 of file SmallDeformationLocalAssemblerFracture.h.

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

References _secondary_data.

◆ preTimestepConcrete()

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

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 83 of file SmallDeformationLocalAssemblerFracture.h.

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

References _integration_method, and _ip_data.

Member Data Documentation

◆ _element

template<typename ShapeFunction, int DisplacementDim>
MeshLib::Element const& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::_element
private

◆ _fracID_to_local

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

◆ _fracture_property

template<typename ShapeFunction, int DisplacementDim>
FractureProperty const* ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::_fracture_property = nullptr
private

◆ _fracture_props

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

◆ _integration_method

◆ _ip_data

◆ _junction_props

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

◆ _process_data

template<typename ShapeFunction, int DisplacementDim>
SmallDeformationProcessData<DisplacementDim>& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::_process_data
private

◆ _secondary_data

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

◆ _shape_matrices

template<typename ShapeFunction, int DisplacementDim>
std::vector<ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices> > ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, DisplacementDim >::_shape_matrices
private

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