24 namespace ThermoMechanics
26 template <
int DisplacementDim>
29 std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
30 std::vector<std::unique_ptr<ParameterLib::ParameterBase>>
const& parameters,
31 unsigned const integration_order,
32 std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
36 bool const use_monolithic_scheme)
37 :
Process(std::move(
name), mesh, std::move(jacobian_assembler), parameters,
38 integration_order, std::move(process_variables),
39 std::move(secondary_variables), use_monolithic_scheme),
40 _process_data(std::move(process_data))
45 _heat_flux = MeshLib::getOrCreateMeshProperty<double>(
49 std::make_unique<IntegrationPointWriter>(
56 std::make_unique<IntegrationPointWriter>(
63 std::make_unique<IntegrationPointWriter>(
70 template <
int DisplacementDim>
76 template <
int DisplacementDim>
79 const int process_id)
const
83 if (_use_monolithic_scheme ||
84 process_id == _process_data.mechanics_process_id)
86 auto const& l = *_local_to_global_index_map;
87 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
88 &l.getGhostIndices(), &this->_sparsity_pattern};
92 auto const& l = *_local_to_global_index_map_single_component;
93 return {l.dofSizeWithoutGhosts(), l.dofSizeWithoutGhosts(),
94 &l.getGhostIndices(), &_sparsity_pattern_with_single_component};
99 template <
int DisplacementDim>
105 if (_use_monolithic_scheme)
107 constructMonolithicProcessDofTable();
110 constructDofTableOfSpecifiedProcessStaggeredScheme(
111 _process_data.mechanics_process_id);
115 std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
116 *_mesh_subset_all_nodes};
117 _local_to_global_index_map_single_component.reset(
119 std::move(all_mesh_subsets_single_component),
123 if (!_use_monolithic_scheme)
125 _sparsity_pattern_with_single_component =
127 *_local_to_global_index_map_single_component, _mesh);
131 template <
int DisplacementDim>
135 unsigned const integration_order)
142 auto add_secondary_variable = [&](std::string
const&
name,
143 int const num_components,
144 auto get_ip_values_function)
146 _secondary_variables.addSecondaryVariable(
150 std::move(get_ip_values_function)));
153 add_secondary_variable(
"sigma",
155 DisplacementDim>::RowsAtCompileTime,
156 &LocalAssemblerInterface::getIntPtSigma);
158 add_secondary_variable(
"epsilon",
160 DisplacementDim>::RowsAtCompileTime,
161 &LocalAssemblerInterface::getIntPtEpsilon);
168 add_secondary_variable);
171 for (
auto const& ip_writer : _integration_point_writer)
174 auto const&
name = ip_writer->name();
179 auto const& mesh_property =
183 if (mesh_property.getMeshItemType() !=
192 if (ip_meta_data.n_components !=
193 mesh_property.getNumberOfGlobalComponents())
196 "Different number of components in meta data ({:d}) than in "
197 "the integration point field data for '{:s}': {:d}.",
198 ip_meta_data.n_components,
name,
199 mesh_property.getNumberOfGlobalComponents());
204 std::size_t position = 0;
205 for (
auto& local_asm : _local_assemblers)
207 std::size_t
const integration_points_read =
208 local_asm->setIPDataInitialConditions(
209 name, &mesh_property[position],
210 ip_meta_data.integration_order);
211 if (integration_points_read == 0)
214 "No integration points read in the integration point "
215 "initial conditions set function.");
217 position += integration_points_read * ip_meta_data.n_components;
224 *_local_to_global_index_map);
227 template <
int DisplacementDim>
230 if (_use_monolithic_scheme)
232 const int process_id_of_thermomechanics = 0;
233 initializeProcessBoundaryConditionsAndSourceTerms(
234 *_local_to_global_index_map, process_id_of_thermomechanics);
240 initializeProcessBoundaryConditionsAndSourceTerms(
241 *_local_to_global_index_map_single_component,
242 _process_data.heat_conduction_process_id);
245 initializeProcessBoundaryConditionsAndSourceTerms(
246 *_local_to_global_index_map, _process_data.mechanics_process_id);
249 template <
int DisplacementDim>
251 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
252 std::vector<GlobalVector*>
const& xdot,
int const process_id,
255 DBUG(
"Assemble ThermoMechanicsProcess.");
257 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
258 dof_table = {std::ref(*_local_to_global_index_map)};
268 template <
int DisplacementDim>
271 const double t,
double const dt, std::vector<GlobalVector*>
const& x,
272 std::vector<GlobalVector*>
const& xdot,
const double dxdot_dx,
273 const double dx_dx,
int const process_id,
GlobalMatrix& M,
276 DBUG(
"AssembleJacobian ThermoMechanicsProcess.");
278 std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
281 if (_use_monolithic_scheme)
284 "Assemble the Jacobian of ThermoMechanics for the monolithic"
286 dof_tables.emplace_back(*_local_to_global_index_map);
291 if (process_id == _process_data.heat_conduction_process_id)
294 "Assemble the Jacobian equations of heat conduction process in "
295 "ThermoMechanics for the staggered scheme.");
300 "Assemble the Jacobian equations of mechanical process in "
301 "ThermoMechanics for the staggered scheme.");
305 if (_process_data.heat_conduction_process_id ==
308 dof_tables.emplace_back(
309 *_local_to_global_index_map_single_component);
310 dof_tables.emplace_back(*_local_to_global_index_map);
314 dof_tables.emplace_back(*_local_to_global_index_map);
315 dof_tables.emplace_back(
316 *_local_to_global_index_map_single_component);
325 dxdot_dx, dx_dx, process_id, M, K, b, Jac);
328 auto copyRhs = [&](
int const variable_id,
auto& output_vector)
330 if (_use_monolithic_scheme)
334 std::negate<double>());
340 std::negate<double>());
343 if (_use_monolithic_scheme ||
344 process_id == _process_data.heat_conduction_process_id)
346 copyRhs(0, *_heat_flux);
348 if (_use_monolithic_scheme ||
349 process_id == _process_data.mechanics_process_id)
351 copyRhs(1, *_nodal_forces);
355 template <
int DisplacementDim>
357 std::vector<GlobalVector*>
const& x,
double const t,
double const dt,
358 const int process_id)
360 DBUG(
"PreTimestep ThermoMechanicsProcess.");
364 assert(process_id < 2);
366 if (process_id == _process_data.mechanics_process_id)
371 *x[process_id], t, dt);
376 template <
int DisplacementDim>
378 std::vector<GlobalVector*>
const& x,
double const t,
double const dt,
379 int const process_id)
386 DBUG(
"PostTimestep ThermoMechanicsProcess.");
387 std::vector<NumLib::LocalToGlobalIndexMap const*> dof_tables;
388 auto const n_processes = x.size();
389 dof_tables.reserve(n_processes);
390 for (std::size_t process_id = 0; process_id < n_processes; ++process_id)
392 dof_tables.push_back(&getDOFTable(process_id));
401 template <
int DisplacementDim>
405 if (_process_data.mechanics_process_id == process_id)
407 return *_local_to_global_index_map;
411 return *_local_to_global_index_map_single_component;
void DBUG(char const *fmt, Args const &... args)
Global vector based on Eigen vector.
bool isAxiallySymmetric() const
unsigned getDimension() const
Returns the dimension of the mesh (determined by the maximum dimension over all elements).
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Properties & getProperties()
bool existsPropertyVector(std::string const &name) const
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)
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 initialize(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
std::vector< std::size_t > const & getActiveElementIDs() const
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
Handles configuration of several secondary variables from the project file.
Local assembler of ThermoMechanics process.
MeshLib::PropertyVector< double > * _nodal_forces
void initializeBoundaryConditions() override
void postTimestepConcreteProcess(std::vector< GlobalVector * > const &x, const double t, const double dt, int const process_id) override
void preTimestepConcreteProcess(std::vector< GlobalVector * > const &x, double const t, double const dt, const int process_id) override
ThermoMechanicsProcess(std::string name, MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase >> const ¶meters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable >>> &&process_variables, ThermoMechanicsProcessData< DisplacementDim > &&process_data, SecondaryVariableCollection &&secondary_variables, bool const use_monolithic_scheme)
void assembleWithJacobianConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, const double dxdot_dx, const double dx_dx, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
bool isLinear() const override
void constructDofTable() override
MathLib::MatrixSpecifications getMatrixSpecifications(const int process_id) const override
std::vector< std::unique_ptr< LocalAssemblerInterface > > _local_assemblers
MeshLib::PropertyVector< double > * _heat_flux
NumLib::LocalToGlobalIndexMap const & getDOFTable(const int process_id) const override
void assembleConcreteProcess(const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
void assembleWithJacobian(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, const double t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, const double dxdot_dx, const double dx_dx, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac)
void assemble(std::size_t const mesh_item_id, LocalAssemblerInterface &local_assembler, std::vector< std::reference_wrapper< NumLib::LocalToGlobalIndexMap >> const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, std::vector< GlobalVector * > const &xdot, int const process_id, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b)
Eigen::Matrix< double, kelvin_vector_dimensions(DisplacementDim), 1, Eigen::ColMajor > KelvinVectorType
@ BY_LOCATION
Ordering data by spatial location.
GlobalSparsityPattern computeSparsityPattern(LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
Computes a sparsity pattern for the given inputs.
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor mapFunction)
IntegrationPointMetaData getIntegrationPointMetaData(MeshLib::Mesh const &mesh, std::string const &name)
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
static void executeMemberOnDereferenced(Method method, Container const &container, Args &&... args)
virtual std::vector< double > getSigma() const =0
virtual std::vector< double > getEpsilon() const =0
virtual std::vector< double > getEpsilonMechanical() const =0