v0.10.0
ThermalElement.hpp
Go to the documentation of this file.
1 /** \file ThermalElement.hpp
2  \ingroup mofem_thermal_elem
3 
4  \brief Operators and data structures for thermal analysis
5 
6  Implementation of thermal element for unsteady and steady case.
7  Radiation and convection blocks implemented by Xuan Meng
8 
9 */
10 
11 /*
12  * This file is part of MoFEM.
13  * MoFEM is free software: you can redistribute it and/or modify it under
14  * the terms of the GNU Lesser General Public License as published by the
15  * Free Software Foundation, either version 3 of the License, or (at your
16  * option) any later version.
17  *
18  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
19  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
21  * License for more details.
22  *
23  * You should have received a copy of the GNU Lesser General Public
24  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
25 
26 #ifndef __THERMAL_ELEMENT_HPP
27 #define __THERMAL_ELEMENT_HPP
28 
29 /** \brief structure grouping operators and data used for thermal problems
30  * \ingroup mofem_thermal_elem
31  *
32  * In order to assemble matrices and right hand vectors, the loops over
33  * elements, entities within the element and finally loop over integration
34  * points are executed.
35  *
36  * Following implementation separate those three types of loops and to each
37  * loop attach operator.
38  *
39  */
41 
42  /// \brief definition of volume element
46 
47  /** \brief it is used to calculate nb. of Gauss integration points
48  *
49  * for more details pleas look
50  * Reference:
51  *
52  * Albert Nijenhuis, Herbert Wilf,
53  * Combinatorial Algorithms for Computers and Calculators,
54  * Second Edition,
55  * Academic Press, 1978,
56  * ISBN: 0-12-519260-6,
57  * LC: QA164.N54.
58  *
59  * More details about algorithm
60  * http://people.sc.fsu.edu/~jburkardt/cpp_src/gm_rule/gm_rule.html
61  **/
62  int getRule(int order) { return 2 * (order - 1); };
63  };
64  MyVolumeFE feRhs; ///< cauclate right hand side for tetrahedral elements
65  MyVolumeFE &getLoopFeRhs() { return feRhs; } ///< get rhs volume element
66  MyVolumeFE feLhs; //< calculate left hand side for tetrahedral elements
67  MyVolumeFE &getLoopFeLhs() { return feLhs; } ///< get lhs volume element
68 
69  /** \brief define surface element
70  *
71  * This element is used to integrate heat fluxes; convection and radiation
72  */
76  int getRule(int order) { return 2 * order; };
77  };
78 
79  MyTriFE feFlux; //< heat flux element
80  MyTriFE &getLoopFeFlux() { return feFlux; } //< get heat flux element
81 
82  MyTriFE feConvectionRhs; //< convection element
85  return feConvectionRhs;
86  } //< get convection element
88 
89  MyTriFE feRadiationRhs; //< radiation element
92  return feRadiationRhs;
93  } //< get radiation element
95 
98  : feRhs(m_field), feLhs(m_field), feFlux(m_field),
99  feConvectionRhs(m_field), feConvectionLhs(m_field),
100  feRadiationRhs(m_field), feRadiationLhs(m_field), mField(m_field) {}
101 
102  /** \brief data for calculation heat conductivity and heat capacity elements
103  * \ingroup mofem_thermal_elem
104  */
105  struct BlockData {
106  // double cOnductivity;
107  MatrixDouble cOnductivity_mat; // This is (3x3) conductivity matrix
108  double cApacity; // rou * c_p == material density multiple heat capacity
109  double initTemp; ///< initial temperature
110  Range tEts; ///< contains elements in block set
111  };
112  std::map<int, BlockData>
113  setOfBlocks; ///< maps block set id with appropriate BlockData
114 
115  /** \brief data for calculation heat flux
116  * \ingroup mofem_thermal_elem
117  */
118  struct FluxData {
119  HeatFluxCubitBcData dAta; ///< for more details look to BCMultiIndices.hpp
120  ///< to see details of HeatFluxCubitBcData
121  Range tRis; ///< surface triangles where hate flux is applied
122  };
123  std::map<int, FluxData>
124  setOfFluxes; ///< maps side set id with appropriate FluxData
125 
126  /** \brief data for convection
127  * \ingroup mofem_thermal_elem
128  */
129  struct ConvectionData {
130  double cOnvection; /*The summation of Convection coefficients*/
131  double tEmperature; /*Ambient temperature of the area contains the black
132  body */
133  Range tRis; ///< those will be on body skin, except this with contact with
134  ///< other body where temperature is applied
135  };
136  std::map<int, ConvectionData>
137  setOfConvection; //< maps block set id with appropriate data
138 
139  /** \brief data for radiation
140  * \ingroup mofem_thermal_elem
141  */
142  struct RadiationData {
143  double sIgma; /* The Stefan-Boltzmann constant*/
144  double eMissivity; /* The surface emissivity coefficients range = [0,1] */
145  // double aBsorption; /* The surface absorption coefficients */
146  double aMbienttEmp; /* The incident radiant heat flow per unit surface area;
147  or the ambient temperature of space*/
148  Range tRis; ///< those will be on body skin, except this with contact with
149  ///< other body where temperature is applied
150  };
151  std::map<int, RadiationData>
152  setOfRadiation; //< maps block set id with appropriate data
153 
154  /** \brief common data used by volume elements
155  * \ingroup mofem_thermal_elem
156  */
157  struct CommonData {
161  inline ublas::matrix_row<MatrixDouble> getGradAtGaussPts(const int gg) {
162  return ublas::matrix_row<MatrixDouble>(gradAtGaussPts, gg);
163  }
164  };
166 
167  /// \brief operator to calculate temperature gradient at Gauss points
170 
172  OpGetGradAtGaussPts(const std::string field_name, CommonData &common_data)
174  field_name, UserDataOperator::OPROW),
175  commonData(common_data) {}
176 
177  /** \brief operator calculating temperature gradients
178  *
179  * temperature gradient is calculated multiplying derivatives of shape
180  * functions by degrees of freedom.
181  */
182  MoFEMErrorCode doWork(int side, EntityType type,
184  };
185 
186  /** \brief operator to calculate temperature and rate of temperature at Gauss
187  * points \ingroup mofem_thermal_elem
188  */
189  template <typename OP>
191 
193  OpGetFieldAtGaussPts(const std::string field_name,
194  VectorDouble &field_at_gauss_pts)
195  : OP::UserDataOperator(field_name, OP::UserDataOperator::OPROW),
196  fieldAtGaussPts(field_at_gauss_pts) {}
197 
198  /** \brief operator calculating temperature and rate of temperature
199  *
200  * temperature temperature or rate of temperature is calculated multiplying
201  * shape functions by degrees of freedom
202  */
203  MoFEMErrorCode doWork(int side, EntityType type,
206  try {
207 
208  if (data.getFieldData().size() == 0)
210  int nb_dofs = data.getFieldData().size();
211  int nb_gauss_pts = data.getN().size1();
212 
213  // initialize
214  fieldAtGaussPts.resize(nb_gauss_pts);
215  if (type == MBVERTEX) {
216  // loop over shape functions on entities always start from
217  // vertices, so if nodal shape functions are processed, vector of
218  // field values is zero at initialization
219  std::fill(fieldAtGaussPts.begin(), fieldAtGaussPts.end(), 0);
220  }
221 
222  for (int gg = 0; gg < nb_gauss_pts; gg++) {
223  fieldAtGaussPts[gg] +=
224  inner_prod(data.getN(gg, nb_dofs), data.getFieldData());
225  }
226 
227  } catch (const std::exception &ex) {
228  std::ostringstream ss;
229  ss << "throw in method: " << ex.what() << std::endl;
230  SETERRQ(PETSC_COMM_SELF, MOFEM_STD_EXCEPTION_THROW, ss.str().c_str());
231  }
232 
234  }
235  };
236 
237  /** \brief operator to calculate temperature at Gauss pts
238  * \ingroup mofem_thermal_elem
239  */
241  : public OpGetFieldAtGaussPts<MoFEM::VolumeElementForcesAndSourcesCore> {
242  OpGetTetTemperatureAtGaussPts(const std::string field_name,
243  CommonData &common_data)
245  field_name, common_data.temperatureAtGaussPts) {}
246  };
247 
248  /** \brief operator to calculate temperature at Gauss pts
249  * \ingroup mofem_thermal_elem
250  */
252  : public OpGetFieldAtGaussPts<MoFEM::FaceElementForcesAndSourcesCore> {
253  OpGetTriTemperatureAtGaussPts(const std::string field_name,
254  CommonData &common_data)
256  field_name, common_data.temperatureAtGaussPts) {}
257  };
258 
259  /** \brief operator to calculate temperature rate at Gauss pts
260  * \ingroup mofem_thermal_elem
261  */
263  : public OpGetFieldAtGaussPts<MoFEM::VolumeElementForcesAndSourcesCore> {
264  OpGetTetRateAtGaussPts(const std::string field_name,
265  CommonData &common_data)
267  field_name, common_data.temperatureRateAtGaussPts) {}
268  };
269 
270  /** \biref operator to calculate right hand side of heat conductivity terms
271  * \ingroup mofem_thermal_elem
272  */
275 
278  bool useTsF;
279  OpThermalRhs(const std::string field_name, BlockData &data,
280  CommonData &common_data)
282  field_name, UserDataOperator::OPROW),
283  dAta(data), commonData(common_data), useTsF(true) {}
284 
285  Vec F;
286  OpThermalRhs(const std::string field_name, Vec _F, BlockData &data,
287  CommonData &common_data)
289  field_name, UserDataOperator::OPROW),
290  dAta(data), commonData(common_data), useTsF(false), F(_F) {}
291 
293 
294  /** \brief calculate thermal conductivity matrix
295  *
296  * F = int diffN^T k gard_T dOmega^2
297  *
298  */
299  MoFEMErrorCode doWork(int side, EntityType type,
301  };
302 
303  /** \biref operator to calculate left hand side of heat conductivity terms
304  * \ingroup mofem_thermal_elem
305  */
308 
311  bool useTsB;
312  OpThermalLhs(const std::string field_name, BlockData &data,
313  CommonData &common_data)
315  field_name, UserDataOperator::OPROWCOL),
316  dAta(data), commonData(common_data), useTsB(true) {}
317 
318  Mat A;
319  OpThermalLhs(const std::string field_name, Mat _A, BlockData &data,
320  CommonData &common_data)
322  field_name, UserDataOperator::OPROWCOL),
323  dAta(data), commonData(common_data), useTsB(false), A(_A) {}
324 
326 
327  /** \brief calculate thermal conductivity matrix
328  *
329  * K = int diffN^T k diffN^T dOmega^2
330  *
331  */
332  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
333  EntityType col_type,
336  };
337 
338  /** \brief operator to calculate right hand side of heat capacity terms
339  * \ingroup mofem_thermal_elem
340  */
343 
346  OpHeatCapacityRhs(const std::string field_name, BlockData &data,
347  CommonData &common_data)
349  field_name, UserDataOperator::OPROW),
350  dAta(data), commonData(common_data) {}
351 
353 
354  /** \brief calculate thermal conductivity matrix
355  *
356  * F = int N^T c (dT/dt) dOmega^2
357  *
358  */
359  MoFEMErrorCode doWork(int side, EntityType type,
361  };
362 
363  /** \brief operator to calculate left hand side of heat capacity terms
364  * \ingroup mofem_thermal_elem
365  */
368 
371  OpHeatCapacityLhs(const std::string field_name, BlockData &data,
372  CommonData &common_data)
374  field_name, UserDataOperator::OPROWCOL),
375  dAta(data), commonData(common_data) {}
376 
378 
379  /** \brief calculate heat capacity matrix
380  *
381  * M = int N^T c N dOmega^2
382  *
383  */
384  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
385  EntityType col_type,
388  };
389 
390  /** \brief operator for calculate heat flux and assemble to right hand side
391  * \ingroup mofem_thermal_elem
392  */
393  struct OpHeatFlux
395 
398  bool useTsF;
399  OpHeatFlux(const std::string field_name, FluxData &data,
400  bool ho_geometry = false)
402  field_name, UserDataOperator::OPROW),
403  dAta(data), hoGeometry(ho_geometry), useTsF(true) {}
404 
405  Vec F;
406  OpHeatFlux(const std::string field_name, Vec _F, FluxData &data,
407  bool ho_geometry = false)
409  field_name, UserDataOperator::OPROW),
410  dAta(data), hoGeometry(ho_geometry), useTsF(false), F(_F) {}
411 
413 
414  /** \brief calculate heat flux
415  *
416  * F = int_S N^T * flux dS
417  *
418  */
419  MoFEMErrorCode doWork(int side, EntityType type,
421  };
422 
423  /**
424  * operator to calculate radiation therms on body surface and assemble to lhs
425  * of equations for the jocabian Matrix of Picard Linearization \ingroup
426  * mofem_thermal_elem
427  */
430  CommonData
431  &commonData; // get the temperature or temperature Rate from CommonData
434  bool useTsB;
435 
436  OpRadiationLhs(const std::string field_name, RadiationData &data,
437  CommonData &common_data, bool ho_geometry = false)
439  field_name, UserDataOperator::OPROWCOL),
440  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
441  useTsB(true) {}
442 
443  Mat A;
444  OpRadiationLhs(const std::string field_name, Mat _A, RadiationData &data,
445  CommonData &common_data, bool ho_geometry = false)
447  field_name, UserDataOperator::OPROWCOL),
448  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
449  useTsB(false), A(_A) {}
450 
452 
453  /** \brief calculate thermal radiation term in the lhs of equations(Tangent
454  * Matrix) for transient Thermal Problem
455  *
456  * K = intS 4* N^T* sIgma* eMissivity* N* T^3 dS (Reference _ see Finite
457  * Element Simulation of Heat Transfer by jean-Michel Bergheau)
458  */
459  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
460  EntityType col_type,
463  };
464 
465  /** \brief operator to calculate radiation therms on body surface and assemble
466  * to rhs of transient equations(Residual Vector) \ingroup mofem_thermal_elem
467  */
470 
471  CommonData
472  &commonData; // get the temperature or temperature Rate from CommonData
475  bool useTsF;
476  OpRadiationRhs(const std::string field_name, RadiationData &data,
477  CommonData &common_data, bool ho_geometry = false)
479  field_name, UserDataOperator::OPROW),
480  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
481  useTsF(true) {}
482 
483  Vec F;
484  OpRadiationRhs(const std::string field_name, Vec _F, RadiationData &data,
485  CommonData &common_data, bool ho_geometry = false)
487  field_name, UserDataOperator::OPROW),
488  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
489  useTsF(false), F(_F) {}
490 
492 
493  /** \brief calculate Transient Radiation condition on the right hand side
494  *residual
495  *
496  * R=int_S N^T * sIgma * eMissivity * (Ta^4 -Ts^4) dS
497  **/
498  MoFEMErrorCode doWork(int side, EntityType type,
500  };
501 
502  /** \brief operator to calculate convection therms on body surface and
503  * assemble to rhs of equations \ingroup mofem_thermal_elem
504  */
507 
508  CommonData
509  &commonData; // get the temperature or temperature Rate from CommonData
512  bool useTsF;
513  OpConvectionRhs(const std::string field_name, ConvectionData &data,
514  CommonData &common_data, bool ho_geometry = false)
516  field_name, UserDataOperator::OPROW),
517  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
518  useTsF(true) {}
519 
520  Vec F;
521  OpConvectionRhs(const std::string field_name, Vec _F, ConvectionData &data,
522  CommonData &common_data, bool ho_geometry = false)
524  field_name, UserDataOperator::OPROW),
525  commonData(common_data), dAta(data), hoGeometry(ho_geometry),
526  useTsF(false), F(_F) {}
527 
529 
530  /** brief calculate Convection condition on the right hand side
531  * R=int_S N^T*alpha*N_f dS **/
532 
533  MoFEMErrorCode doWork(int side, EntityType type,
535  };
536 
537  /// \biref operator to calculate convection therms on body surface and
538  /// assemble to lhs of equations
541 
544  bool useTsB;
545 
546  OpConvectionLhs(const std::string field_name, ConvectionData &data,
547  bool ho_geometry = false)
549  field_name, UserDataOperator::OPROWCOL),
550  dAta(data), hoGeometry(ho_geometry), useTsB(true) {}
551 
552  Mat A;
553  OpConvectionLhs(const std::string field_name, Mat _A, ConvectionData &data,
554  bool ho_geometry = false)
556  field_name, UserDataOperator::OPROWCOL),
557  dAta(data), hoGeometry(ho_geometry), useTsB(false), A(_A) {}
558 
560  /** \brief calculate thermal convection term in the lhs of equations
561  *
562  * K = intS N^T alpha N dS
563  */
564  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
565  EntityType col_type,
568  };
569 
570  /** \brief this calass is to control time stepping
571  * \ingroup mofem_thermal_elem
572  *
573  * It is used to save data for temperature rate vector to MoFEM field.
574  */
575  struct UpdateAndControl : public FEMethod {
576 
578  const std::string tempName;
579  const std::string rateName;
580 
581  UpdateAndControl(MoFEM::Interface &m_field, const std::string temp_name,
582  const std::string rate_name)
583  : mField(m_field), tempName(temp_name), rateName(rate_name) {}
584 
587  };
588 
589  /** \brief TS monitore it records temperature at time steps
590  * \ingroup mofem_thermal_elem
591  */
592  struct TimeSeriesMonitor : public FEMethod {
593 
595  const std::string seriesName;
596  const std::string tempName;
598 
599  TimeSeriesMonitor(MoFEM::Interface &m_field, const std::string series_name,
600  const std::string temp_name)
601  : mField(m_field), seriesName(series_name), tempName(temp_name) {
602  mask.set();
603  }
604 
606  };
607 
608  /** \brief add thermal element on tets
609  * \ingroup mofem_thermal_elem
610  *
611  * It get data from block set and define element in moab
612  *w
613  * \param field name
614  * \param name of mesh nodal positions (if not defined nodal coordinates are
615  *used)
616  */
618  const std::string field_name,
619  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
620 
621  /** \brief add heat flux element
622  * \ingroup mofem_thermal_elem
623  *
624  * It get data from heat flux set and define element in moab. Alternatively
625  * uses block set with name HEAT_FLUX.
626  *
627  * \param field name
628  * \param name of mesh nodal positions (if not defined nodal coordinates are
629  * used)
630  */
632  const std::string field_name,
633  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
634 
635  /** \brief add convection element
636  * \ingroup mofem_thermal_elem
637  *
638  * It get data from convection set and define element in moab. Alternatively
639  * uses block set with name CONVECTION.
640  *
641  * \param field name
642  * \param name of mesh nodal positions (if not defined nodal coordinates are
643  * used)
644  */
646  const std::string field_name,
647  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
648 
649  /** \brief add Non-linear Radiation element
650  * \ingroup mofem_thermal_elem
651  *
652  * It get data from Radiation set and define element in moab. Alternatively
653  * uses block set with name RADIATION.
654  *
655  * \param field name
656  * \param name of mesh nodal positions (if not defined nodal coordinates are
657  * used)
658  */
660  const std::string field_name,
661  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
662 
663  /** \brief this function is used in case of stationary problem to set elements
664  * for rhs \ingroup mofem_thermal_elem
665  */
666  MoFEMErrorCode setThermalFiniteElementRhsOperators(string field_name, Vec &F);
667 
668  /** \brief this function is used in case of stationary heat conductivity
669  * problem for lhs \ingroup mofem_thermal_elem
670  */
671  MoFEMErrorCode setThermalFiniteElementLhsOperators(string field_name, Mat A);
672 
673  /** \brief this function is used in case of stationary problem for heat flux
674  * terms \ingroup mofem_thermal_elem
675  */
677  string field_name, Vec &F,
678  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
679 
680  /* \brief linear Steady convection terms in lhs
681  */
683  string field_name, Vec &F,
684  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
685 
686  /* \brief linear Steady convection terms in rhs
687  */
689  string field_name, Mat A,
690  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
691 
692  /** \brief set up operators for unsteady heat flux; convection; radiation
693  * problem \ingroup mofem_thermal_elem
694  */
696  string field_name, string rate_name,
697  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
698 
699  /** \brief set up operators for unsteady heat flux; convection; radiation
700  * problem \ingroup mofem_thermal_elem
701  */
703  TsCtx &ts_ctx, string field_name, string rate_name,
704  const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
705 };
706 
707 #endif //__THERMAL_ELEMENT_HPP
708 
709 /**
710  * \defgroup mofem_thermal_elem Thermal element
711  * \ingroup user_modules
712  **/
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
OpGetFieldAtGaussPts(const std::string field_name, VectorDouble &field_at_gauss_pts)
MoFEMErrorCode setThermalConvectionFiniteElementRhsOperators(string field_name, Vec &F, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
OpConvectionRhs(const std::string field_name, ConvectionData &data, CommonData &common_data, bool ho_geometry=false)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
calculate thermal conductivity matrix
Deprecated interface functions.
MoFEMErrorCode setThermalFiniteElementLhsOperators(string field_name, Mat A)
this function is used in case of stationary heat conductivity problem for lhs
ublas::matrix_row< MatrixDouble > getGradAtGaussPts(const int gg)
VolumeElementForcesAndSourcesCoreSwitch< 0 > VolumeElementForcesAndSourcesCore
Volume finite element default.
data for calculation heat flux
int getRule(int order)
it is used to calculate nb. of Gauss integration points
OpThermalLhs(const std::string field_name, Mat _A, BlockData &data, CommonData &common_data)
TimeSeriesMonitor(MoFEM::Interface &m_field, const std::string series_name, const std::string temp_name)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
calculate thermal conductivity matrix
MyVolumeFE & getLoopFeLhs()
get lhs volume element
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
calculate thermal radiation term in the lhs of equations(Tangent Matrix) for transient Thermal Proble...
CommonData commonData
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
OpRadiationRhs(const std::string field_name, Vec _F, RadiationData &data, CommonData &common_data, bool ho_geometry=false)
OpThermalRhs(const std::string field_name, BlockData &data, CommonData &common_data)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
define surface element
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
MoFEMErrorCode addThermalFluxElement(const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
add heat flux elementIt get data from heat flux set and define element in moab. Alternatively uses bl...
operator to calculate temperature at Gauss pts
UpdateAndControl(MoFEM::Interface &m_field, const std::string temp_name, const std::string rate_name)
MoFEMErrorCode setThermalFluxFiniteElementRhsOperators(string field_name, Vec &F, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
this function is used in case of stationary problem for heat flux terms
OpRadiationLhs(const std::string field_name, Mat _A, RadiationData &data, CommonData &common_data, bool ho_geometry=false)
common data used by volume elements
data for calculation heat conductivity and heat capacity elements
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
std::map< int, ConvectionData > setOfConvection
OpConvectionLhs(const std::string field_name, Mat _A, ConvectionData &data, bool ho_geometry=false)
TS monitore it records temperature at time steps.
MyVolumeFE(MoFEM::Interface &m_field)
OpGetTetRateAtGaussPts(const std::string field_name, CommonData &common_data)
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
MyTriFE & getLoopFeConvectionLhs()
operator to calculate temperature rate at Gauss pts
MoFEMErrorCode addThermalConvectionElement(const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
add convection elementIt get data from convection set and define element in moab. Alternatively uses ...
operator to calculate left hand side of heat capacity terms
MoFEMErrorCode addThermalRadiationElement(const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
add Non-linear Radiation elementIt get data from Radiation set and define element in moab....
operator to calculate right hand side of heat capacity terms
MoFEM::Interface & mField
MoFEMErrorCode setTimeSteppingProblem(string field_name, string rate_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
set up operators for unsteady heat flux; convection; radiation problem
double initTemp
initial temperature
OpHeatCapacityLhs(const std::string field_name, BlockData &data, CommonData &common_data)
MoFEMErrorCode setThermalConvectionFiniteElementLhsOperators(string field_name, Mat A, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
ForcesAndSourcesCore::UserDataOperator UserDataOperator
MyTriFE(MoFEM::Interface &m_field)
DataForcesAndSourcesCore::EntData EntData
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
calculate heat capacity matrix
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
calculate thermal convection term in the lhs of equations
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
operator calculating temperature and rate of temperature
MyVolumeFE feRhs
cauclate right hand side for tetrahedral elements
definition of volume element
Range tEts
contains elements in block set
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
operator calculating temperature gradients
operator to calculate convection therms on body surface and assemble to rhs of equations
OpGetTetTemperatureAtGaussPts(const std::string field_name, CommonData &common_data)
MyTriFE feRadiationRhs
Range tRis
surface triangles where hate flux is applied
OpConvectionRhs(const std::string field_name, Vec _F, ConvectionData &data, CommonData &common_data, bool ho_geometry=false)
MyVolumeFE & getLoopFeRhs()
get rhs volume element
std::map< int, RadiationData > setOfRadiation
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
calculate Transient Radiation condition on the right hand side residual
MyTriFE feRadiationLhs
OpRadiationRhs(const std::string field_name, RadiationData &data, CommonData &common_data, bool ho_geometry=false)
operator to calculate temperature and rate of temperature at Gauss points
constexpr int order
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
calculate heat flux
HeatFluxCubitBcData dAta
OpGetTriTemperatureAtGaussPts(const std::string field_name, CommonData &common_data)
MyTriFE feConvectionRhs
operator to calculate temperature gradient at Gauss points
operator to calculate radiation therms on body surface and assemble to rhs of transient equations(Res...
operator for calculate heat flux and assemble to right hand side
MoFEMErrorCode setThermalFiniteElementRhsOperators(string field_name, Vec &F)
this function is used in case of stationary problem to set elements for rhs
MyTriFE & getLoopFeRadiationRhs()
OpGetGradAtGaussPts(const std::string field_name, CommonData &common_data)
OpHeatFlux(const std::string field_name, FluxData &data, bool ho_geometry=false)
MyTriFE & getLoopFeConvectionRhs()
MyTriFE & getLoopFeRadiationLhs()
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
FaceElementForcesAndSourcesCoreSwitch< 0 > FaceElementForcesAndSourcesCore
Face finite element default.
MoFEMErrorCode addThermalElements(const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
add thermal element on tetsIt get data from block set and define element in moab w
OpThermalRhs(const std::string field_name, Vec _F, BlockData &data, CommonData &common_data)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
calculate thermal conductivity matrix
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...
#define _F(n)
Definition: quad.c:25
ThermalElement(MoFEM::Interface &m_field)
MyVolumeFE feLhs
MyTriFE feConvectionLhs
VectorDouble temperatureRateAtGaussPts
OpHeatCapacityRhs(const std::string field_name, BlockData &data, CommonData &common_data)
operator to calculate temperature at Gauss pts
OpRadiationLhs(const std::string field_name, RadiationData &data, CommonData &common_data, bool ho_geometry=false)
OpHeatFlux(const std::string field_name, Vec _F, FluxData &data, bool ho_geometry=false)
this calass is to control time steppingIt is used to save data for temperature rate vector to MoFEM f...
OpThermalLhs(const std::string field_name, BlockData &data, CommonData &common_data)
structure grouping operators and data used for thermal problemsIn order to assemble matrices and righ...
MyTriFE & getLoopFeFlux()
OpConvectionLhs(const std::string field_name, ConvectionData &data, bool ho_geometry=false)