OGS
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim > Class Template Reference

Detailed Description

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

Definition at line 34 of file SmallDeformationLocalAssemblerFracture.h.

#include <SmallDeformationLocalAssemblerFracture.h>

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

Public Types

using ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, DisplacementDim >
 
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 = typename HMatricesType::NodalForceVectorType
 
using ForceVectorType = typename HMatricesType::ForceVectorType
 

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, bool const is_axially_symmetric, unsigned const integration_order, 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. More...
 
std::vector< double > const & getIntPtSigmaXX (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtSigmaYY (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtSigmaZZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtSigmaXY (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtSigmaXZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtSigmaYZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonXX (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonYY (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonZZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonXY (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonXZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilonYZ (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 &, const double, const double, std::vector< double > &, std::vector< double > &, 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
 
void setInitialConditions (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, 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_xdot, 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_xdot, const double dxdot_dx, const double dx_dx, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, 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_dot, 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, double const t, double const dt)
 
void postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id)
 
virtual std::vector< double > interpolateNodalValuesToIntegrationPoints (std::vector< double > const &)
 
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. More...
 
- Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual ~ExtrapolatableElement ()=default
 

Static Private Member Functions

static std::vector< double > const & getIntPtSigma (std::vector< double > &cache, std::size_t const)
 
static std::vector< double > const & getIntPtEpsilon (std::vector< double > &cache, std::size_t 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
 
FractureProperty const * _fracture_property = nullptr
 
std::vector< IntegrationPointDataFracture< HMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointDataFracture< HMatricesType, DisplacementDim > > > _ip_data
 
IntegrationMethod _integration_method
 
std::vector< ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > _shape_matrices
 
MeshLib::Element const & _element
 
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
 

Additional Inherited Members

Member Typedef Documentation

◆ ForceVectorType

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

Definition at line 51 of file SmallDeformationLocalAssemblerFracture.h.

◆ HMatricesType

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

Definition at line 43 of file SmallDeformationLocalAssemblerFracture.h.

◆ HMatrixType

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

Definition at line 45 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalDisplacementVectorType

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

Definition at line 48 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalForceVectorType

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

Definition at line 47 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalMatrixType

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

Definition at line 40 of file SmallDeformationLocalAssemblerFracture.h.

◆ NodalVectorType

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

Definition at line 41 of file SmallDeformationLocalAssemblerFracture.h.

◆ ShapeMatrices

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

Definition at line 42 of file SmallDeformationLocalAssemblerFracture.h.

◆ ShapeMatricesType

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

Definition at line 38 of file SmallDeformationLocalAssemblerFracture.h.

◆ StiffnessMatrixType

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

Definition at line 46 of file SmallDeformationLocalAssemblerFracture.h.

Constructor & Destructor Documentation

◆ SmallDeformationLocalAssemblerFracture() [1/3]

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

◆ SmallDeformationLocalAssemblerFracture() [2/3]

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

◆ SmallDeformationLocalAssemblerFracture() [3/3]

template<typename ShapeFunction , typename IntegrationMethod , int DisplacementDim>
ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, 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,
bool const  is_axially_symmetric,
unsigned const  integration_order,
SmallDeformationProcessData< DisplacementDim > &  process_data 
)

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

40  n_variables * ShapeFunction::NPOINTS * DisplacementDim,
41  dofIndex_to_localIndex),
42  _process_data(process_data),
43  _integration_method(integration_order),
46  DisplacementDim>(e, is_axially_symmetric,
48  _element(e)
49 {
50  assert(_element.getDimension() == DisplacementDim - 1);
51 
52  unsigned const n_integration_points =
53  _integration_method.getNumberOfPoints();
54 
55  _ip_data.reserve(n_integration_points);
56  _secondary_data.N.resize(n_integration_points);
57 
58  auto mat_id = (*_process_data._mesh_prop_materialIDs)[e.getID()];
59  auto frac_id = _process_data._map_materialID_to_fractureID[mat_id];
60  _fracture_property = &_process_data.fracture_properties[frac_id];
61  for (auto fid : process_data._vec_ele_connected_fractureIDs[e.getID()])
62  {
63  _fracID_to_local.insert({fid, _fracture_props.size()});
64  _fracture_props.push_back(&_process_data.fracture_properties[fid]);
65  }
66 
67  for (auto jid : process_data._vec_ele_connected_junctionIDs[e.getID()])
68  {
69  _junction_props.push_back(&_process_data.junction_properties[jid]);
70  }
71 
73  x_position.setElementID(_element.getID());
74  for (unsigned ip = 0; ip < n_integration_points; ip++)
75  {
76  x_position.setIntegrationPoint(ip);
77 
78  _ip_data.emplace_back(*_process_data._fracture_model);
79  auto const& sm = _shape_matrices[ip];
80  auto& ip_data = _ip_data[ip];
81  ip_data.integration_weight =
82  _integration_method.getWeightedPoint(ip).getWeight() *
83  sm.integralMeasure * sm.detJ;
84  ip_data.h_matrices.setZero(DisplacementDim,
85  ShapeFunction::NPOINTS * DisplacementDim);
86 
87  computeHMatrix<DisplacementDim, ShapeFunction::NPOINTS,
89  HMatrixType>(sm.N, ip_data.h_matrices);
90 
91  // Initialize current time step values
92  ip_data.w.setZero(DisplacementDim);
93  ip_data.sigma.setZero(DisplacementDim);
94 
95  // Previous time step values are not initialized and are set later.
96  ip_data.sigma_prev.resize(DisplacementDim);
97  ip_data.w_prev.resize(DisplacementDim);
98 
99  ip_data.C.resize(DisplacementDim, DisplacementDim);
100 
101  ip_data.aperture0 = _fracture_property->aperture0(0, x_position)[0];
102  ip_data.aperture_prev = ip_data.aperture0;
103 
104  _secondary_data.N[ip] = sm.N;
105  }
106 }
virtual constexpr unsigned getDimension() const =0
Get dimension of the mesh element.
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:82
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
std::vector< IntegrationPointDataFracture< HMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointDataFracture< HMatricesType, DisplacementDim > > > _ip_data
std::vector< ShapeMatrices, Eigen::aligned_allocator< typename ShapeMatricesType::ShapeMatrices > > _shape_matrices
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)
void computeHMatrix(N_Type const &N, HMatrixType &H)
Fills a H-matrix based on given shape function.
Definition: HMatrixUtils.h:53
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
Definition: SecondaryData.h:28
ParameterLib::Parameter< double > const & aperture0
Initial aperture.

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_element, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_fracID_to_local, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_fracture_property, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_fracture_props, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_integration_method, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_ip_data, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_junction_props, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_process_data, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_secondary_data, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_shape_matrices, ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::_vec_ele_connected_fractureIDs, ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::_vec_ele_connected_junctionIDs, ProcessLib::LIE::FractureProperty::aperture0, ProcessLib::computeHMatrix(), MeshLib::Element::getDimension(), MeshLib::Element::getID(), ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N, ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

Member Function Documentation

◆ assemble()

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

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

References OGS_FATAL.

◆ assembleWithJacobian()

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

116 {
117  auto const N_DOF_PER_VAR = ShapeFunction::NPOINTS * DisplacementDim;
118  auto const n_fractures = _fracture_props.size();
119  auto const n_junctions = _junction_props.size();
120  auto const n_enrich_var = n_fractures + n_junctions;
121 
122  //--------------------------------------------------------------------------------------
123  // prepare sub vectors, matrices for regular displacement (u) and
124  // displacement jumps (g)
125  //
126  // example with two fractures with one intersection:
127  // b = |b(g1)|
128  // |b(g2)|
129  //
130  // J = |J(g1,g1) J(g1,g2)|
131  // |J(g2,g1) J(g2,g2)|
132  //--------------------------------------------------------------------------------------
133 
134  using BlockVectorType =
135  typename Eigen::VectorXd::FixedSegmentReturnType<N_DOF_PER_VAR>::Type;
136  using BlockMatrixType =
137  Eigen::Block<Eigen::MatrixXd, N_DOF_PER_VAR, N_DOF_PER_VAR>;
138 
139  std::vector<BlockVectorType> vec_local_b_g;
140  for (unsigned i = 0; i < n_enrich_var; i++)
141  {
142  vec_local_b_g.push_back(
143  local_b.segment<N_DOF_PER_VAR>(N_DOF_PER_VAR * i));
144  }
145  std::vector<std::vector<BlockMatrixType>> vec_local_J_gg(n_enrich_var);
146  for (unsigned i = 0; i < n_enrich_var; i++)
147  {
148  for (unsigned j = 0; j < n_enrich_var; j++)
149  {
150  auto sub_gg = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(
151  N_DOF_PER_VAR * i, N_DOF_PER_VAR * j);
152  vec_local_J_gg[i].push_back(sub_gg);
153  }
154  }
155 
156  std::vector<Eigen::VectorXd> vec_nodal_g;
157  for (unsigned i = 0; i < n_enrich_var; i++)
158  {
159  auto sub = const_cast<Eigen::VectorXd&>(local_u).segment<N_DOF_PER_VAR>(
160  N_DOF_PER_VAR * i);
161  vec_nodal_g.push_back(sub);
162  }
163 
164  //------------------------------------------------
165  // integration
166  //------------------------------------------------
167  // the index of a normal (normal to a fracture plane) component
168  // in a displacement vector
169  int const index_normal = DisplacementDim - 1;
170  auto const& R = _fracture_property->R;
171 
172  unsigned const n_integration_points =
173  _integration_method.getNumberOfPoints();
174 
176  x_position.setElementID(_element.getID());
177 
178  for (unsigned ip = 0; ip < n_integration_points; ip++)
179  {
180  x_position.setIntegrationPoint(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& N = _secondary_data.N[ip];
193 
194  Eigen::Vector3d const ip_physical_coords(
195  computePhysicalCoordinates(_element, N).getCoords());
196  std::vector<double> const levelsets(duGlobalEnrichments(
198  _fracID_to_local, ip_physical_coords));
199 
200  // du = du^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x)
201  // * [u]_i)
202  Eigen::VectorXd nodal_gap(N_DOF_PER_VAR);
203  nodal_gap.setZero();
204  for (unsigned i = 0; i < n_enrich_var; i++)
205  {
206  nodal_gap += levelsets[i] * vec_nodal_g[i];
207  }
208 
209  // displacement jumps
210  w.noalias() = R * H * nodal_gap;
211 
212  // total aperture
213  ip_data.aperture = ip_data.aperture0 + w[index_normal];
214 
215  // local C, local stress
216  mat.computeConstitutiveRelation(
217  t, x_position, ip_data.aperture0,
218  Eigen::Matrix<double, DisplacementDim, 1>::Zero(), // TODO (naumov)
219  // Replace with
220  // initial
221  // stress values
222  w_prev, w, sigma_prev, sigma, C, state);
223 
224  // r_[u] += H^T*Stress
225  for (unsigned i = 0; i < n_enrich_var; i++)
226  {
227  vec_local_b_g[i].noalias() -= levelsets[i] * H.transpose() *
228  R.transpose() * sigma *
229  integration_weight;
230  }
231 
232  // J_[u][u] += H^T*C*H
233  for (unsigned i = 0; i < n_enrich_var; i++)
234  {
235  for (unsigned j = 0; j < n_enrich_var; j++)
236  {
237  // J_[u][u] += (levelset * B)^T * C * (levelset * B)
238  vec_local_J_gg[i][j].noalias() +=
239  (levelsets[i] * H.transpose() * R.transpose()) * C *
240  (levelsets[j] * R * H) * integration_weight;
241  }
242  }
243  }
244 }
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)
MathLib::Point3d computePhysicalCoordinates(MeshLib::Element const &e, Eigen::MatrixBase< Derived > const &shape)
Definition: Utils.h:24
Eigen::MatrixXd R
Rotation matrix from global to local coordinates.

References ProcessLib::LIE::computePhysicalCoordinates(), ProcessLib::LIE::duGlobalEnrichments(), ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

◆ computeSecondaryVariableConcreteWithVector()

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

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

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

252 {
253  auto const N_DOF_PER_VAR = ShapeFunction::NPOINTS * DisplacementDim;
254  auto const n_fractures = _fracture_props.size();
255  auto const n_junctions = _junction_props.size();
256  auto const n_enrich_var = n_fractures + n_junctions;
257 
258  auto const& R = _fracture_property->R;
259 
260  // the index of a normal (normal to a fracture plane) component
261  // in a displacement vector
262  int const index_normal = DisplacementDim - 1;
263 
264  unsigned const n_integration_points =
265  _integration_method.getNumberOfPoints();
266 
268  x_position.setElementID(_element.getID());
269 
270  std::vector<Eigen::VectorXd> vec_nodal_g;
271  for (unsigned i = 0; i < n_enrich_var; i++)
272  {
273  auto sub = const_cast<Eigen::VectorXd&>(local_u).segment<N_DOF_PER_VAR>(
274  N_DOF_PER_VAR * i);
275  vec_nodal_g.push_back(sub);
276  }
277 
278  for (unsigned ip = 0; ip < n_integration_points; ip++)
279  {
280  x_position.setIntegrationPoint(ip);
281 
282  auto& ip_data = _ip_data[ip];
283  auto const& H = ip_data.h_matrices;
284  auto& mat = ip_data.fracture_material;
285  auto& sigma = ip_data.sigma;
286  auto const& sigma_prev = ip_data.sigma_prev;
287  auto& w = ip_data.w;
288  auto const& w_prev = ip_data.w_prev;
289  auto& C = ip_data.C;
290  auto& state = *ip_data.material_state_variables;
291  auto& b_m = ip_data.aperture;
292  auto& N = _secondary_data.N[ip];
293 
294  Eigen::Vector3d const ip_physical_coords(
295  computePhysicalCoordinates(_element, N).getCoords());
296  std::vector<double> const levelsets(duGlobalEnrichments(
298  _fracID_to_local, ip_physical_coords));
299 
300  // du = du^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x)
301  // * [u]_i)
302  Eigen::VectorXd nodal_gap(N_DOF_PER_VAR);
303  nodal_gap.setZero();
304  for (unsigned i = 0; i < n_enrich_var; i++)
305  {
306  nodal_gap += levelsets[i] * vec_nodal_g[i];
307  }
308 
309  // displacement jumps in local coordinates
310  w.noalias() = R * H * nodal_gap;
311 
312  // aperture
313  b_m = ip_data.aperture0 + w[index_normal];
314  if (b_m < 0.0)
315  {
316  OGS_FATAL(
317  "Element {:d}, gp {:d}: Fracture aperture is {:g}, but it must "
318  "be non-negative.",
319  _element.getID(), ip, b_m);
320  }
321 
322  // local C, local stress
323  mat.computeConstitutiveRelation(
324  t, x_position, ip_data.aperture0,
325  Eigen::Matrix<double, DisplacementDim, 1>::Zero(), // TODO (naumov)
326  // Replace with
327  // initial
328  // stress values
329  w_prev, w, sigma_prev, sigma, C, state);
330  }
331 
332  double ele_b = 0;
333  typename HMatricesType::ForceVectorType ele_sigma =
334  HMatricesType::ForceVectorType::Zero(DisplacementDim);
335  typename HMatricesType::ForceVectorType ele_w =
336  HMatricesType::ForceVectorType::Zero(DisplacementDim);
337 
338  for (unsigned ip = 0; ip < n_integration_points; ip++)
339  {
340  auto& ip_data = _ip_data[ip];
341  ele_b += ip_data.aperture;
342  ele_w += ip_data.w;
343  ele_sigma += ip_data.sigma;
344  }
345  ele_b /= n_integration_points;
346  ele_w /= n_integration_points;
347  ele_sigma /= n_integration_points;
348  (*_process_data._mesh_prop_b)[_element.getID()] = ele_b;
349  (*_process_data._mesh_prop_w_n)[_element.getID()] = ele_w[index_normal];
350  (*_process_data._mesh_prop_w_s)[_element.getID()] = ele_w[0];
351  (*_process_data._mesh_prop_fracture_stress_normal)[_element.getID()] =
352  ele_sigma[index_normal];
353  (*_process_data._mesh_prop_fracture_stress_shear)[_element.getID()] =
354  ele_sigma[0];
355 }
VectorType< DisplacementDim > ForceVectorType
Definition: HMatrixUtils.h:45

References ProcessLib::LIE::computePhysicalCoordinates(), ProcessLib::LIE::duGlobalEnrichments(), OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), and ParameterLib::SpatialPosition::setIntegrationPoint().

◆ getIntPtEpsilon()

◆ getIntPtEpsilonXX()

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

◆ getIntPtEpsilonXY()

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

◆ getIntPtEpsilonXZ()

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

◆ getIntPtEpsilonYY()

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

◆ getIntPtEpsilonYZ()

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

◆ getIntPtEpsilonZZ()

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

◆ getIntPtSigma()

◆ getIntPtSigmaXX()

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

◆ getIntPtSigmaXY()

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

◆ getIntPtSigmaXZ()

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

◆ getIntPtSigmaYY()

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

◆ getIntPtSigmaYZ()

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

◆ getIntPtSigmaZZ()

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

◆ getShapeMatrix()

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

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 100 of file SmallDeformationLocalAssemblerFracture.h.

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

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_secondary_data, and ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N.

◆ preTimestepConcrete()

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

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 84 of file SmallDeformationLocalAssemblerFracture.h.

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

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_integration_method, and ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_ip_data.

Member Data Documentation

◆ _element

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

◆ _fracID_to_local

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

◆ _fracture_property

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

◆ _fracture_props

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

◆ _integration_method

◆ _ip_data

◆ _junction_props

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

◆ _process_data

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

◆ _secondary_data

◆ _shape_matrices

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

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