OGS
ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, int DisplacementDim>
class ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >

Definition at line 134 of file PhaseFieldFEM.h.

#include <PhaseFieldFEM.h>

Inheritance diagram for ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]
Collaboration diagram for ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >:
[legend]

Public Types

using ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, DisplacementDim >
 
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
 
using BMatricesType = BMatrixPolicyType< ShapeFunction, DisplacementDim >
 
using NodalForceVectorType = typename BMatricesType::NodalForceVectorType
 
using DeformationVector = typename ShapeMatricesType::template VectorType< displacement_size >
 
using PhaseFieldVector = typename ShapeMatricesType::template VectorType< phasefield_size >
 
using PhaseFieldMatrix = typename ShapeMatricesType::template MatrixType< phasefield_size, phasefield_size >
 
using NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
 
using NodalVectorType = typename ShapeMatricesType::NodalVectorType
 
using IpData = IntegrationPointData< BMatricesType, ShapeMatricesType, DisplacementDim >
 

Public Member Functions

 PhaseFieldLocalAssembler (PhaseFieldLocalAssembler const &)=delete
 
 PhaseFieldLocalAssembler (PhaseFieldLocalAssembler &&)=delete
 
 PhaseFieldLocalAssembler (MeshLib::Element const &e, std::size_t const, NumLib::GenericIntegrationMethod const &integration_method, bool const is_axially_symmetric, PhaseFieldProcessData< 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 assembleWithJacobianForStaggeredScheme (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, std::vector< double > &local_Jac_data) override
 
void initializeConcrete () override
 
void postTimestepConcrete (Eigen::VectorXd const &, double const, double const) override
 
void computeCrackIntegral (std::size_t mesh_item_id, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, GlobalVector const &x, double const t, double &crack_volume, CoupledSolutionsForStaggeredScheme const *const cpl_xs) override
 
void computeEnergy (std::size_t mesh_item_id, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, GlobalVector const &x, double const t, double &elastic_energy, double &surface_energy, double &pressure_work, CoupledSolutionsForStaggeredScheme const *const cpl_xs) override
 
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix (const unsigned integration_point) const override
 Provides the shape matrix at the given integration point. More...
 
- 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 assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, 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 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 (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
std::vector< double > const & getIntPtEpsilon (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &cache) const override
 
void assembleWithJacobianPhaseFieldEquations (double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
 
void assembleWithJacobianForDeformationEquations (double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
 

Private Attributes

PhaseFieldProcessData< DisplacementDim > & _process_data
 
std::vector< IpData, Eigen::aligned_allocator< IpData > > _ip_data
 
NumLib::GenericIntegrationMethod const & _integration_method
 
MeshLib::Element const & _element
 
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
 
bool const _is_axially_symmetric
 

Static Private Attributes

static constexpr int displacement_index = 0
 
static constexpr int displacement_size
 
static constexpr int phasefield_index
 
static constexpr int phasefield_size = ShapeFunction::NPOINTS
 
static const int phase_process_id = 1
 
static const int mechanics_process_id = 0
 

Member Typedef Documentation

◆ BMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::BMatricesType = BMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 150 of file PhaseFieldFEM.h.

◆ DeformationVector

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::DeformationVector = typename ShapeMatricesType::template VectorType<displacement_size>

Definition at line 154 of file PhaseFieldFEM.h.

◆ IpData

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::IpData = IntegrationPointData<BMatricesType, ShapeMatricesType, DisplacementDim>

Definition at line 163 of file PhaseFieldFEM.h.

◆ NodalForceVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::NodalForceVectorType = typename BMatricesType::NodalForceVectorType

Definition at line 152 of file PhaseFieldFEM.h.

◆ NodalMatrixType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::NodalMatrixType = typename ShapeMatricesType::NodalMatrixType

Definition at line 161 of file PhaseFieldFEM.h.

◆ NodalVectorType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType

Definition at line 162 of file PhaseFieldFEM.h.

◆ PhaseFieldMatrix

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::PhaseFieldMatrix = typename ShapeMatricesType::template MatrixType<phasefield_size, phasefield_size>

Definition at line 158 of file PhaseFieldFEM.h.

◆ PhaseFieldVector

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::PhaseFieldVector = typename ShapeMatricesType::template VectorType<phasefield_size>

Definition at line 156 of file PhaseFieldFEM.h.

◆ ShapeMatrices

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

Definition at line 149 of file PhaseFieldFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunction , int DisplacementDim>
using ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, DisplacementDim>

Definition at line 145 of file PhaseFieldFEM.h.

Constructor & Destructor Documentation

◆ PhaseFieldLocalAssembler() [1/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::PhaseFieldLocalAssembler ( PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim > const &  )
delete

◆ PhaseFieldLocalAssembler() [2/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::PhaseFieldLocalAssembler ( PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim > &&  )
delete

◆ PhaseFieldLocalAssembler() [3/3]

template<typename ShapeFunction , int DisplacementDim>
ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::PhaseFieldLocalAssembler ( MeshLib::Element const &  e,
std::size_t const  ,
NumLib::GenericIntegrationMethod const &  integration_method,
bool const  is_axially_symmetric,
PhaseFieldProcessData< DisplacementDim > &  process_data 
)
inline

Definition at line 169 of file PhaseFieldFEM.h.

175  : _process_data(process_data),
176  _integration_method(integration_method),
177  _element(e),
178  _is_axially_symmetric(is_axially_symmetric)
179  {
180  unsigned const n_integration_points =
182 
183  _ip_data.reserve(n_integration_points);
184  _secondary_data.N.resize(n_integration_points);
185 
186  auto& solid_material =
188  _process_data.solid_materials,
189  _process_data.material_ids,
190  e.getID());
192  DisplacementDim> const*>(&solid_material))
193  {
194  OGS_FATAL(
195  "For phasefield process only linear elastic solid material "
196  "support is implemented.");
197  }
198 
199  auto const shape_matrices =
201  DisplacementDim>(e, is_axially_symmetric,
203 
205  x_position.setElementID(_element.getID());
206 
207  for (unsigned ip = 0; ip < n_integration_points; ip++)
208  {
209  _ip_data.emplace_back(solid_material);
210  auto& ip_data = _ip_data[ip];
211  ip_data.integration_weight =
213  shape_matrices[ip].integralMeasure * shape_matrices[ip].detJ;
214 
215  static const int kelvin_vector_size =
217  DisplacementDim);
218  ip_data.eps.setZero(kelvin_vector_size);
219  ip_data.eps_prev.resize(kelvin_vector_size);
220  ip_data.C_tensile.setZero(kelvin_vector_size, kelvin_vector_size);
221  ip_data.C_compressive.setZero(kelvin_vector_size,
222  kelvin_vector_size);
223  ip_data.sigma_tensile.setZero(kelvin_vector_size);
224  ip_data.sigma_compressive.setZero(kelvin_vector_size);
225  ip_data.sigma.setZero(kelvin_vector_size);
226  ip_data.strain_energy_tensile = 0.0;
227  ip_data.elastic_energy = 0.0;
228 
229  ip_data.N = shape_matrices[ip].N;
230  ip_data.dNdx = shape_matrices[ip].dNdx;
231 
232  _secondary_data.N[ip] = shape_matrices[ip].N;
233  }
234  }
#define OGS_FATAL(...)
Definition: Error.h:26
double getWeight() const
Definition: WeightedPoint.h:80
virtual std::size_t getID() const final
Returns the ID of the element.
Definition: Element.h:89
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
void setElementID(std::size_t element_id)
PhaseFieldProcessData< DisplacementDim > & _process_data
SecondaryData< typename ShapeMatrices::ShapeType > _secondary_data
NumLib::GenericIntegrationMethod const & _integration_method
ShapeMatrixPolicyType< ShapeFunction, DisplacementDim > ShapeMatricesType
std::vector< IpData, Eigen::aligned_allocator< IpData > > _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:24
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)
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N
Definition: SecondaryData.h:27

References ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_element, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_integration_method, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_process_data, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, MeshLib::Element::getID(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::initShapeMatrices(), MathLib::KelvinVector::kelvin_vector_dimensions(), ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N, OGS_FATAL, MaterialLib::Solids::selectSolidConstitutiveRelation(), and ParameterLib::SpatialPosition::setElementID().

Member Function Documentation

◆ assemble()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< 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 236 of file PhaseFieldFEM.h.

242  {
243  OGS_FATAL(
244  "PhaseFieldLocalAssembler: assembly without jacobian is not "
245  "implemented.");
246  }

References OGS_FATAL.

◆ assembleWithJacobianForDeformationEquations()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianForDeformationEquations ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
std::vector< double > &  local_b_data,
std::vector< double > &  local_Jac_data 
)
private

Definition at line 43 of file PhaseFieldFEM-impl.h.

47 {
48  using DeformationMatrix =
49  typename ShapeMatricesType::template MatrixType<displacement_size,
51  auto const d = local_x.template segment<phasefield_size>(phasefield_index);
52  auto const u =
53  local_x.template segment<displacement_size>(displacement_index);
54 
55  auto local_Jac = MathLib::createZeroedMatrix<DeformationMatrix>(
56  local_Jac_data, displacement_size, displacement_size);
57 
58  auto local_rhs = MathLib::createZeroedVector<DeformationVector>(
59  local_b_data, displacement_size);
60 
62  x_position.setElementID(_element.getID());
63 
64  auto local_pressure = 0.0;
65  if (_process_data.crack_pressure)
66  {
67  local_pressure = _process_data.unity_pressure;
68  }
69 
70  int const n_integration_points = _integration_method.getNumberOfPoints();
71  for (int ip = 0; ip < n_integration_points; ip++)
72  {
73  x_position.setIntegrationPoint(ip);
74  auto const& w = _ip_data[ip].integration_weight;
75  auto const& N = _ip_data[ip].N;
76  auto const& dNdx = _ip_data[ip].dNdx;
77 
78  auto const x_coord =
79  NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
80  _element, N);
81 
82  auto const& B =
83  LinearBMatrix::computeBMatrix<DisplacementDim,
84  ShapeFunction::NPOINTS,
86  dNdx, N, x_coord, _is_axially_symmetric);
87 
88  auto& eps = _ip_data[ip].eps;
89  eps.noalias() = B * u;
90  double const k = _process_data.residual_stiffness(t, x_position)[0];
91  double const ls = _process_data.crack_length_scale(t, x_position)[0];
92  double const d_ip = N.dot(d);
93  double const degradation =
94  _process_data.degradation_derivative->degradation(d_ip, k, ls);
95  _ip_data[ip].updateConstitutiveRelation(
96  t, x_position, dt, u, degradation,
97  _process_data.energy_split_model);
98 
99  auto& sigma = _ip_data[ip].sigma;
100  auto const& C_tensile = _ip_data[ip].C_tensile;
101  auto const& C_compressive = _ip_data[ip].C_compressive;
102 
103  typename ShapeMatricesType::template MatrixType<DisplacementDim,
105  N_u = ShapeMatricesType::template MatrixType<
106  DisplacementDim, displacement_size>::Zero(DisplacementDim,
108 
109  for (int i = 0; i < DisplacementDim; ++i)
110  {
111  N_u.template block<1, displacement_size / DisplacementDim>(
112  i, i * displacement_size / DisplacementDim)
113  .noalias() = N;
114  }
115 
116  auto const rho_sr = _process_data.solid_density(t, x_position)[0];
117  auto const& b = _process_data.specific_body_force;
118 
119  local_rhs.noalias() -=
120  (B.transpose() * sigma - N_u.transpose() * rho_sr * b -
121  local_pressure * N_u.transpose() * dNdx * d) *
122  w;
123  local_Jac.noalias() +=
124  B.transpose() * (degradation * C_tensile + C_compressive) * B * w;
125  }
126 }
void setIntegrationPoint(unsigned integration_point)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
Definition: BMatrixPolicy.h:55
static const double u
static const double t
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(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MathLib::t, and MathLib::u.

◆ assembleWithJacobianForStaggeredScheme()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianForStaggeredScheme ( 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,
std::vector< double > &  local_Jac_data 
)
overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 21 of file PhaseFieldFEM-impl.h.

28 {
29  // For the equations with phase field.
30  if (process_id == phase_process_id)
31  {
32  assembleWithJacobianPhaseFieldEquations(t, dt, local_x, local_b_data,
33  local_Jac_data);
34  return;
35  }
36 
37  // For the equations with deformation
38  assembleWithJacobianForDeformationEquations(t, dt, local_x, local_b_data,
39  local_Jac_data);
40 }
void assembleWithJacobianForDeformationEquations(double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
void assembleWithJacobianPhaseFieldEquations(double const t, double const dt, Eigen::VectorXd const &local_x, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)

References MathLib::t.

◆ assembleWithJacobianPhaseFieldEquations()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::assembleWithJacobianPhaseFieldEquations ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
std::vector< double > &  local_b_data,
std::vector< double > &  local_Jac_data 
)
private

Definition at line 129 of file PhaseFieldFEM-impl.h.

134 {
135  auto const d = local_x.template segment<phasefield_size>(phasefield_index);
136  auto const u =
137  local_x.template segment<displacement_size>(displacement_index);
138 
139  auto local_Jac = MathLib::createZeroedMatrix<PhaseFieldMatrix>(
140  local_Jac_data, phasefield_size, phasefield_size);
141  auto local_rhs = MathLib::createZeroedVector<PhaseFieldVector>(
142  local_b_data, phasefield_size);
143 
145  x_position.setElementID(_element.getID());
146 
147  auto local_pressure = 0.0;
148  if (_process_data.crack_pressure)
149  {
150  local_pressure = _process_data.unity_pressure;
151  }
152  else if (_process_data.hydro_crack)
153  {
154  local_pressure = _process_data.pressure;
155  }
156 
157  double const k = _process_data.residual_stiffness(t, x_position)[0];
158  double const ls = _process_data.crack_length_scale(t, x_position)[0];
159  double const gc = _process_data.crack_resistance(t, x_position)[0];
160 
161  int const n_integration_points = _integration_method.getNumberOfPoints();
162  for (int ip = 0; ip < n_integration_points; ip++)
163  {
164  x_position.setIntegrationPoint(ip);
165  auto const& w = _ip_data[ip].integration_weight;
166  auto const& N = _ip_data[ip].N;
167  auto const& dNdx = _ip_data[ip].dNdx;
168 
169  double const d_ip = N.dot(d);
170  double const degradation =
171  _process_data.degradation_derivative->degradation(d_ip, k, ls);
172  double const degradation_df1 =
173  _process_data.degradation_derivative->degradation_df1(d_ip, ls);
174  double const degradation_df2 =
175  _process_data.degradation_derivative->degradation_df2(d_ip, ls);
176  auto const x_coord =
177  NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
178  _element, N);
179  auto const& B =
180  LinearBMatrix::computeBMatrix<DisplacementDim,
181  ShapeFunction::NPOINTS,
182  typename BMatricesType::BMatrixType>(
183  dNdx, N, x_coord, _is_axially_symmetric);
184 
185  auto& eps = _ip_data[ip].eps;
186  eps.noalias() = B * u;
187  _ip_data[ip].updateConstitutiveRelation(
188  t, x_position, dt, u, degradation,
189  _process_data.energy_split_model);
190 
191  auto const& strain_energy_tensile = _ip_data[ip].strain_energy_tensile;
192 
193  auto& ip_data = _ip_data[ip];
194  ip_data.strain_energy_tensile = strain_energy_tensile;
195 
196  typename ShapeMatricesType::template MatrixType<DisplacementDim,
198  N_u = ShapeMatricesType::template MatrixType<
199  DisplacementDim, displacement_size>::Zero(DisplacementDim,
201 
202  for (int i = 0; i < DisplacementDim; ++i)
203  {
204  N_u.template block<1, displacement_size / DisplacementDim>(
205  i, i * displacement_size / DisplacementDim)
206  .noalias() = N;
207  }
208 
209  local_Jac.noalias() +=
210  (N.transpose() * N * degradation_df2 * strain_energy_tensile) * w;
211 
212  local_rhs.noalias() -=
213  (N.transpose() * degradation_df1 * strain_energy_tensile -
214  local_pressure * dNdx.transpose() * N_u * u) *
215  w;
216 
217  switch (_process_data.phasefield_model)
218  {
220  {
221  auto const local_Jac_AT1 =
222  (gc * 0.75 * dNdx.transpose() * dNdx * ls * w).eval();
223  local_Jac.noalias() += local_Jac_AT1;
224 
225  local_rhs.noalias() -=
226  gc * (-0.375 * N.transpose() / ls) * w + local_Jac_AT1 * d;
227  break;
228  }
230  {
231  auto const local_Jac_AT2 =
232  (gc *
233  (N.transpose() * N / ls + dNdx.transpose() * dNdx * ls) *
234  w)
235  .eval();
236  local_Jac.noalias() += local_Jac_AT2;
237 
238  local_rhs.noalias() -=
239  local_Jac_AT2 * d - gc * (N.transpose() / ls) * w;
240  break;
241  }
243  {
244  auto const local_Jac_COHESIVE =
245  (2.0 / boost::math::double_constants::pi * gc *
246  (-N.transpose() * N / ls + dNdx.transpose() * dNdx * ls) *
247  w)
248  .eval();
249 
250  local_Jac.noalias() += local_Jac_COHESIVE;
251 
252  local_rhs.noalias() -= local_Jac_COHESIVE * d;
253  break;
254  }
255  }
256  }
257 }

References ProcessLib::PhaseField::AT1, ProcessLib::PhaseField::AT2, ProcessLib::PhaseField::COHESIVE, ProcessLib::LinearBMatrix::computeBMatrix(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MathLib::t, and MathLib::u.

◆ computeCrackIntegral()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::computeCrackIntegral ( std::size_t  mesh_item_id,
std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &  dof_tables,
GlobalVector const &  x,
double const  t,
double &  crack_volume,
CoupledSolutionsForStaggeredScheme const *const  cpl_xs 
)
overridevirtual

Implements ProcessLib::PhaseField::PhaseFieldLocalAssemblerInterface.

Definition at line 260 of file PhaseFieldFEM-impl.h.

267 {
268  assert(cpl_xs != nullptr);
269 
270  std::vector<std::vector<GlobalIndexType>> indices_of_processes;
271  indices_of_processes.reserve(dof_tables.size());
272  std::transform(dof_tables.begin(), dof_tables.end(),
273  std::back_inserter(indices_of_processes),
274  [&](NumLib::LocalToGlobalIndexMap const& dof_table)
275  { return NumLib::getIndices(mesh_item_id, dof_table); });
276 
277  auto local_coupled_xs =
278  getCoupledLocalSolutions(cpl_xs->coupled_xs, indices_of_processes);
279  assert(local_coupled_xs.size() == displacement_size + phasefield_size);
280 
281  auto const d = Eigen::Map<PhaseFieldVector const>(
282  &local_coupled_xs[phasefield_index], phasefield_size);
283  auto const u = Eigen::Map<DeformationVector const>(
284  &local_coupled_xs[displacement_index], displacement_size);
285 
287  x_position.setElementID(_element.getID());
288 
289  int const n_integration_points = _integration_method.getNumberOfPoints();
290  for (int ip = 0; ip < n_integration_points; ip++)
291  {
292  x_position.setIntegrationPoint(ip);
293  auto const& w = _ip_data[ip].integration_weight;
294  auto const& N = _ip_data[ip].N;
295  auto const& dNdx = _ip_data[ip].dNdx;
296 
297  typename ShapeMatricesType::template MatrixType<DisplacementDim,
299  N_u = ShapeMatricesType::template MatrixType<
300  DisplacementDim, displacement_size>::Zero(DisplacementDim,
302 
303  for (int i = 0; i < DisplacementDim; ++i)
304  {
305  N_u.template block<1, displacement_size / DisplacementDim>(
306  i, i * displacement_size / DisplacementDim)
307  .noalias() = N;
308  }
309 
310  crack_volume += (N_u * u).dot(dNdx * d) * w;
311  }
312 }
std::vector< double > getCoupledLocalSolutions(std::vector< GlobalVector * > const &global_solutions, std::vector< std::vector< GlobalIndexType >> const &indices)

References ProcessLib::CoupledSolutionsForStaggeredScheme::coupled_xs, ProcessLib::getCoupledLocalSolutions(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), and MathLib::u.

◆ computeEnergy()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::computeEnergy ( std::size_t  mesh_item_id,
std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &  dof_tables,
GlobalVector const &  x,
double const  t,
double &  elastic_energy,
double &  surface_energy,
double &  pressure_work,
CoupledSolutionsForStaggeredScheme const *const  cpl_xs 
)
overridevirtual

Implements ProcessLib::PhaseField::PhaseFieldLocalAssemblerInterface.

Definition at line 315 of file PhaseFieldFEM-impl.h.

322 {
323  assert(cpl_xs != nullptr);
324 
325  std::vector<std::vector<GlobalIndexType>> indices_of_processes;
326  indices_of_processes.reserve(dof_tables.size());
327  std::transform(dof_tables.begin(), dof_tables.end(),
328  std::back_inserter(indices_of_processes),
329  [&](NumLib::LocalToGlobalIndexMap const& dof_table)
330  { return NumLib::getIndices(mesh_item_id, dof_table); });
331 
332  auto const local_coupled_xs =
333  getCoupledLocalSolutions(cpl_xs->coupled_xs, indices_of_processes);
334  assert(local_coupled_xs.size() == displacement_size + phasefield_size);
335 
336  auto const d = Eigen::Map<PhaseFieldVector const>(
337  &local_coupled_xs[phasefield_index], phasefield_size);
338  auto const u = Eigen::Map<DeformationVector const>(
339  &local_coupled_xs[displacement_index], displacement_size);
340 
342  x_position.setElementID(_element.getID());
343 
344  double element_elastic_energy = 0.0;
345  double element_surface_energy = 0.0;
346  double element_pressure_work = 0.0;
347 
348  double const gc = _process_data.crack_resistance(t, x_position)[0];
349  double const ls = _process_data.crack_length_scale(t, x_position)[0];
350  int const n_integration_points = _integration_method.getNumberOfPoints();
351  for (int ip = 0; ip < n_integration_points; ip++)
352  {
353  x_position.setIntegrationPoint(ip);
354  auto const& w = _ip_data[ip].integration_weight;
355  auto const& N = _ip_data[ip].N;
356  auto const& dNdx = _ip_data[ip].dNdx;
357  auto pressure_ip = _process_data.pressure;
358 
359  typename ShapeMatricesType::template MatrixType<DisplacementDim,
361  N_u = ShapeMatricesType::template MatrixType<
362  DisplacementDim, displacement_size>::Zero(DisplacementDim,
364 
365  for (int i = 0; i < DisplacementDim; ++i)
366  {
367  N_u.template block<1, displacement_size / DisplacementDim>(
368  i, i * displacement_size / DisplacementDim)
369  .noalias() = N;
370  }
371 
372  element_elastic_energy += _ip_data[ip].elastic_energy * w;
373 
374  double const d_ip = N.dot(d);
375  switch (_process_data.phasefield_model)
376  {
378  {
379  element_surface_energy +=
380  gc * 0.375 *
381  ((1 - d_ip) / ls + (dNdx * d).dot((dNdx * d)) * ls) * w;
382 
383  break;
384  }
386  {
387  element_surface_energy += 0.5 * gc *
388  ((1 - d_ip) * (1 - d_ip) / ls +
389  (dNdx * d).dot((dNdx * d)) * ls) *
390  w;
391  break;
392  }
394  {
395  element_surface_energy +=
396  gc / boost::math::double_constants::pi *
397  ((1 - d_ip * d_ip) / ls + (dNdx * d).dot((dNdx * d)) * ls) *
398  w;
399  break;
400  }
401  }
402 
403  if (_process_data.crack_pressure)
404  {
405  element_pressure_work += pressure_ip * (N_u * u).dot(dNdx * d) * w;
406  }
407  }
408 
409 #ifdef USE_PETSC
410  int const n_all_nodes = indices_of_processes[1].size();
411  int const n_regular_nodes = std::count_if(
412  begin(indices_of_processes[1]), end(indices_of_processes[1]),
413  [](GlobalIndexType const& index) { return index >= 0; });
414  if (n_all_nodes != n_regular_nodes)
415  {
416  element_elastic_energy *=
417  static_cast<double>(n_regular_nodes) / n_all_nodes;
418  element_surface_energy *=
419  static_cast<double>(n_regular_nodes) / n_all_nodes;
420  element_pressure_work *=
421  static_cast<double>(n_regular_nodes) / n_all_nodes;
422  }
423 #endif // USE_PETSC
424  elastic_energy += element_elastic_energy;
425  surface_energy += element_surface_energy;
426  pressure_work += element_pressure_work;
427 }
GlobalMatrix::IndexType GlobalIndexType

References ProcessLib::PhaseField::AT1, ProcessLib::PhaseField::AT2, ProcessLib::PhaseField::COHESIVE, ProcessLib::CoupledSolutionsForStaggeredScheme::coupled_xs, ProcessLib::getCoupledLocalSolutions(), ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MathLib::t, and MathLib::u.

◆ getIntPtEpsilon()

template<typename ShapeFunction , int DisplacementDim>
std::vector<double> const& ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::getIntPtEpsilon ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &  cache 
) const
inlineoverrideprivatevirtual

◆ getIntPtSigma()

template<typename ShapeFunction , int DisplacementDim>
std::vector<double> const& ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::getIntPtSigma ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &  cache 
) const
inlineoverrideprivatevirtual

◆ getShapeMatrix()

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

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 296 of file PhaseFieldFEM.h.

298  {
299  auto const& N = _secondary_data.N[integration_point];
300 
301  // assumes N is stored contiguously in memory
302  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
303  }

References ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data, and ProcessLib::HeatTransportBHE::SecondaryData< ShapeMatrixType >::N.

◆ initializeConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::initializeConcrete ( )
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 255 of file PhaseFieldFEM.h.

256  {
257  unsigned const n_integration_points =
259 
260  for (unsigned ip = 0; ip < n_integration_points; ip++)
261  {
262  _ip_data[ip].pushBackState();
263  }
264  }

References ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_integration_method, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, and NumLib::GenericIntegrationMethod::getNumberOfPoints().

◆ postTimestepConcrete()

template<typename ShapeFunction , int DisplacementDim>
void ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::postTimestepConcrete ( Eigen::VectorXd const &  ,
double const  ,
double const   
)
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 266 of file PhaseFieldFEM.h.

269  {
270  unsigned const n_integration_points =
272 
273  for (unsigned ip = 0; ip < n_integration_points; ip++)
274  {
275  _ip_data[ip].pushBackState();
276  }
277  }

References ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_integration_method, ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_ip_data, and NumLib::GenericIntegrationMethod::getNumberOfPoints().

Member Data Documentation

◆ _element

template<typename ShapeFunction , int DisplacementDim>
MeshLib::Element const& ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_element
private

◆ _integration_method

◆ _ip_data

◆ _is_axially_symmetric

template<typename ShapeFunction , int DisplacementDim>
bool const ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_is_axially_symmetric
private

Definition at line 341 of file PhaseFieldFEM.h.

◆ _process_data

template<typename ShapeFunction , int DisplacementDim>
PhaseFieldProcessData<DisplacementDim>& ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_process_data
private

◆ _secondary_data

template<typename ShapeFunction , int DisplacementDim>
SecondaryData<typename ShapeMatrices::ShapeType> ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::_secondary_data
private

◆ displacement_index

template<typename ShapeFunction , int DisplacementDim>
constexpr int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::displacement_index = 0
staticconstexprprivate

Definition at line 137 of file PhaseFieldFEM.h.

◆ displacement_size

template<typename ShapeFunction , int DisplacementDim>
constexpr int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::displacement_size
staticconstexprprivate
Initial value:
=
ShapeFunction::NPOINTS * DisplacementDim

Definition at line 138 of file PhaseFieldFEM.h.

◆ mechanics_process_id

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::mechanics_process_id = 0
staticprivate

Definition at line 344 of file PhaseFieldFEM.h.

◆ phase_process_id

template<typename ShapeFunction , int DisplacementDim>
const int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::phase_process_id = 1
staticprivate

Definition at line 343 of file PhaseFieldFEM.h.

◆ phasefield_index

template<typename ShapeFunction , int DisplacementDim>
constexpr int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::phasefield_index
staticconstexprprivate
Initial value:

Definition at line 140 of file PhaseFieldFEM.h.

◆ phasefield_size

template<typename ShapeFunction , int DisplacementDim>
constexpr int ProcessLib::PhaseField::PhaseFieldLocalAssembler< ShapeFunction, DisplacementDim >::phasefield_size = ShapeFunction::NPOINTS
staticconstexprprivate

Definition at line 142 of file PhaseFieldFEM.h.


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