OGS
ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, typename IntegrationMethod, int GlobalDim>
class ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >

Definition at line 58 of file HeatConductionFEM.h.

#include <HeatConductionFEM.h>

Inheritance diagram for ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >:
[legend]
Collaboration diagram for ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >:
[legend]

Public Member Functions

 LocalAssemblerData (MeshLib::Element const &element, std::size_t const local_matrix_size, bool is_axially_symmetric, unsigned const integration_order, 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_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
 
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 & getIntPtHeatFluxX (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const override
 
std::vector< double > const & getIntPtHeatFluxY (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) const override
 
std::vector< double > const & getIntPtHeatFluxZ (const double, std::vector< GlobalVector * > const &, std::vector< NumLib::LocalToGlobalIndexMap const * > const &, std::vector< double > &) 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 Types

using ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, GlobalDim >
 
using ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
 
using LocalAssemblerTraits = ProcessLib::LocalAssemblerTraits< ShapeMatricesType, ShapeFunction::NPOINTS, NUM_NODAL_DOF, GlobalDim >
 
using NodalMatrixType = typename LocalAssemblerTraits::LocalMatrix
 
using NodalVectorType = typename LocalAssemblerTraits::LocalVector
 
using GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
 

Private Attributes

MeshLib::Element const & _element
 
HeatConductionProcessData const & _process_data
 
IntegrationMethod const _integration_method
 
std::vector< ShapeMatrices, Eigen::aligned_allocator< ShapeMatrices > > _shape_matrices
 
std::vector< std::vector< double > > _heat_fluxes
 

Member Typedef Documentation

◆ GlobalDimVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
private

Definition at line 68 of file HeatConductionFEM.h.

◆ LocalAssemblerTraits

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalAssemblerTraits = ProcessLib::LocalAssemblerTraits< ShapeMatricesType, ShapeFunction::NPOINTS, NUM_NODAL_DOF, GlobalDim>
private

Definition at line 63 of file HeatConductionFEM.h.

◆ NodalMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::NodalMatrixType = typename LocalAssemblerTraits::LocalMatrix
private

Definition at line 66 of file HeatConductionFEM.h.

◆ NodalVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::NodalVectorType = typename LocalAssemblerTraits::LocalVector
private

Definition at line 67 of file HeatConductionFEM.h.

◆ ShapeMatrices

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
private

Definition at line 61 of file HeatConductionFEM.h.

◆ ShapeMatricesType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
using ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, GlobalDim>
private

Definition at line 60 of file HeatConductionFEM.h.

Constructor & Destructor Documentation

◆ LocalAssemblerData()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalAssemblerData ( MeshLib::Element const &  element,
std::size_t const  local_matrix_size,
bool  is_axially_symmetric,
unsigned const  integration_order,
HeatConductionProcessData const &  process_data 
)
inline

The thermal_conductivity factor is directly integrated into the local element matrix.

Definition at line 73 of file HeatConductionFEM.h.

78  : _element(element),
79  _process_data(process_data),
80  _integration_method(integration_order),
83  GlobalDim>(
84  element, is_axially_symmetric, _integration_method)),
86  GlobalDim,
87  std::vector<double>(_integration_method.getNumberOfPoints()))
88  {
89  // This assertion is valid only if all nodal d.o.f. use the same shape
90  // matrices.
91  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
92  (void)local_matrix_size;
93  }
std::vector< std::vector< double > > _heat_fluxes
HeatConductionProcessData const & _process_data
ShapeMatrixPolicyType< ShapeFunction, GlobalDim > ShapeMatricesType
std::vector< ShapeMatrices, Eigen::aligned_allocator< ShapeMatrices > > _shape_matrices
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)

References ProcessLib::HeatConduction::NUM_NODAL_DOF.

Member Function Documentation

◆ assemble()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::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 > &   
)
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 95 of file HeatConductionFEM.h.

101  {
102  auto const local_matrix_size = local_x.size();
103  // This assertion is valid only if all nodal d.o.f. use the same shape
104  // matrices.
105  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
106 
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);
111 
112  unsigned const n_integration_points =
113  _integration_method.getNumberOfPoints();
114 
116  pos.setElementID(_element.getID());
117 
118  auto const& medium =
119  *_process_data.media_map->getMedium(_element.getID());
121  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
122  medium
123  .property(
125  .template value<double>(vars, pos, t, dt);
126 
127  for (unsigned ip = 0; ip < n_integration_points; ip++)
128  {
129  pos.setIntegrationPoint(ip);
130  auto const& sm = _shape_matrices[ip];
131  auto const& wp = _integration_method.getWeightedPoint(ip);
132  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
133  medium
134  .property(
136  .value(vars, pos, t, dt));
137  auto const heat_capacity =
138  medium
139  .property(
141  .template value<double>(vars, pos, t, dt);
142  auto const density =
143  medium
144  .property(
146  .template value<double>(vars, pos, t, dt);
147 
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() *
152  sm.integralMeasure;
153  }
155  {
156  local_M = local_M.colwise().sum().eval().asDiagonal();
157  }
158  }
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)
std::array< VariableType, static_cast< int >(Variable::number_of_variables)> VariableArray
Definition: VariableType.h:108
std::unique_ptr< MaterialPropertyLib::MaterialSpatialDistributionMap > media_map
bool const mass_lumping
If set mass lumping will be applied to the equation.

References ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_shape_matrices, MaterialPropertyLib::density, MeshLib::Element::getID(), MaterialPropertyLib::heat_capacity, ProcessLib::HeatConduction::HeatConductionProcessData::mass_lumping, ProcessLib::HeatConduction::HeatConductionProcessData::media_map, ProcessLib::HeatConduction::NUM_NODAL_DOF, MaterialPropertyLib::reference_temperature, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::temperature, and MaterialPropertyLib::thermal_conductivity.

◆ assembleWithJacobian()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::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 
)
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 160 of file HeatConductionFEM.h.

168  {
169  auto const local_matrix_size = local_x.size();
170  // This assertion is valid only if all nodal d.o.f. use the same shape
171  // matrices.
172  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
173 
174  auto x = Eigen::Map<NodalVectorType const>(local_x.data(),
175  local_matrix_size);
176 
177  auto x_dot = Eigen::Map<NodalVectorType const>(local_xdot.data(),
178  local_matrix_size);
179 
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);
184 
185  NodalMatrixType laplace =
186  NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
188  NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
189 
190  unsigned const n_integration_points =
191  _integration_method.getNumberOfPoints();
192 
194  pos.setElementID(_element.getID());
195 
196  auto const& medium =
197  *_process_data.media_map->getMedium(_element.getID());
199  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
200  medium
201  .property(
203  .template value<double>(vars, pos, t, dt);
204 
205  for (unsigned ip = 0; ip < n_integration_points; ip++)
206  {
207  pos.setIntegrationPoint(ip);
208  auto const& sm = _shape_matrices[ip];
209  double const w =
210  _integration_method.getWeightedPoint(ip).getWeight() * sm.detJ *
211  sm.integralMeasure;
212 
213  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
214  medium
215  .property(
217  .value(vars, pos, t, dt));
218  auto const heat_capacity =
219  medium
220  .property(
222  .template value<double>(vars, pos, t, dt);
223  auto const density =
224  medium
225  .property(
227  .template value<double>(vars, pos, t, dt);
228 
229  laplace.noalias() += sm.dNdx.transpose() * k * sm.dNdx * w;
230  storage.noalias() +=
231  sm.N.transpose() * density * heat_capacity * sm.N * w;
232  }
234  {
235  storage = storage.colwise().sum().eval().asDiagonal();
236  }
237 
238  local_Jac.noalias() += laplace + storage / dt;
239  local_rhs.noalias() -= laplace * x + storage * x_dot;
240  }
typename LocalAssemblerTraits::LocalMatrix NodalMatrixType

References ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_shape_matrices, MaterialPropertyLib::density, MeshLib::Element::getID(), MaterialPropertyLib::heat_capacity, ProcessLib::HeatConduction::HeatConductionProcessData::mass_lumping, ProcessLib::HeatConduction::HeatConductionProcessData::media_map, ProcessLib::HeatConduction::NUM_NODAL_DOF, MaterialPropertyLib::reference_temperature, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::storage, MaterialPropertyLib::temperature, and MaterialPropertyLib::thermal_conductivity.

◆ computeSecondaryVariableConcrete()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
void ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::computeSecondaryVariableConcrete ( double const  t,
double const  dt,
Eigen::VectorXd const &  local_x,
Eigen::VectorXd const &   
)
inlineoverridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 242 of file HeatConductionFEM.h.

245  {
246  unsigned const n_integration_points =
247  _integration_method.getNumberOfPoints();
248 
250  pos.setElementID(_element.getID());
251 
252  auto const& medium =
253  *_process_data.media_map->getMedium(_element.getID());
255  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
256  medium
257  .property(
259  .template value<double>(vars, pos, t, dt);
260 
261  for (unsigned ip = 0; ip < n_integration_points; ip++)
262  {
263  pos.setIntegrationPoint(ip);
264  auto const& sm = _shape_matrices[ip];
265  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
266  medium
267  .property(
269  .value(vars, pos, t, dt));
270  // heat flux only computed for output.
271  GlobalDimVectorType const heat_flux = -k * sm.dNdx * local_x;
272 
273  for (int d = 0; d < GlobalDim; ++d)
274  {
275  _heat_fluxes[d][ip] = heat_flux[d];
276  }
277  }
278  }
typename ShapeMatricesType::GlobalDimVectorType GlobalDimVectorType

References ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_element, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_heat_fluxes, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data, ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_shape_matrices, MeshLib::Element::getID(), ProcessLib::HeatConduction::HeatConductionProcessData::media_map, MaterialPropertyLib::reference_temperature, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), MaterialPropertyLib::temperature, and MaterialPropertyLib::thermal_conductivity.

◆ getIntPtHeatFluxX()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
std::vector<double> const& ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtHeatFluxX ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &   
) const
inlineoverridevirtual

◆ getIntPtHeatFluxY()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
std::vector<double> const& ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtHeatFluxY ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &   
) const
inlineoverridevirtual

◆ getIntPtHeatFluxZ()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
std::vector<double> const& ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtHeatFluxZ ( const double  ,
std::vector< GlobalVector * > const &  ,
std::vector< NumLib::LocalToGlobalIndexMap const * > const &  ,
std::vector< double > &   
) const
inlineoverridevirtual

◆ getShapeMatrix()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>
Eigen::Map<const Eigen::RowVectorXd> ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getShapeMatrix ( const unsigned  integration_point) const
inlineoverridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 280 of file HeatConductionFEM.h.

282  {
283  auto const& N = _shape_matrices[integration_point].N;
284 
285  // assumes N is stored contiguously in memory
286  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
287  }

References ProcessLib::HeatConduction::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_shape_matrices.

Member Data Documentation

◆ _element

◆ _heat_fluxes

◆ _integration_method

◆ _process_data

◆ _shape_matrices


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