template<typename ShapeFunction, int GlobalDim>
class ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >
Definition at line 49 of file HeatConductionFEM.h.
|
| LocalAssemblerData (MeshLib::Element const &element, std::size_t const local_matrix_size, NumLib::GenericIntegrationMethod const &integration_method, bool is_axially_symmetric, HeatConductionProcessData const &process_data) |
|
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 |
|
void | assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_x_prev, std::vector< double > &local_rhs_data, std::vector< double > &local_Jac_data) override |
|
Eigen::Map< const Eigen::RowVectorXd > | getShapeMatrix (const unsigned integration_point) const override |
| Provides the shape matrix at the given integration point.
|
|
std::vector< double > const & | getIntPtHeatFlux (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override |
|
virtual | ~LocalAssemblerInterface ()=default |
|
virtual void | setInitialConditions (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, 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_x_prev, 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_x_prev, int const process_id, 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_prev, 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, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, int const process_id) |
|
void | postNonLinearSolver (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &x_prev, double const t, double const dt, 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.
|
|
virtual std::optional< VectorSegment > | getVectorDeformationSegment () const |
|
template<typename ShapeFunction , int GlobalDim>
void ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >::assembleWithJacobian |
( |
double const | t, |
|
|
double const | dt, |
|
|
std::vector< double > const & | local_x, |
|
|
std::vector< double > const & | local_x_prev, |
|
|
std::vector< double > & | local_rhs_data, |
|
|
std::vector< double > & | local_Jac_data ) |
|
inlineoverridevirtual |
Reimplemented from ProcessLib::LocalAssemblerInterface.
Definition at line 152 of file HeatConductionFEM.h.
157 {
158 auto const local_matrix_size = local_x.size();
159
160
161 assert(local_matrix_size == ShapeFunction::NPOINTS *
NUM_NODAL_DOF);
162
163 auto x = Eigen::Map<NodalVectorType const>(local_x.data(),
164 local_matrix_size);
165
166 auto x_prev = Eigen::Map<NodalVectorType const>(local_x_prev.data(),
167 local_matrix_size);
168
170 local_Jac_data, local_matrix_size, local_matrix_size);
172 local_rhs_data, local_matrix_size);
173
175 NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
177 NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
178
179 unsigned const n_integration_points =
181
182 auto const& medium =
185
186 for (unsigned ip = 0; ip < n_integration_points; ip++)
187 {
194 sm.N))};
195
196 double const w =
198 sm.integralMeasure;
199
200
201
202 double T_int_pt = 0.0;
205
207 medium
208 .property(
210 .value(vars, pos, t, dt));
212 medium
215 .template value<double>(vars, pos, t, dt);
218 .template value<double>(vars, pos, t, dt);
219
220 laplace.noalias() += sm.dNdx.transpose() * k * sm.dNdx * w;
223 }
225 {
227 }
228
229 local_Jac.noalias() += laplace +
storage / dt;
230 local_rhs.noalias() -= laplace * x +
storage * (x - x_prev) / dt;
231 }
constexpr double getWeight() const
typename LocalAssemblerTraits::LocalMatrix NodalMatrixType
Eigen::Map< Vector > createZeroedVector(std::vector< double > &data, Eigen::VectorXd::Index size)
References ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >::_element, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >::_integration_method, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >::_process_data, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, GlobalDim >::_shape_matrices, MathLib::createZeroedMatrix(), MathLib::createZeroedVector(), MaterialPropertyLib::density, MaterialPropertyLib::formEigenTensor(), MeshLib::Element::getID(), MaterialPropertyLib::MaterialSpatialDistributionMap::getMedium(), NumLib::GenericIntegrationMethod::getNumberOfPoints(), MathLib::WeightedPoint::getWeight(), NumLib::GenericIntegrationMethod::getWeightedPoint(), NumLib::interpolateCoordinates(), ProcessLib::HeatConduction::HeatConductionProcessData::mass_lumping, ProcessLib::HeatConduction::HeatConductionProcessData::media_map, ProcessLib::HeatConduction::NUM_NODAL_DOF, NumLib::detail::shapeFunctionInterpolate(), MaterialPropertyLib::VariableArray::temperature, and MaterialPropertyLib::thermal_conductivity.