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

Detailed Description

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

Definition at line 37 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

#include <SmallDeformationLocalAssemblerMatrixNearFracture.h>

Inheritance diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< 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 BMatricesType = BMatrixPolicyType< ShapeFunction, DisplacementDim >
 
using BMatrixType = typename BMatricesType::BMatrixType
 
using StiffnessMatrixType = typename BMatricesType::StiffnessMatrixType
 
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
 
using NodalDisplacementVectorType = typename BMatricesType::NodalForceVectorType
 

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, 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 (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. 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
 

Private Member Functions

std::vector< double > const & getIntPtSigma (std::vector< double > &cache, std::size_t const component) const
 
std::vector< double > const & getIntPtEpsilon (std::vector< double > &cache, std::size_t const component) 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< IntegrationPointDataMatrix< ShapeMatricesType, BMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointDataMatrix< ShapeMatricesType, BMatricesType, DisplacementDim > > > _ip_data
 
IntegrationMethod _integration_method
 
MeshLib::Element const & _element
 
bool const _is_axially_symmetric
 
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
 

Additional Inherited Members

Member Typedef Documentation

◆ BMatricesType

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

◆ BMatrixType

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

◆ NodalDisplacementVectorType

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

◆ NodalForceVectorType

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

◆ NodalMatrixType

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

◆ NodalVectorType

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

◆ ShapeMatrices

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

◆ ShapeMatricesType

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

◆ StiffnessMatrixType

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

Constructor & Destructor Documentation

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [1/3]

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

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [2/3]

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

◆ SmallDeformationLocalAssemblerMatrixNearFracture() [3/3]

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

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

55  n_variables * ShapeFunction::NPOINTS * DisplacementDim,
56  dofIndex_to_localIndex),
57  _process_data(process_data),
58  _integration_method(integration_order),
59  _element(e),
60  _is_axially_symmetric(is_axially_symmetric)
61 {
62  std::vector<ShapeMatrices, Eigen::aligned_allocator<
64  shape_matrices =
66  DisplacementDim>(e, is_axially_symmetric,
68 
69  unsigned const n_integration_points =
70  _integration_method.getNumberOfPoints();
71 
72  _ip_data.reserve(n_integration_points);
73  _secondary_data.N.resize(n_integration_points);
74 
76  _process_data.solid_materials, _process_data.material_ids, e.getID());
77 
78  for (unsigned ip = 0; ip < n_integration_points; ip++)
79  {
80  _ip_data.emplace_back(solid_material);
81  auto& ip_data = _ip_data[ip];
82  auto const& sm = shape_matrices[ip];
83  ip_data.N = sm.N;
84  ip_data.dNdx = sm.dNdx;
85  ip_data.integration_weight =
86  _integration_method.getWeightedPoint(ip).getWeight() *
87  sm.integralMeasure * sm.detJ;
88 
89  // Initialize current time step values
90  static const int kelvin_vector_size =
92  ip_data._sigma.setZero(kelvin_vector_size);
93  ip_data._eps.setZero(kelvin_vector_size);
94 
95  // Previous time step values are not initialized and are set later.
96  ip_data._sigma_prev.resize(kelvin_vector_size);
97  ip_data._eps_prev.resize(kelvin_vector_size);
98 
99  ip_data._C.resize(kelvin_vector_size, kelvin_vector_size);
100 
101  _secondary_data.N[ip] = sm.N;
102  }
103 
104  for (auto fid : process_data._vec_ele_connected_fractureIDs[e.getID()])
105  {
106  _fracID_to_local.insert({fid, _fracture_props.size()});
107  _fracture_props.push_back(&_process_data.fracture_properties[fid]);
108  }
109 
110  for (auto jid : process_data._vec_ele_connected_junctionIDs[e.getID()])
111  {
112  _junction_props.push_back(&_process_data.junction_properties[jid]);
113  }
114 }
std::vector< IntegrationPointDataMatrix< ShapeMatricesType, BMatricesType, DisplacementDim >, Eigen::aligned_allocator< IntegrationPointDataMatrix< ShapeMatricesType, BMatricesType, DisplacementDim > > > _ip_data
MechanicsBase< DisplacementDim > & selectSolidConstitutiveRelation(std::map< int, std::unique_ptr< MechanicsBase< DisplacementDim >>> 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.
Definition: KelvinVector.h:23
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)
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
Definition: SecondaryData.h:28

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_fracID_to_local, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_fracture_props, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_integration_method, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_ip_data, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_junction_props, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_process_data, ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_secondary_data, ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::_vec_ele_connected_fractureIDs, ProcessLib::LIE::SmallDeformation::SmallDeformationProcessData< DisplacementDim >::_vec_ele_connected_junctionIDs, MeshLib::Element::getID(), NumLib::initShapeMatrices(), MathLib::KelvinVector::kelvin_vector_dimensions(), ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N, and MaterialLib::Solids::selectSolidConstitutiveRelation().

Member Function Documentation

◆ assemble()

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

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

References OGS_FATAL.

◆ assembleWithJacobian()

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

124 {
125  assert(_element.getDimension() == DisplacementDim);
126 
127  auto const N_DOF_PER_VAR = ShapeFunction::NPOINTS * DisplacementDim;
128  auto const n_fractures = _fracture_props.size();
129  auto const n_junctions = _junction_props.size();
130  auto const n_enrich_var = n_fractures + n_junctions;
131 
132  using BlockVectorType =
133  typename Eigen::VectorXd::FixedSegmentReturnType<N_DOF_PER_VAR>::Type;
134  using BlockMatrixType =
135  Eigen::Block<Eigen::MatrixXd, N_DOF_PER_VAR, N_DOF_PER_VAR>;
136 
137  //--------------------------------------------------------------------------------------
138  // prepare sub vectors, matrices for regular displacement (u) and
139  // displacement jumps (g)
140  //
141  // example with two fractures with one intersection:
142  // |b(u)|
143  // b = |b(g1)|
144  // |b(g2)|
145  // |b(j1)|
146  //
147  // |J(u,u) J(u,g1) J(u,g2) J(u,j1) |
148  // J = |J(g1,u) J(g1,g1) J(g1,g2) J(g1,j1)|
149  // |J(g2,u) J(g2,g1) J(g2,g2) J(g2,j1)|
150  // |J(j1,u) J(j1,g1) J(j1,g2) J(j1,j1)|
151  //--------------------------------------------------------------------------------------
152  auto local_b_u = local_b.segment<N_DOF_PER_VAR>(0);
153  std::vector<BlockVectorType> vec_local_b_g;
154  for (unsigned i = 0; i < n_enrich_var; i++)
155  {
156  vec_local_b_g.push_back(
157  local_b.segment<N_DOF_PER_VAR>(N_DOF_PER_VAR * (i + 1)));
158  }
159 
160  auto local_J_uu = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(0, 0);
161  std::vector<BlockMatrixType> vec_local_J_ug;
162  std::vector<BlockMatrixType> vec_local_J_gu;
163  std::vector<std::vector<BlockMatrixType>> vec_local_J_gg(n_enrich_var);
164  for (unsigned i = 0; i < n_enrich_var; i++)
165  {
166  auto sub_ug = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(
167  0, N_DOF_PER_VAR * (i + 1));
168  vec_local_J_ug.push_back(sub_ug);
169 
170  auto sub_gu = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(
171  N_DOF_PER_VAR * (i + 1), 0);
172  vec_local_J_gu.push_back(sub_gu);
173 
174  for (unsigned j = 0; j < n_enrich_var; j++)
175  {
176  auto sub_gg = local_J.block<N_DOF_PER_VAR, N_DOF_PER_VAR>(
177  N_DOF_PER_VAR * (i + 1), N_DOF_PER_VAR * (j + 1));
178  vec_local_J_gg[i].push_back(sub_gg);
179  }
180  }
181 
182  auto const nodal_u = local_u.segment<N_DOF_PER_VAR>(0);
183  std::vector<BlockVectorType> vec_nodal_g;
184  for (unsigned i = 0; i < n_enrich_var; i++)
185  {
186  auto sub = const_cast<Eigen::VectorXd&>(local_u).segment<N_DOF_PER_VAR>(
187  N_DOF_PER_VAR * (i + 1));
188  vec_nodal_g.push_back(sub);
189  }
190 
191  //------------------------------------------------
192  // integration
193  //------------------------------------------------
194  unsigned const n_integration_points =
195  _integration_method.getNumberOfPoints();
196 
197  MPL::VariableArray variables;
198  MPL::VariableArray variables_prev;
200  x_position.setElementID(_element.getID());
201 
202  for (unsigned ip = 0; ip < n_integration_points; ip++)
203  {
204  x_position.setIntegrationPoint(ip);
205 
206  auto& ip_data = _ip_data[ip];
207  auto const& w = _ip_data[ip].integration_weight;
208 
209  auto const& N = ip_data.N;
210  auto const& dNdx = ip_data.dNdx;
211 
212  // levelset functions
213  Eigen::Vector3d const ip_physical_coords(
214  computePhysicalCoordinates(_element, N).getCoords());
215  std::vector<double> const levelsets(
217  _fracID_to_local, ip_physical_coords));
218 
219  // u = u^hat + sum_i(enrich^br_i(x) * [u]_i) + sum_i(enrich^junc_i(x) *
220  // [u]_i)
221  NodalDisplacementVectorType nodal_total_u = nodal_u;
222  for (unsigned i = 0; i < n_enrich_var; i++)
223  {
224  nodal_total_u += levelsets[i] * vec_nodal_g[i];
225  }
226 
227  auto const x_coord =
228  NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
229  _element, N);
230  auto const B =
231  LinearBMatrix::computeBMatrix<DisplacementDim,
232  ShapeFunction::NPOINTS,
233  typename BMatricesType::BMatrixType>(
234  dNdx, N, x_coord, _is_axially_symmetric);
235 
236  // strain, stress
237  auto const& eps_prev = ip_data._eps_prev;
238  auto const& sigma_prev = ip_data._sigma_prev;
239 
240  auto& eps = ip_data._eps;
241  auto& sigma = ip_data._sigma;
242  auto& state = ip_data._material_state_variables;
243 
244  eps.noalias() = B * nodal_total_u;
245 
246  variables[static_cast<int>(
249  eps);
250 
251  variables_prev[static_cast<int>(MaterialPropertyLib::Variable::stress)]
253  sigma_prev);
254  variables_prev[static_cast<int>(
257  eps_prev);
258  variables_prev[static_cast<int>(
260  .emplace<double>(_process_data._reference_temperature);
261 
262  auto&& solution = _ip_data[ip]._solid_material.integrateStress(
263  variables_prev, variables, t, x_position, dt, *state);
264 
265  if (!solution)
266  {
267  OGS_FATAL("Computation of local constitutive relation failed.");
268  }
269 
271  std::tie(sigma, state, C) = std::move(*solution);
272 
273  // r_u = B^T * Sigma = B^T * C * B * (u+phi*[u])
274  // r_[u] = (phi*B)^T * Sigma = (phi*B)^T * C * B * (u+phi*[u])
275  local_b_u.noalias() -= B.transpose() * sigma * w;
276  for (unsigned i = 0; i < n_enrich_var; i++)
277  {
278  vec_local_b_g[i].noalias() -=
279  levelsets[i] * B.transpose() * sigma * w;
280  }
281 
282  // J_uu += B^T * C * B
283  local_J_uu.noalias() += B.transpose() * C * B * w;
284 
285  for (unsigned i = 0; i < n_enrich_var; i++)
286  {
287  // J_u[u] += B^T * C * (levelset * B)
288  vec_local_J_ug[i].noalias() +=
289  B.transpose() * C * (levelsets[i] * B) * w;
290 
291  // J_[u]u += (levelset * B)^T * C * B
292  vec_local_J_gu[i].noalias() +=
293  (levelsets[i] * B.transpose()) * C * B * w;
294 
295  for (unsigned j = 0; j < n_enrich_var; j++)
296  {
297  // J_[u][u] += (levelset * B)^T * C * (levelset * B)
298  vec_local_J_gg[i][j].noalias() +=
299  (levelsets[i] * B.transpose()) * C * (levelsets[j] * B) * w;
300  }
301  }
302  }
303 }
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)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
Definition: BMatrixPolicy.h:55
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Definition: VariableType.h:108
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:48
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), kelvin_vector_dimensions(DisplacementDim), Eigen::RowMajor > KelvinMatrixType
Definition: KelvinVector.h:56
MathLib::Point3d 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)
BMatrixType computeBMatrix(DNDX_Type const &dNdx, N_Type const &N, const double radius, const bool is_axially_symmetric)
Fills a B-matrix based on given shape function dN/dx values.
Definition: LinearBMatrix.h:42

References ProcessLib::LinearBMatrix::computeBMatrix(), ProcessLib::LIE::computePhysicalCoordinates(), MaterialPropertyLib::mechanical_strain, OGS_FATAL, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::stress, MaterialPropertyLib::temperature, and ProcessLib::LIE::uGlobalEnrichments().

◆ computeSecondaryVariableConcreteWithVector()

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

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

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

313 {
314  // Compute average value per element
315  const int n = DisplacementDim == 2 ? 4 : 6;
316  Eigen::VectorXd ele_stress = Eigen::VectorXd::Zero(n);
317  Eigen::VectorXd ele_strain = Eigen::VectorXd::Zero(n);
318 
319  unsigned const n_integration_points =
320  _integration_method.getNumberOfPoints();
321  for (unsigned ip = 0; ip < n_integration_points; ip++)
322  {
323  auto& ip_data = _ip_data[ip];
324 
325  ele_stress += ip_data._sigma;
326  ele_strain += ip_data._eps;
327  }
328  ele_stress /= n_integration_points;
329  ele_strain /= n_integration_points;
330 
331  (*_process_data._mesh_prop_stress_xx)[_element.getID()] = ele_stress[0];
332  (*_process_data._mesh_prop_stress_yy)[_element.getID()] = ele_stress[1];
333  (*_process_data._mesh_prop_stress_zz)[_element.getID()] = ele_stress[2];
334  (*_process_data._mesh_prop_stress_xy)[_element.getID()] = ele_stress[3];
335  if (DisplacementDim == 3)
336  {
337  (*_process_data._mesh_prop_stress_yz)[_element.getID()] = ele_stress[4];
338  (*_process_data._mesh_prop_stress_xz)[_element.getID()] = ele_stress[5];
339  }
340 
341  (*_process_data._mesh_prop_strain_xx)[_element.getID()] = ele_strain[0];
342  (*_process_data._mesh_prop_strain_yy)[_element.getID()] = ele_strain[1];
343  (*_process_data._mesh_prop_strain_zz)[_element.getID()] = ele_strain[2];
344  (*_process_data._mesh_prop_strain_xy)[_element.getID()] = ele_strain[3];
345  if (DisplacementDim == 3)
346  {
347  (*_process_data._mesh_prop_strain_yz)[_element.getID()] = ele_strain[4];
348  (*_process_data._mesh_prop_strain_xz)[_element.getID()] = ele_strain[5];
349  }
350 }

◆ getIntPtEpsilon()

template<typename ShapeFunction , typename IntegrationMethod , int DisplacementDim>
std::vector<double> const& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilon ( std::vector< double > &  cache,
std::size_t const  component 
) const
inlineprivate

Definition at line 244 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

246  {
247  cache.clear();
248  cache.reserve(_ip_data.size());
249 
250  for (auto const& ip_data : _ip_data)
251  {
252  if (component < 3)
253  { // xx, yy, zz components
254  cache.push_back(ip_data._eps[component]);
255  }
256  else
257  { // mixed xy, yz, xz components
258  cache.push_back(ip_data._eps[component] / std::sqrt(2));
259  }
260  }
261 
262  return cache;
263  }

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_ip_data.

Referenced by ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonXX(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonXY(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonXZ(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonYY(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonYZ(), and ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonZZ().

◆ getIntPtEpsilonXX()

template<typename ShapeFunction , typename IntegrationMethod , int DisplacementDim>
std::vector<double> const& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtEpsilonZZ ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &  cache 
) const
inlineoverridevirtual

◆ getIntPtSigma()

template<typename ShapeFunction , typename IntegrationMethod , int DisplacementDim>
std::vector<double> const& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigma ( std::vector< double > &  cache,
std::size_t const  component 
) const
inlineprivate

Definition at line 223 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

225  {
226  cache.clear();
227  cache.reserve(_ip_data.size());
228 
229  for (auto const& ip_data : _ip_data)
230  {
231  if (component < 3)
232  { // xx, yy, zz components
233  cache.push_back(ip_data._sigma[component]);
234  }
235  else
236  { // mixed xy, yz, xz components
237  cache.push_back(ip_data._sigma[component] / std::sqrt(2));
238  }
239  }
240 
241  return cache;
242  }

References ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::_ip_data.

Referenced by ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaXX(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaXY(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaXZ(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaYY(), ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaYZ(), and ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< ShapeFunction, IntegrationMethod, DisplacementDim >::getIntPtSigmaZZ().

◆ getIntPtSigmaXX()

template<typename ShapeFunction , typename IntegrationMethod , int DisplacementDim>
std::vector<double> const& ProcessLib::LIE::SmallDeformation::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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::SmallDeformationLocalAssemblerMatrixNearFracture< 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 101 of file SmallDeformationLocalAssemblerMatrixNearFracture.h.

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

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

◆ preTimestepConcrete()

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

Member Data Documentation

◆ _element

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

◆ _fracID_to_local

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

◆ _fracture_props

◆ _integration_method

◆ _ip_data

◆ _is_axially_symmetric

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

◆ _junction_props

◆ _process_data

◆ _secondary_data


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