29 namespace HeatConduction
40 std::vector<GlobalVector*>
const& x,
41 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
42 std::vector<double>& cache)
const = 0;
46 std::vector<GlobalVector*>
const& x,
47 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
48 std::vector<double>& cache)
const = 0;
52 std::vector<GlobalVector*>
const& x,
53 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& dof_table,
54 std::vector<double>& cache)
const = 0;
57 template <
typename ShapeFunction,
typename IntegrationMethod,
int GlobalDim>
74 std::size_t
const local_matrix_size,
75 bool is_axially_symmetric,
76 unsigned const integration_order,
91 assert(local_matrix_size == ShapeFunction::NPOINTS *
NUM_NODAL_DOF);
92 (void)local_matrix_size;
96 std::vector<double>
const& local_x,
97 std::vector<double>
const& ,
98 std::vector<double>& local_M_data,
99 std::vector<double>& local_K_data,
100 std::vector<double>& )
override
102 auto const local_matrix_size = local_x.size();
105 assert(local_matrix_size == ShapeFunction::NPOINTS *
NUM_NODAL_DOF);
107 auto local_M = MathLib::createZeroedMatrix<NodalMatrixType>(
108 local_M_data, local_matrix_size, local_matrix_size);
109 auto local_K = MathLib::createZeroedMatrix<NodalMatrixType>(
110 local_K_data, local_matrix_size, local_matrix_size);
112 unsigned const n_integration_points =
125 .template value<double>(vars, pos, t, dt);
127 for (
unsigned ip = 0; ip < n_integration_points; ip++)
132 auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
136 .value(vars, pos, t, dt));
141 .template value<double>(vars, pos, t, dt);
146 .template value<double>(vars, pos, t, dt);
148 local_K.noalias() += sm.dNdx.transpose() * k * sm.dNdx * sm.detJ *
149 wp.getWeight() * sm.integralMeasure;
150 local_M.noalias() += sm.N.transpose() * density *
heat_capacity *
151 sm.N * sm.detJ * wp.getWeight() *
156 local_M = local_M.colwise().sum().eval().asDiagonal();
161 std::vector<double>
const& local_x,
162 std::vector<double>
const& local_xdot,
163 const double ,
const double ,
164 std::vector<double>& ,
165 std::vector<double>& ,
166 std::vector<double>& local_rhs_data,
167 std::vector<double>& local_Jac_data)
override
169 auto const local_matrix_size = local_x.size();
172 assert(local_matrix_size == ShapeFunction::NPOINTS *
NUM_NODAL_DOF);
174 auto x = Eigen::Map<NodalVectorType const>(local_x.data(),
177 auto x_dot = Eigen::Map<NodalVectorType const>(local_xdot.data(),
180 auto local_Jac = MathLib::createZeroedMatrix<NodalMatrixType>(
181 local_Jac_data, local_matrix_size, local_matrix_size);
182 auto local_rhs = MathLib::createZeroedVector<NodalVectorType>(
183 local_rhs_data, local_matrix_size);
186 NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
188 NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
190 unsigned const n_integration_points =
203 .template value<double>(vars, pos, t, dt);
205 for (
unsigned ip = 0; ip < n_integration_points; ip++)
213 auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
217 .value(vars, pos, t, dt));
222 .template value<double>(vars, pos, t, dt);
227 .template value<double>(vars, pos, t, dt);
229 laplace.noalias() += sm.dNdx.transpose() * k * sm.dNdx * w;
238 local_Jac.noalias() += laplace +
storage / dt;
239 local_rhs.noalias() -= laplace * x +
storage * x_dot;
243 double const t,
double const dt, Eigen::VectorXd
const& local_x,
244 Eigen::VectorXd
const& )
override
246 unsigned const n_integration_points =
259 .template value<double>(vars, pos, t, dt);
261 for (
unsigned ip = 0; ip < n_integration_points; ip++)
265 auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
269 .value(vars, pos, t, dt));
273 for (
int d = 0; d < GlobalDim; ++d)
281 const unsigned integration_point)
const override
286 return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
291 std::vector<GlobalVector*>
const& ,
292 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
293 std::vector<double>& )
const override
301 std::vector<GlobalVector*>
const& ,
302 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
303 std::vector<double>& )
const override
311 std::vector<GlobalVector*>
const& ,
312 std::vector<NumLib::LocalToGlobalIndexMap const*>
const& ,
313 std::vector<double>& )
const override
324 std::vector<ShapeMatrices, Eigen::aligned_allocator<ShapeMatrices>>
virtual std::size_t getID() const final
Returns the ID of the element.
void setElementID(std::size_t element_id)
void setIntegrationPoint(unsigned integration_point)
virtual std::vector< double > const & getIntPtHeatFluxY(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtHeatFluxX(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
virtual std::vector< double > const & getIntPtHeatFluxZ(const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const =0
std::vector< std::vector< double > > _heat_fluxes
LocalAssemblerData(MeshLib::Element const &element, std::size_t const local_matrix_size, bool is_axially_symmetric, unsigned const integration_order, HeatConductionProcessData const &process_data)
std::vector< double > const & getIntPtHeatFluxY(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const override
typename LocalAssemblerTraits::LocalMatrix NodalMatrixType
void assemble(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &) override
HeatConductionProcessData const & _process_data
Eigen::Map< const Eigen::RowVectorXd > getShapeMatrix(const unsigned integration_point) const override
Provides the shape matrix at the given integration point.
MeshLib::Element const & _element
ShapeMatrixPolicyType< ShapeFunction, GlobalDim > ShapeMatricesType
std::vector< double > const & getIntPtHeatFluxZ(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const override
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType
typename ShapeMatricesType::ShapeMatrices ShapeMatrices
typename LocalAssemblerTraits::LocalVector NodalVectorType
void assembleWithJacobian(double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double, const double, std::vector< double > &, std::vector< double > &, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override
void computeSecondaryVariableConcrete(double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &) override
std::vector< ShapeMatrices, Eigen::aligned_allocator< ShapeMatrices > > _shape_matrices
IntegrationMethod const _integration_method
std::vector< double > const & getIntPtHeatFluxX(const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const override
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
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)
const unsigned NUM_NODAL_DOF
NumLib::ShapeMatrices< NodalRowVectorType, DimNodalMatrixType, DimMatrixType, GlobalDimNodalMatrixType > ShapeMatrices
VectorType< GlobalDim > GlobalDimVectorType
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map
bool const mass_lumping
If set mass lumping will be applied to the equation.
Vector< NNodes *NodalDOF > LocalVector
Matrix< NNodes *NodalDOF, NNodes *NodalDOF > LocalMatrix