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 59 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, 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, 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 69 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 64 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 67 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 68 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 62 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 61 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 74 of file HeatConductionFEM.h.

79  : _element(element),
80  _process_data(process_data),
81  _integration_method(integration_order),
84  GlobalDim>(
85  element, is_axially_symmetric, _integration_method)),
87  GlobalDim,
88  std::vector<double>(_integration_method.getNumberOfPoints()))
89  {
90  // This assertion is valid only if all nodal d.o.f. use the same shape
91  // matrices.
92  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
93  (void)local_matrix_size;
94  }
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 96 of file HeatConductionFEM.h.

102  {
103  auto const local_matrix_size = local_x.size();
104  // This assertion is valid only if all nodal d.o.f. use the same shape
105  // matrices.
106  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
107 
108  auto local_M = MathLib::createZeroedMatrix<NodalMatrixType>(
109  local_M_data, local_matrix_size, local_matrix_size);
110  auto local_K = MathLib::createZeroedMatrix<NodalMatrixType>(
111  local_K_data, local_matrix_size, local_matrix_size);
112 
113  unsigned const n_integration_points =
114  _integration_method.getNumberOfPoints();
115 
117  pos.setElementID(_element.getID());
118 
119  auto const& medium =
120  *_process_data.media_map->getMedium(_element.getID());
122 
123  for (unsigned ip = 0; ip < n_integration_points; ip++)
124  {
125  pos.setIntegrationPoint(ip);
126  auto const& sm = _shape_matrices[ip];
127  auto const& wp = _integration_method.getWeightedPoint(ip);
128 
129  // get the local temperature and put it in the variable array for
130  // access in MPL
131  double T_int_pt = 0.0;
132  NumLib::shapeFunctionInterpolate(local_x, sm.N, T_int_pt);
133  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
134  T_int_pt;
135 
136  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
137  medium
138  .property(
140  .value(vars, pos, t, dt));
141  auto const specific_heat_capacity =
142  medium
145  .template value<double>(vars, pos, t, dt);
146  auto const density =
147  medium
148  .property(
150  .template value<double>(vars, pos, t, dt);
151 
152  local_K.noalias() += sm.dNdx.transpose() * k * sm.dNdx * sm.detJ *
153  wp.getWeight() * sm.integralMeasure;
154  local_M.noalias() += sm.N.transpose() * density *
155  specific_heat_capacity * sm.N * sm.detJ *
156  wp.getWeight() * sm.integralMeasure;
157  }
159  {
160  local_M = local_M.colwise().sum().eval().asDiagonal();
161  }
162  }
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
static const double t
void shapeFunctionInterpolate(const NodalValues &nodal_values, const ShapeMatrix &shape_matrix_N, double &interpolated_value, ScalarTypes &... interpolated_values)
Definition: Interpolation.h:79
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(), ProcessLib::HeatConduction::HeatConductionProcessData::mass_lumping, ProcessLib::HeatConduction::HeatConductionProcessData::media_map, ProcessLib::HeatConduction::NUM_NODAL_DOF, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), MaterialPropertyLib::specific_heat_capacity, MathLib::t, 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,
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 164 of file HeatConductionFEM.h.

171  {
172  auto const local_matrix_size = local_x.size();
173  // This assertion is valid only if all nodal d.o.f. use the same shape
174  // matrices.
175  assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
176 
177  auto x = Eigen::Map<NodalVectorType const>(local_x.data(),
178  local_matrix_size);
179 
180  auto x_dot = Eigen::Map<NodalVectorType const>(local_xdot.data(),
181  local_matrix_size);
182 
183  auto local_Jac = MathLib::createZeroedMatrix<NodalMatrixType>(
184  local_Jac_data, local_matrix_size, local_matrix_size);
185  auto local_rhs = MathLib::createZeroedVector<NodalVectorType>(
186  local_rhs_data, local_matrix_size);
187 
188  NodalMatrixType laplace =
189  NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
191  NodalMatrixType::Zero(local_matrix_size, local_matrix_size);
192 
193  unsigned const n_integration_points =
194  _integration_method.getNumberOfPoints();
195 
197  pos.setElementID(_element.getID());
198 
199  auto const& medium =
200  *_process_data.media_map->getMedium(_element.getID());
202 
203  for (unsigned ip = 0; ip < n_integration_points; ip++)
204  {
205  pos.setIntegrationPoint(ip);
206  auto const& sm = _shape_matrices[ip];
207  double const w =
208  _integration_method.getWeightedPoint(ip).getWeight() * sm.detJ *
209  sm.integralMeasure;
210 
211  // get the local temperature and put it in the variable array for
212  // access in MPL
213  double T_int_pt = 0.0;
214  NumLib::shapeFunctionInterpolate(local_x, sm.N, T_int_pt);
215  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
216  T_int_pt;
217 
218  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
219  medium
220  .property(
222  .value(vars, pos, t, dt));
223  auto const specific_heat_capacity =
224  medium
227  .template value<double>(vars, pos, t, dt);
228  auto const density =
229  medium
230  .property(
232  .template value<double>(vars, pos, t, dt);
233 
234  laplace.noalias() += sm.dNdx.transpose() * k * sm.dNdx * w;
235  storage.noalias() +=
236  sm.N.transpose() * density * specific_heat_capacity * sm.N * w;
237  }
239  {
240  storage = storage.colwise().sum().eval().asDiagonal();
241  }
242 
243  local_Jac.noalias() += laplace + storage / dt;
244  local_rhs.noalias() -= laplace * x + storage * x_dot;
245  }
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(), ProcessLib::HeatConduction::HeatConductionProcessData::mass_lumping, ProcessLib::HeatConduction::HeatConductionProcessData::media_map, ProcessLib::HeatConduction::NUM_NODAL_DOF, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), MaterialPropertyLib::specific_heat_capacity, MaterialPropertyLib::storage, MathLib::t, 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 247 of file HeatConductionFEM.h.

250  {
251  unsigned const n_integration_points =
252  _integration_method.getNumberOfPoints();
253 
255  pos.setElementID(_element.getID());
256 
257  auto const& medium =
258  *_process_data.media_map->getMedium(_element.getID());
260 
261  for (unsigned ip = 0; ip < n_integration_points; ip++)
262  {
263  pos.setIntegrationPoint(ip);
264  auto const& sm = _shape_matrices[ip];
265  // get the local temperature and put it in the variable array for
266  // access in MPL
267  double T_int_pt = 0.0;
268  NumLib::shapeFunctionInterpolate(local_x, sm.N, T_int_pt);
269  vars[static_cast<int>(MaterialPropertyLib::Variable::temperature)] =
270  T_int_pt;
271 
272  auto const k = MaterialPropertyLib::formEigenTensor<GlobalDim>(
273  medium
274  .property(
276  .value(vars, pos, t, dt));
277  // heat flux only computed for output.
278  GlobalDimVectorType const heat_flux = -k * sm.dNdx * local_x;
279 
280  for (int d = 0; d < GlobalDim; ++d)
281  {
282  _heat_fluxes[d][ip] = heat_flux[d];
283  }
284  }
285  }
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, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), MathLib::t, 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 287 of file HeatConductionFEM.h.

289  {
290  auto const& N = _shape_matrices[integration_point].N;
291 
292  // assumes N is stored contiguously in memory
293  return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
294  }

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: