v0.8.23
Files | Classes | Functions
Thermal element
Collaboration diagram for Thermal element:

Files

file  ThermalElement.cpp
 
file  ThermalElement.hpp
 Operators and data structures for thermal analysis.
 
file  ThermalStressElement.hpp
 Implementation of thermal stresses element.
 
file  thermal_steady.cpp
 Example of steady thermal analysis.
 
file  thermal_unsteady.cpp
 Example of thermal unsteady analyze.
 

Classes

struct  MixTransport::MixTransportElement::BlockData
 data for evaluation of het conductivity and heat capacity elements More...
 
struct  ThermalElement
 structure grouping operators and data used for thermal problemsIn order to assemble matrices and right hand vectors, the loops over elements, entities within the element and finally loop over integration points are executed. More...
 
struct  ThermalElement::TimeSeriesMonitor
 TS monitore it records temperature at time steps. More...
 
struct  ThermalElement::UpdateAndControl
 this calass is to control time steppingIt is used to save data for temperature rate vector to MoFEM field. More...
 
struct  ThermalElement::OpConvectionLhs
 
struct  ThermalElement::OpConvectionRhs
 operator to calculate convection therms on body surface and assemble to rhs of equations More...
 
struct  ThermalElement::OpRadiationRhs
 operator to calculate radiation therms on body surface and assemble to rhs of transient equations(Residual Vector) More...
 
struct  ThermalElement::OpRadiationLhs
 
struct  ThermalElement::OpHeatFlux
 operator for calculate heat flux and assemble to right hand side More...
 
struct  ThermalElement::OpHeatCapacityLhs
 operator to calculate left hand side of heat capacity terms More...
 
struct  ThermalElement::OpHeatCapacityRhs
 operator to calculate right hand side of heat capacity terms More...
 
struct  ThermalElement::OpThermalLhs
 
struct  ThermalElement::OpThermalRhs
 
struct  ThermalElement::OpGetTetRateAtGaussPts
 operator to calculate temperature rate at Gauss pts More...
 
struct  ThermalElement::OpGetTriTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetTetTemperatureAtGaussPts
 operator to calculate temperature at Gauss pts More...
 
struct  ThermalElement::OpGetFieldAtGaussPts< OP >
 operator to calculate temperature and rate of temperature at Gauss points More...
 
struct  ThermalElement::OpGetGradAtGaussPts
 operator to calculate temperature gradient at Gauss points More...
 
struct  ThermalElement::CommonData
 common data used by volume elements More...
 
struct  ThermalElement::RadiationData
 data for radiation More...
 
struct  ThermalElement::ConvectionData
 data for convection More...
 
struct  ThermalElement::FluxData
 data for calculation heat flux More...
 
struct  ThermalElement::BlockData
 data for calculation heat conductivity and heat capacity elements More...
 
struct  ThermalElement::MyTriFE
 define surface element More...
 
struct  ThermalElement::MyVolumeFE
 definition of volume element More...
 
struct  ThermalStressElement
 Implentation of thermal stress element. More...
 
struct  ThermalStressElement::OpThermalStressRhs
 
struct  ThermalStressElement::OpGetTemperatureAtGaussPts
 
struct  ThermalStressElement::CommonData
 
struct  ThermalStressElement::BlockData
 
struct  ThermalStressElement::MyVolumeFE
 

Functions

MoFEMErrorCode ThermalElement::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 More...
 
MoFEMErrorCode ThermalElement::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 block set with name HEAT_FLUX. More...
 
MoFEMErrorCode ThermalElement::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 block set with name CONVECTION. More...
 
MoFEMErrorCode ThermalElement::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. Alternatively uses block set with name RADIATION. More...
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementRhsOperators (string field_name, Vec &F)
 this function is used in case of stationary problem to set elements for rhs More...
 
MoFEMErrorCode ThermalElement::setThermalFiniteElementLhsOperators (string field_name, Mat A)
 this function is used in case of stationary heat conductivity problem for lhs More...
 
MoFEMErrorCode ThermalElement::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 More...
 
MoFEMErrorCode ThermalElement::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 More...
 
MoFEMErrorCode ThermalElement::setTimeSteppingProblem (TsCtx &ts_ctx, 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 More...
 

Detailed Description

Function Documentation

◆ addThermalConvectionElement()

MoFEMErrorCode ThermalElement::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 block set with name CONVECTION.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 572 of file ThermalElement.cpp.

573  {
575 
576  CHKERR mField.add_finite_element("THERMAL_CONVECTION_FE", MF_ZERO);
577  CHKERR mField.modify_finite_element_add_field_row("THERMAL_CONVECTION_FE",
578  field_name);
579  CHKERR mField.modify_finite_element_add_field_col("THERMAL_CONVECTION_FE",
580  field_name);
581  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
582  field_name);
583  if (mField.check_field(mesh_nodals_positions)) {
584  CHKERR mField.modify_finite_element_add_field_data("THERMAL_CONVECTION_FE",
585  mesh_nodals_positions);
586  }
587 
588  // this is alternative method for setting boundary conditions, to bypass bu
589  // in cubit file reader. not elegant, but good enough
591  if (it->getName().compare(0, 10, "CONVECTION") == 0) {
592 
593  std::vector<double> data;
594  CHKERR it->getAttributes(data);
595  if (data.size() != 2) {
596  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
597  }
598  setOfConvection[it->getMeshsetId()].cOnvection = data[0];
599  setOfConvection[it->getMeshsetId()].tEmperature = data[1];
600  CHKERR mField.get_moab().get_entities_by_type(
601  it->meshset, MBTRI, setOfConvection[it->getMeshsetId()].tRis, true);
603  setOfConvection[it->getMeshsetId()].tRis, MBTRI,
604  "THERMAL_CONVECTION_FE");
605  }
606  }
607 
609 }
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:500
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:507
std::map< int, ConvectionData > setOfConvection
MoFEM::Interface & mField
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:595
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use

◆ addThermalElements()

MoFEMErrorCode ThermalElement::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

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 480 of file ThermalElement.cpp.

481  {
483 
484  CHKERR mField.add_finite_element("THERMAL_FE", MF_ZERO);
485  CHKERR mField.modify_finite_element_add_field_row("THERMAL_FE", field_name);
486  CHKERR mField.modify_finite_element_add_field_col("THERMAL_FE", field_name);
487  CHKERR mField.modify_finite_element_add_field_data("THERMAL_FE", field_name);
488  if (mField.check_field(mesh_nodals_positions)) {
490  mesh_nodals_positions);
491  }
492 
493  // takes skin of block of entities
494  // Skinner skin(&mField.get_moab());
495  // loop over all blocksets and get data which name is FluidPressure
497  mField, BLOCKSET | MAT_THERMALSET, it)) {
498 
499  Mat_Thermal temp_data;
500  ierr = it->getAttributeDataStructure(temp_data);
501 
502  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.resize(
503  3, 3); //(3X3) conductivity matrix
504  setOfBlocks[it->getMeshsetId()].cOnductivity_mat.clear();
505  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(0, 0) =
506  temp_data.data.Conductivity;
507  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(1, 1) =
508  temp_data.data.Conductivity;
509  setOfBlocks[it->getMeshsetId()].cOnductivity_mat(2, 2) =
510  temp_data.data.Conductivity;
511  // setOfBlocks[it->getMeshsetId()].cOnductivity =
512  // temp_data.data.Conductivity;
513 
514  setOfBlocks[it->getMeshsetId()].cApacity = temp_data.data.HeatCapacity;
515  CHKERR mField.get_moab().get_entities_by_type(
516  it->meshset, MBTET, setOfBlocks[it->getMeshsetId()].tEts, true);
518  setOfBlocks[it->getMeshsetId()].tEts, MBTET, "THERMAL_FE");
519  }
520 
522 }
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MoFEM::Interface & mField
Thermal material data structure.
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
block name is "MAT_THERMAL"
Definition: definitions.h:223
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define CHKERR
Inline error check.
Definition: definitions.h:595
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ addThermalFluxElement()

MoFEMErrorCode ThermalElement::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 block set with name HEAT_FLUX.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 525 of file ThermalElement.cpp.

526  {
528 
529  CHKERR mField.add_finite_element("THERMAL_FLUX_FE", MF_ZERO);
531  field_name);
533  field_name);
535  field_name);
536  if (mField.check_field(mesh_nodals_positions)) {
538  mesh_nodals_positions);
539  }
540 
542  it)) {
543  CHKERR it->getBcDataStructure(setOfFluxes[it->getMeshsetId()].dAta);
544  CHKERR mField.get_moab().get_entities_by_type(
545  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
547  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
548  }
549 
550  // this is alternative method for setting boundary conditions, to bypass bu
551  // in cubit file reader. not elegant, but good enough
553  if (it->getName().compare(0, 9, "HEAT_FLUX") == 0) {
554  std::vector<double> data;
555  CHKERR it->getAttributes(data);
556  if (data.size() != 1) {
557  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
558  }
559  strcpy(setOfFluxes[it->getMeshsetId()].dAta.data.name, "HeatFlu");
560  setOfFluxes[it->getMeshsetId()].dAta.data.flag1 = 1;
561  setOfFluxes[it->getMeshsetId()].dAta.data.value1 = data[0];
562  CHKERR mField.get_moab().get_entities_by_type(
563  it->meshset, MBTRI, setOfFluxes[it->getMeshsetId()].tRis, true);
565  setOfFluxes[it->getMeshsetId()].tRis, MBTRI, "THERMAL_FLUX_FE");
566  }
567  }
568 
570 }
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MoFEM::Interface & mField
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:595
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ addThermalRadiationElement()

MoFEMErrorCode ThermalElement::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. Alternatively uses block set with name RADIATION.

Parameters
fieldname
nameof mesh nodal positions (if not defined nodal coordinates are used)

Definition at line 611 of file ThermalElement.cpp.

612  {
614 
615  CHKERR mField.add_finite_element("THERMAL_RADIATION_FE", MF_ZERO);
616  CHKERR mField.modify_finite_element_add_field_row("THERMAL_RADIATION_FE",
617  field_name);
618  CHKERR mField.modify_finite_element_add_field_col("THERMAL_RADIATION_FE",
619  field_name);
620  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
621  field_name);
622  if (mField.check_field(mesh_nodals_positions)) {
623  CHKERR mField.modify_finite_element_add_field_data("THERMAL_RADIATION_FE",
624  mesh_nodals_positions);
625  }
626 
627  // this is alternative method for setting boundary conditions, to bypass bu
628  // in cubit file reader. not elegant, but good enough
630  if (it->getName().compare(0, 9, "RADIATION") == 0) {
631  std::vector<double> data;
632  ierr = it->getAttributes(data);
633  if (data.size() != 3) {
634  SETERRQ(PETSC_COMM_SELF, 1, "Data inconsistency");
635  }
636  setOfRadiation[it->getMeshsetId()].sIgma = data[0];
637  setOfRadiation[it->getMeshsetId()].eMissivity = data[1];
638  setOfRadiation[it->getMeshsetId()].aMbienttEmp = data[2];
639  CHKERR mField.get_moab().get_entities_by_type(
640  it->meshset, MBTRI, setOfRadiation[it->getMeshsetId()].tRis, true);
642  setOfRadiation[it->getMeshsetId()].tRis, MBTRI,
643  "THERMAL_RADIATION_FE");
644  }
645  }
646 
648 }
virtual MoFEMErrorCode modify_finite_element_add_field_col(const std::string &fe_name, const std::string &name_row)=0
set field col which finite element use
virtual MoFEMErrorCode add_ents_to_finite_element_by_type(const EntityHandle entities, const EntityType type, const std::string &name, const bool recursive=true)=0
add entities to finite element
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MoFEM::Interface & mField
std::map< int, RadiationData > setOfRadiation
virtual MoFEMErrorCode add_finite_element(const std::string &fe_name, enum MoFEMTypes bh=MF_EXCL, int verb=DEFAULT_VERBOSITY)=0
add finite element
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
#define CHKERR
Inline error check.
Definition: definitions.h:595
virtual MoFEMErrorCode modify_finite_element_add_field_data(const std::string &fe_name, const std::string &name_filed)=0
set finite element field data
virtual bool check_field(const std::string &name) const =0
check if field is in database
virtual MoFEMErrorCode modify_finite_element_add_field_row(const std::string &fe_name, const std::string &name_row)=0
set field row which finite element use
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ setThermalFiniteElementLhsOperators()

MoFEMErrorCode ThermalElement::setThermalFiniteElementLhsOperators ( string  field_name,
Mat  A 
)

this function is used in case of stationary heat conductivity problem for lhs

Definition at line 665 of file ThermalElement.cpp.

665  {
667  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
668  for (; sit != setOfBlocks.end(); sit++) {
669  // add finite elemen
670  feLhs.getOpPtrVector().push_back(
671  new OpThermalLhs(field_name, A, sit->second, commonData));
672  }
674 }
CommonData commonData
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
MyVolumeFE feLhs

◆ setThermalFiniteElementRhsOperators()

MoFEMErrorCode ThermalElement::setThermalFiniteElementRhsOperators ( string  field_name,
Vec &  F 
)

this function is used in case of stationary problem to set elements for rhs

Definition at line 651 of file ThermalElement.cpp.

651  {
653  std::map<int, BlockData>::iterator sit = setOfBlocks.begin();
654  feRhs.getOpPtrVector().push_back(
655  new OpGetGradAtGaussPts(field_name, commonData));
656  for (; sit != setOfBlocks.end(); sit++) {
657  // add finite element
658  feRhs.getOpPtrVector().push_back(
659  new OpThermalRhs(field_name, F, sit->second, commonData));
660  }
662 }
CommonData commonData
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MyVolumeFE feRhs
cauclate right hand side for tetrahedral elements
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ setThermalFluxFiniteElementRhsOperators()

MoFEMErrorCode ThermalElement::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

Definition at line 676 of file ThermalElement.cpp.

677  {
679  bool hoGeometry = false;
680  if (mField.check_field(mesh_nodals_positions)) {
681  hoGeometry = true;
682  }
683  std::map<int, FluxData>::iterator sit = setOfFluxes.begin();
684  for (; sit != setOfFluxes.end(); sit++) {
685  // add finite element
686  feFlux.getOpPtrVector().push_back(
687  new OpHeatFlux(field_name, F, sit->second, hoGeometry));
688  }
690 }
std::map< int, FluxData > setOfFluxes
maps side set id with appropriate FluxData
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
MoFEM::Interface & mField
boost::ptr_vector< UserDataOperator > & getOpPtrVector()
Use to push back operator for row operator.
virtual bool check_field(const std::string &name) const =0
check if field is in database
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406

◆ setTimeSteppingProblem() [1/2]

MoFEMErrorCode ThermalElement::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

◆ setTimeSteppingProblem() [2/2]

MoFEMErrorCode ThermalElement::setTimeSteppingProblem ( TsCtx &  ts_ctx,
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