OGS 6.2.0-97-g4a610c866
ThermoMechanicsProcess.cpp
Go to the documentation of this file.
1 
10 #include "ThermoMechanicsProcess.h"
11 
12 #include <cassert>
13 
14 #include "BaseLib/Functional.h"
17 
18 #include "ThermoMechanicsFEM.h"
19 
20 namespace ProcessLib
21 {
22 namespace ThermoMechanics
23 {
24 template <int DisplacementDim>
26  MeshLib::Mesh& mesh,
27  std::unique_ptr<ProcessLib::AbstractJacobianAssembler>&& jacobian_assembler,
28  std::vector<std::unique_ptr<ParameterLib::ParameterBase>> const& parameters,
29  unsigned const integration_order,
30  std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>&&
31  process_variables,
33  SecondaryVariableCollection&& secondary_variables,
34  NumLib::NamedFunctionCaller&& named_function_caller,
35  bool const use_monolithic_scheme)
36  : Process(mesh, std::move(jacobian_assembler), parameters,
37  integration_order, std::move(process_variables),
38  std::move(secondary_variables), std::move(named_function_caller),
39  use_monolithic_scheme),
40  _process_data(std::move(process_data))
41 {
42  _nodal_forces = MeshLib::getOrCreateMeshProperty<double>(
43  mesh, "NodalForces", MeshLib::MeshItemType::Node, DisplacementDim);
44 
45  _heat_flux = MeshLib::getOrCreateMeshProperty<double>(
46  mesh, "HeatFlux", MeshLib::MeshItemType::Node, 1);
47 
48  _integration_point_writer.emplace_back(
49  std::make_unique<KelvinVectorIntegrationPointWriter>(
50  "sigma_ip",
51  static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
52  2 /*integration order*/, [this]() {
53  // Result containing integration point data for each local
54  // assembler.
55  std::vector<std::vector<double>> result;
56  result.resize(_local_assemblers.size());
57 
58  for (std::size_t i = 0; i < _local_assemblers.size(); ++i)
59  {
60  auto const& local_asm = *_local_assemblers[i];
61 
62  result[i] = local_asm.getSigma();
63  }
64 
65  return result;
66  }));
67 
68  _integration_point_writer.emplace_back(
69  std::make_unique<KelvinVectorIntegrationPointWriter>(
70  "epsilon_ip",
71  static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
72  2 /*integration order*/, [this]() {
73  // Result containing integration point data for each local
74  // assembler.
75  std::vector<std::vector<double>> result;
76  result.resize(_local_assemblers.size());
77 
78  for (std::size_t i = 0; i < _local_assemblers.size(); ++i)
79  {
80  auto const& local_asm = *_local_assemblers[i];
81 
82  result[i] = local_asm.getEpsilon();
83  }
84 
85  return result;
86  }));
87 
88  _integration_point_writer.emplace_back(
89  std::make_unique<KelvinVectorIntegrationPointWriter>(
90  "epsilon_m_ip",
91  static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
92  2 /*integration order*/, [this]() {
93  // Result containing integration point data for each local
94  // assembler.
95  std::vector<std::vector<double>> result;
96  result.resize(_local_assemblers.size());
97 
98  for (std::size_t i = 0; i < _local_assemblers.size(); ++i)
99  {
100  auto const& local_asm = *_local_assemblers[i];
101 
102  result[i] = local_asm.getEpsilonMechanical();
103  }
104 
105  return result;
106  }));
107 }
108 
109 template <int DisplacementDim>
111 {
112  return false;
113 }
114 
115 template <int DisplacementDim>
117  NumLib::LocalToGlobalIndexMap const& dof_table,
118  MeshLib::Mesh const& mesh,
119  unsigned const integration_order)
120 {
122  DisplacementDim, ThermoMechanicsLocalAssembler>(
123  mesh.getElements(), dof_table, _local_assemblers,
124  mesh.isAxiallySymmetric(), integration_order, _process_data);
125 
126  // TODO move the two data members somewhere else.
127  // for extrapolation of secondary variables
128  std::vector<MeshLib::MeshSubset> all_mesh_subsets_single_component{
132  std::move(all_mesh_subsets_single_component),
133  // by location order is needed for output
135 
137  "sigma",
140  DisplacementDim>::RowsAtCompileTime,
141  getExtrapolator(), _local_assemblers,
143 
145  "epsilon",
148  DisplacementDim>::RowsAtCompileTime,
149  getExtrapolator(), _local_assemblers,
151 
152  // Set initial conditions for integration point data.
153  for (auto const& ip_writer : _integration_point_writer)
154  {
155  // Find the mesh property with integration point writer's name.
156  auto const& name = ip_writer->name();
157  if (!mesh.getProperties().existsPropertyVector<double>(name))
158  {
159  continue;
160  }
161  auto const& mesh_property =
162  *mesh.getProperties().template getPropertyVector<double>(name);
163 
164  // The mesh property must be defined on integration points.
165  if (mesh_property.getMeshItemType() !=
167  {
168  continue;
169  }
170 
171  auto const ip_meta_data = getIntegrationPointMetaData(mesh, name);
172 
173  // Check the number of components.
174  if (ip_meta_data.n_components != mesh_property.getNumberOfComponents())
175  {
176  OGS_FATAL(
177  "Different number of components in meta data (%d) than in "
178  "the integration point field data for '%s': %d.",
179  ip_meta_data.n_components, name.c_str(),
180  mesh_property.getNumberOfComponents());
181  }
182 
183  // Now we have a properly named vtk's field data array and the
184  // corresponding meta data.
185  std::size_t position = 0;
186  for (auto& local_asm : _local_assemblers)
187  {
188  std::size_t const integration_points_read =
189  local_asm->setIPDataInitialConditions(
190  name, &mesh_property[position],
191  ip_meta_data.integration_order);
192  if (integration_points_read == 0)
193  {
194  OGS_FATAL(
195  "No integration points read in the integration point "
196  "initial conditions set function.");
197  }
198  position += integration_points_read * ip_meta_data.n_components;
199  }
200  }
201 }
202 
203 template <int DisplacementDim>
205  const double t, GlobalVector const& x, GlobalMatrix& M, GlobalMatrix& K,
206  GlobalVector& b)
207 {
208  DBUG("Assemble ThermoMechanicsProcess.");
209 
210  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
211  dof_table = {std::ref(*_local_to_global_index_map)};
212  const int process_id =
214  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
215 
216  // Call global assembler for each local assembly item.
219  pv.getActiveElementIDs(), dof_table, t, x, M, K, b,
221 }
222 
223 template <int DisplacementDim>
225  assembleWithJacobianConcreteProcess(const double t, GlobalVector const& x,
226  GlobalVector const& xdot,
227  const double dxdot_dx,
228  const double dx_dx, GlobalMatrix& M,
229  GlobalMatrix& K, GlobalVector& b,
230  GlobalMatrix& Jac)
231 {
232  DBUG("AssembleJacobian ThermoMechanicsProcess.");
233 
234  std::vector<std::reference_wrapper<NumLib::LocalToGlobalIndexMap>>
235  dof_table = {std::ref(*_local_to_global_index_map)};
236  const int process_id =
238  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
239 
240  // Call global assembler for each local assembly item.
243  _local_assemblers, pv.getActiveElementIDs(), dof_table, t, x,
244  xdot, dxdot_dx, dx_dx, M, K, b, Jac, _coupled_solutions);
245 
246  // TODO (naumov): Refactor the copy rhs part. This is copy from HM.
247  auto copyRhs = [&](int const variable_id, auto& output_vector) {
249  {
250  transformVariableFromGlobalVector(b, variable_id, dof_table[0],
251  output_vector,
252  std::negate<double>());
253  }
254  else
255  {
257  b, 0, dof_table[_coupled_solutions->process_id], output_vector,
258  std::negate<double>());
259  }
260  };
262  {
263  copyRhs(0, *_heat_flux);
264  }
266  {
267  copyRhs(1, *_nodal_forces);
268  }
269 }
270 
271 template <int DisplacementDim>
273  GlobalVector const& x, double const t, double const dt,
274  const int process_id)
275 {
276  DBUG("PreTimestep ThermoMechanicsProcess.");
277 
278  _process_data.dt = dt;
279  _process_data.t = t;
280 
281  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
282 
286  dt);
287 }
288 
289 template <int DisplacementDim>
291  GlobalVector const& x, const double /*t*/, const double /*delta_t*/,
292  int const process_id)
293 {
294  DBUG("PostTimestep ThermoMechanicsProcess.");
295 
296  ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];
297 
302 }
303 
304 template class ThermoMechanicsProcess<2>;
305 template class ThermoMechanicsProcess<3>;
306 
307 } // namespace ThermoMechanics
308 } // namespace ProcessLib
ThermoMechanicsProcess(MeshLib::Mesh &mesh, std::unique_ptr< ProcessLib::AbstractJacobianAssembler > &&jacobian_assembler, std::vector< std::unique_ptr< ParameterLib::ParameterBase >> const &parameters, unsigned const integration_order, std::vector< std::vector< std::reference_wrapper< ProcessVariable >>> &&process_variables, ThermoMechanicsProcessData< DisplacementDim > &&process_data, SecondaryVariableCollection &&secondary_variables, NumLib::NamedFunctionCaller &&named_function_caller, bool const use_monolithic_scheme)
IntegrationPointMetaData getIntegrationPointMetaData(MeshLib::Mesh const &mesh, std::string const &name)
void assembleWithJacobianConcreteProcess(const double t, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac) override
std::unique_ptr< MeshLib::MeshSubset const > _mesh_subset_all_nodes
Definition: Process.h:263
static void executeSelectedMemberOnDereferenced(Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
virtual std::vector< double > const & getIntPtSigma(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const =0
std::vector< std::reference_wrapper< ProcessVariable > > const & getProcessVariables(const int process_id) const
Definition: Process.h:124
SecondaryVariableFunctions makeExtrapolator(const unsigned num_components, NumLib::Extrapolator &extrapolator, LocalAssemblerCollection const &local_assemblers, typename NumLib::ExtrapolatableLocalAssemblerCollection< LocalAssemblerCollection >::IntegrationPointValuesMethod integration_point_values_method)
void assembleConcreteProcess(const double t, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b) override
std::vector< std::unique_ptr< ThermoMechanicsLocalAssemblerInterface > > _local_assemblers
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map_single_component
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, GlobalVector const &x, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
MeshLib::Properties & getProperties()
Definition: Mesh.h:134
void createLocalAssemblers(std::vector< MeshLib::Element *> const &mesh_elements, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< std::unique_ptr< LocalAssemblerInterface >> &local_assemblers, ExtraCtorArgs &&... extra_ctor_args)
NumLib::Extrapolator & getExtrapolator() const
Definition: Process.h:151
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, GlobalVector const &x, GlobalVector const &xdot, const double dxdot_dx, const double dx_dx, GlobalMatrix &M, GlobalMatrix &K, GlobalVector &b, GlobalMatrix &Jac, CoupledSolutionsForStaggeredScheme const *const cpl_xs)
bool isAxiallySymmetric() const
Definition: Mesh.h:137
std::unique_ptr< NumLib::LocalToGlobalIndexMap > _local_to_global_index_map
Definition: Process.h:265
Builds expression trees of named functions dynamically at runtime.
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)
static void executeSelectedMemberDereferenced(Object &object, Method method, Container const &container, std::vector< std::size_t > const &active_container_ids, Args &&... args)
SecondaryVariableCollection _secondary_variables
Definition: Process.h:267
void postTimestepConcreteProcess(GlobalVector const &x, const double t, const double delta_t, int const process_id) override
ThermoMechanicsProcessData< DisplacementDim > _process_data
void preTimestepConcreteProcess(GlobalVector const &x, double const t, double const dt, const int process_id) override
virtual std::vector< double > const & getIntPtEpsilon(const double, GlobalVector const &, NumLib::LocalToGlobalIndexMap const &, std::vector< double > &cache) const =0
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
bool existsPropertyVector(std::string const &name) const
Definition: Properties.h:79
virtual void postTimestep(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x)
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)
Definition: DOFTableUtil.h:59
CoupledSolutionsForStaggeredScheme * _coupled_solutions
Definition: Process.h:280
Handles configuration of several secondary variables from the project file.
void initializeConcreteProcess(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh, unsigned const integration_order) override
Process specific initialization called by initialize().
Ordering data by spatial location.
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
void addSecondaryVariable(std::string const &internal_name, SecondaryVariableFunctions &&fcts)
std::vector< std::unique_ptr< IntegrationPointWriter > > _integration_point_writer
Definition: Process.h:291
const bool _use_monolithic_scheme
Definition: Process.h:276
VectorMatrixAssembler _global_assembler
Definition: Process.h:274
Eigen::Matrix< double, KelvinVectorDimensions< DisplacementDim >::value, 1, Eigen::ColMajor > KelvinVectorType
Definition: KelvinVector.h:49
std::vector< std::size_t > & getActiveElementIDs() const