35namespace SmallDeformation
39template <
typename BMatricesType,
typename ShapeMatricesType,
44 typename ShapeMatricesType::NodalRowVectorType
N;
45 typename ShapeMatricesType::GlobalDimNodalMatrixType
dNdx;
52template <
typename ShapeMatrixType>
55 std::vector<ShapeMatrixType, Eigen::aligned_allocator<ShapeMatrixType>>
N;
58template <
typename ShapeFunction,
int DisplacementDim>
93 bool const is_axially_symmetric,
96 e, integration_method, is_axially_symmetric, process_data)
98 unsigned const n_integration_points =
101 _ip_data.resize(n_integration_points);
104 auto const shape_matrices =
109 for (
unsigned ip = 0; ip < n_integration_points; ip++)
112 auto const& sm = shape_matrices[ip];
115 sm.integralMeasure * sm.detJ;
118 ip_data.dNdx = sm.dNdx;
126 unsigned const n_integration_points =
128 for (
unsigned ip = 0; ip < n_integration_points; ip++)
146 double>::quiet_NaN() ,
153 t, x_position, *material_state.material_state_variables);
155 material_state.pushBackState();
158 for (
unsigned ip = 0; ip < n_integration_points; ip++)
166 Eigen::Ref<Eigen::VectorXd const>
const& u,
167 Eigen::Ref<Eigen::VectorXd const>
const& u_prev,
183 auto const& N = ip_data.
N;
184 auto const& dNdx = ip_data.
dNdx;
186 NumLib::interpolateXCoordinate<ShapeFunction, ShapeMatricesType>(
190 ShapeFunction::NPOINTS,
196 ? (*this->
process_data_.reference_temperature)(t, x_position)[0]
197 : std::numeric_limits<double>::quiet_NaN();
205 CS.
eval(models, t, dt, x_position,
207 T_ref, B * u, B * u_prev,
208 current_state, prev_state, material_state, tmp, output_data,
215 std::vector<double>
const& ,
216 std::vector<double>
const& ,
217 std::vector<double>& ,
218 std::vector<double>& ,
219 std::vector<double>& )
override
222 "SmallDeformationLocalAssembler: assembly without jacobian is not "
227 std::vector<double>
const& local_x,
228 std::vector<double>
const& local_x_prev,
229 std::vector<double>& ,
230 std::vector<double>& ,
231 std::vector<double>& local_b_data,
232 std::vector<double>& local_Jac_data)
override
234 auto const local_matrix_size = local_x.size();
236 auto local_Jac = MathLib::createZeroedMatrix<StiffnessMatrixType>(
237 local_Jac_data, local_matrix_size, local_matrix_size);
239 auto local_b = MathLib::createZeroedVector<NodalDisplacementVectorType>(
240 local_b_data, local_matrix_size);
243 auto [u_prev] =
localDOF(local_x_prev);
245 unsigned const n_integration_points =
256 for (
unsigned ip = 0; ip < n_integration_points; ip++)
259 auto const& w =
_ip_data[ip].integration_weight;
261 auto const& dNdx =
_ip_data[ip].dNdx;
268 DisplacementDim, ShapeFunction::NPOINTS,
273 u, u_prev, x_position, t, dt,
_ip_data[ip],
279 auto const& b = *CD.volumetric_body_force;
280 auto const& C = CD.s_mech_data.stiffness_tensor;
283 (B.transpose() * sigma -
N_u_op(N).transpose() * b) * w;
284 local_Jac.noalias() += B.transpose() * C * B * w;
289 Eigen::VectorXd
const& local_x_prev,
290 double const t,
double const dt,
293 unsigned const n_integration_points =
305 for (
unsigned ip = 0; ip < n_integration_points; ip++)
310 local_x, local_x_prev, x_position, t, dt,
_ip_data[ip],
318 for (
unsigned ip = 0; ip < n_integration_points; ip++)
325 std::vector<double>
const& local_x,
326 std::vector<double>& nodal_values)
override
339 const unsigned integration_point)
const override
344 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
348 static constexpr auto localDOF(std::vector<double>
const& x)
355 std::vector<IpData, Eigen::aligned_allocator<IpData>>
_ip_data;
std::size_t getID() const
Returns the ID of the element.
MathLib::WeightedPoint const & getWeightedPoint(unsigned const igp) const
unsigned getNumberOfPoints() const
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
MatrixType< _kelvin_vector_size, _number_of_dof > BMatrixType
MatrixType< _number_of_dof, _number_of_dof > StiffnessMatrixType
VectorType< _number_of_dof > NodalForceVectorType
Rhs residual.
VectorTypeFixedSize< MathLib::VectorizedTensor::size(DisplacementDim)> GradientVectorType
MatrixType< MathLib::VectorizedTensor::size(DisplacementDim), _number_of_dof > GradientMatrixType
Eigen::Matrix< double, Eigen::MatrixBase< Derived >::RowsAtCompileTime, 1 > symmetricTensorToKelvinVector(Eigen::MatrixBase< Derived > const &v)
constexpr Eigen::CwiseNullaryOp< EigenBlockMatrixViewFunctor< D, M >, typename EigenBlockMatrixViewFunctor< D, M >::Matrix > eigenBlockMatrixView(const Eigen::MatrixBase< M > &matrix)
double interpolateXCoordinate(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
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)
auto localDOF(ElementDOFVector const &x)
std::array< double, 3 > interpolateCoordinates(MeshLib::Element const &e, typename ShapeMatricesType::ShapeMatrices::ShapeType const &N)
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.
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
MatrixType< ShapeFunction::NPOINTS, ShapeFunction::NPOINTS > NodalMatrixType
VectorType< ShapeFunction::NPOINTS > NodalVectorType
RowVectorType< ShapeFunction::NPOINTS > NodalRowVectorType
std::vector< ShapeMatrixType, Eigen::aligned_allocator< ShapeMatrixType > > N